[R] using mapply to avoid loops

[David Winsemius]

On Oct 14, 2009, at 12:58 AM, Stephen Samaha wrote:

> 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)
>
>

(It would be more courteous to offer the error messages you are  
currently keeping secret.)

Could it be because you are offering a mix of vectors and scalars to  
mapply without properly segregating them? See the help page for mapply  
and the moreArgs argument.

-- 
David Winsemius, MD
Heritage Laboratories
West Hartford, CT