[R] Failure to run mcsamp() in package arm

[Michael Kubovy]

More problems. If I run
sim(fm1 <- lmer(Reaction ~ Days + (Days|Subject), sleepstudy))
from the lmer() help page.

I get the error
Error in mvrnorm(n.sims, bhat[j, ], V.beta) :
	'Sigma' is not positive definite

On Mar 7, 2007, at 1:30 PM, Michael Kubovy wrote:

> Dear r-helpers,
>
> I can run the examples on the mcsamp help page. For example:
> ****************************************
>> M1 <- lmer (y1 ~ x + (1|group))
>> (M1.sim <- mcsamp (M1))
> fit using lmer,
> 3 chains, each with 1000 iterations (first 500 discarded)
> n.sims = 1500 iterations saved
>                            mean  sd 2.5%  25%  50%  75% 97.5% Rhat  
> n.eff
> beta.(Intercept)           0.1 0.7 -1.2 -0.3  0.1  0.5   1.4  1.0   
> 1500
> beta.x                     2.5 0.4  1.7  2.2  2.5  2.7   3.2  1.0   
> 1500
> sigma.y                    3.8 0.3  3.3  3.6  3.7  3.9   4.3   
> 1.0    61
> sigma.grop.(In)            1.5 0.8  0.0  1.0  1.4  1.9   3.3   
> 1.4    12
> eta.group.(Intercept)[1]   0.0 1.0 -2.1 -0.5  0.0  0.6   2.0  1.0   
> 1500
> eta.group.(Intercept)[2]   1.0 1.1 -0.9  0.2  0.9  1.7   3.4   
> 1.0    59
> eta.group.(Intercept)[3]  -1.3 1.2 -4.0 -2.0 -1.3 -0.4   0.5   
> 1.0    66
> eta.group.(Intercept)[4]   1.3 1.1 -0.6  0.4  1.1  2.0   3.7   
> 1.1    43
> eta.group.(Intercept)[5]  -0.7 1.0 -3.0 -1.4 -0.6  0.0   1.2  1.0    
> 120
> eta.group.(Intercept)[6]   1.5 1.2 -0.3  0.6  1.4  2.2   4.0   
> 1.0    49
> eta.group.(Intercept)[7]   0.3 1.0 -1.7 -0.3  0.1  0.8   2.5  1.0    
> 440
> eta.group.(Intercept)[8]  -1.6 1.2 -4.0 -2.4 -1.5 -0.6   0.3   
> 1.1    41
> eta.group.(Intercept)[9]   0.4 1.0 -1.6 -0.2  0.2  0.9   2.7  1.0    
> 180
> eta.group.(Intercept)[10] -1.0 1.1 -3.3 -1.6 -0.9 -0.2   0.8   
> 1.0    86
>
> For each parameter, n.eff is a crude measure of effective sample size,
> and Rhat is the potential scale reduction factor (at convergence,
> Rhat=1).
> ****************************************
> But when I try to do this with my own data I get an error:
> ****************************************
>> display(e7.lmer2)
> lmer(formula = baLO ~ I(baRatio - 0.985) + delta + (1 + I(baRatio -
> 0.985) + delta | subject), data = e7)
>                     coef.est coef.se
> (Intercept)        -0.19     0.06
> I(baRatio - 0.985) -4.95     0.74
> delta               0.41     0.06
> Error terms:
> Groups   Name               Std.Dev. Corr
> subject  (Intercept)        0.13
>            I(baRatio - 0.985) 2.57      0.45
>            delta              0.22     -0.12 -0.94
> Residual                    0.39
> number of obs: 494, groups: subject, 13
> deviance = 551.4
>
>> e7.sim <- mcsamp(e7.lmer2)
> Error in as.bugs.array(sims, program = "lmer", n.iter = n.iter,
> n.burnin = n.burnin,  :
> 	error in parameter sigma. in parameters.to.save
> ****************************************
> I would appreciate a pointer to what the problem might be.
> _____________________________
> Professor Michael Kubovy
> University of Virginia
> Department of Psychology
> USPS:     P.O.Box 400400    Charlottesville, VA 22904-4400
> Parcels:    Room 102        Gilmer Hall
>          McCormick Road    Charlottesville, VA 22903
> Office:    B011    +1-434-982-4729
> Lab:        B019    +1-434-982-4751
> Fax:        +1-434-982-4766
> WWW:    http://www.people.virginia.edu/~mk9y/
>
> ______________________________________________
> R-help at stat.math.ethz.ch 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.