R/int.eprior.R

int.eprior <- function(sdat, g.hat, d.hat) {
    g.star <- d.star <- NULL
    r <- nrow(sdat)
    for (i in 1:r) {
        g <- g.hat[-i]
        d <- d.hat[-i]
        x <- sdat[i, !is.na(sdat[i, ])]
        n <- length(x)
        j <- numeric(n) + 1
        dat <- matrix(as.numeric(x), length(g), n, 
            byrow = T)
        resid2 <- (dat - g)^2
        sum2 <- resid2 %*% j
        LH <- 1/(2 * pi * d)^(n/2) * exp(-sum2/(2 * 
            d))
        LH[LH == "NaN"] = 0
        g.star <- c(g.star, sum(g * LH)/sum(LH))
        d.star <- c(d.star, sum(d * LH)/sum(LH))
        # if(i%%1000==0){cat(i,'\n')}
    }
    adjust <- rbind(g.star, d.star)
    rownames(adjust) <- c("g.star", "d.star")
    adjust
}
BowenNCSU/xMSanalyzer documentation built on May 24, 2019, 7:36 a.m.