individual log-likelihood for the multivariate normal distribution

  • Mauricio Garnier-Villarreal


    I am working to calculate the WAIC. For this I need to calculate the log-likelihood for each row of the data. This can be done outside of JAGS in R, but I am working on calculting in the same JAGS code.

    I can reproduce the log-likelihood for univariate models, but when I try to calculate it for the multivariate normal model, the log-likelihood is very different from the deviance/-2 from JAGS. I am comparing my calculation with the deviance of JAGS, but it hasnt work

    Here is the code that I have been using, could somebody help me how to get this calculation right?

            for(i in 1:nRows){
                Y[i, 1:nVars] ~ dmnorm(Y.hat[i, 1:nVars], TauY[1:nVars,1:nVars])
                Y.hat[i,1:nVars] <- b[id[i],1,1:nVars] + age[i]*b[id[i],2,1:nVars] + (b[id[i], 3, 1:nVars]*((age[i] - b[id[i],4,1:nVars])*step(age[i]-b[id[i],4,1:nVars])))
                for(j in 1:nVars){
                dif1[i,j] <- (Y[i,j]-Y.hat[i,j])
            dift <- t(dif1)
            for(i in 1:nRows){
                log_lik[i] <- (-(nVars/2)*log(2*pi)) - (0.5*(logdet(SigmaY[,]))) - (0.5*dif1[i,]%*%inverse(SigmaY[,])%*%dift[,i]) 
        for(p in 1:nVars){
        for(j in 1:N){
        b[j, 1:params,p] ~dmnorm(mub[1:params,p], Taub[1:params, 1:params,p])
        mub[1:params,p] ~ dmnorm(pm[1:params,p], prec[1:params, 1:params])
        Taub[1:params, 1:params,p] ~ dwish(R[1:params,1:params], params)
        Sigma2b[1:params,1:params,p] <- inverse(Taub[1:params,1:params,p])
        for (i in 1 : params) {Stdb[i, p] <- sqrt(Sigma2b[i, i, p])} 
        log_lik0 <- sum(log_lik[1:nRows]) 
        dev <- -2*sum(log_lik[1:nRows]) 
        for (j in 1:nVars) {stdy[j] <- sqrt(SigmaY[j,j])}
        TauY[1:nVars,1:nVars] ~ dwish(W[1:nVars,1:nVars], nVars)
        SigmaY[1:nVars,1:nVars] <- inverse(TauY[1:nVars,1:nVars])



Get latest updates about Open Source Projects, Conferences and News.

Sign up for the SourceForge newsletter:

No, thanks