Hello, I would like to use mapply to avoid using a loop but for some reason, I can't seem to get it to work. I've included copies of my code below. The first set of code uses a loop (and it works fine), and the second set of code attempts to use mapply but I get a "subscript out of bounds" error. Any guidance would be greatly appreciated. Xj, Yj, and Wj are also lists, and s2, TAU, and GAMMA are scalars. Thank You. # THIS WORKS USING THE LOOP for (j in 1:J) { V.tilde.j <- solve((1/s2)*t(Xj[[j]])%*%Xj[[j]] + solve(TAU)) # Not singular case: if(round(det(t(Xj[[j]])%*%Xj[[j]]),8)!=0) { Beta.hat.j <- solve(t(Xj[[j]])%*%Xj[[j]])%*%t(Xj[[j]])%*%Yj[[j]] V.j <- s2*solve(t(Xj[[j]])%*%Xj[[j]]) Lambda.j <- solve(solve(V.j) + solve(TAU))%*%solve(V.j) Beta.tilde.j <- Lambda.j%*%Beta.hat.j + (diag(P) - Lambda.j)%*%Wj[[j]]%*%GAMMA } # Singular case else { Beta.tilde.j <- V.tilde.j%*%((1/s2)*t(Xj[[j]])%*%Yj[[j]] + solve(TAU)%*%Wj[[j]]%*%GAMMA) } BETA.Js[[j]] <- t(rmnorm(1, mean=as.vector(Beta.tilde.j), V.tilde.j)) } # THIS DOESN'T WORK USING MAPPLY update.betas <- function(s2,Xj,Yj,TAU,Wj,GAMMA) { V.tilde.j <- solve((1/s2)*t(Xj[[j]])%*%Xj[[j]] + solve(TAU)) # Not singular case: if(round(det(t(Xj[[j]])%*%Xj[[j]]),8)!=0) { Beta.hat.j <- solve(t(Xj[[j]])%*%Xj[[j]])%*%t(Xj[[j]])%*%Yj[[j]] V.j <- s2*solve(t(Xj[[j]])%*%Xj[[j]]) Lambda.j <- solve(solve(V.j) + solve(TAU))%*%solve(V.j) Beta.tilde.j <- Lambda.j%*%Beta.hat.j + (diag(P) - Lambda.j)%*%Wj[[j]]%*%GAMMA } # Singular case else { Beta.tilde.j <- V.tilde.j%*%((1/s2)*t(Xj[[j]])%*%Yj[[j]] + solve(TAU)%*%Wj[[j]]%*%GAMMA) } BETA.Js[[j]] <- t(rmnorm(1, mean=as.vector(Beta.tilde.j), V.tilde.j)) return(Beta.tilde.j) } BETA.Js <- mapply(update.betas,s2,Xj,Yj,TAU,Wj,GAMMA,SIMPLIFY=FALSE)
[[alternative HTML version deleted]]
______________________________________________ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.