[R-sig-ME] Regarding mcmcsamp in the most recent alpha release

[Douglas Bates]

The first thing to notice is that the estimated correlation of the
random effects is -1.000 which, as Simon indicates, calls into
question the model that you are fitting.  Also, does Trial vary within
Patient?  You can only fit a random Trial:Patient interaction if Trial
varies within Patient.

You only have 10 different Patients so it is unlikely that you will be
able to estimate many variance parameters for the random effects. It
may not be obvious but the model that you have fit is somewhat
complicated.  I would start with simpler models.  First fit a model
with the random effect term (1|Patient).  The next model to attempt to
fit is (1|Patient/Trial), which is what I think that Simon meant to
write.  This provides a random effect for Patient and a random effect
for the Patient:Trial interaction.  It is easier to fit this model
than to fit the other version of a "random interaction" which is
(1+Trial|Patient).

Compare the AIC for those models with the AIC for the model that you
fit.  (To be more confident of the results you should fit all three
models with method = "ML")

That error message can occur even when the model is well-defined.  It
can be caused by the MCMC sampler getting stuck in regions with low,
but very flat, posterior probability. I have seen it on samples like

mcmcsamp(fm1 <- lmer(Reaction ~ Days + (Days|Subject), sleepstudy), 10000)

but it is a low probability event for such models and I haven't found
a seed where I can reproduce the problem repeatably.


On Sat, Mar 29, 2008 at 8:21 AM, Simon Blomberg <s.blomberg1 at uq.edu.au> wrote:
> Why do you have a random coefficient for Trial? Perhaps (1|Trial/Patient) might work?

>  -----Original Message-----
>  From: r-sig-mixed-models-bounces at r-project.org on behalf of Dave Warren
>  Sent: Sat 29/03/2008 2:39 PM
>  To: r-sig-mixed-models at r-project.org
>  Subject: [R-sig-ME] Regarding mcmcsamp  in the most recent alpha release
>
>  Hi all,
>
>     I've been trying out the most recent alpha release and keep running
>  into the following exception when I try using mcmcsamp.
>
>   > path.lmer = lmer( log( PathLength ) ~ Trial * Status + (Trial |
>  Patient), data = path_agg )
>   > summary( path.lmer )
>  Linear mixed model fit by REML
>  Formula: log(PathLength) ~ Trial * Status + (Trial | Patient)
>    Data: path_agg
>    AIC   BIC logLik deviance REMLdev
>   640.4 674.1 -312.2    594.9   624.4
>  Random effects:
>   Groups   Name        Variance   Std.Dev.  Corr
>   Patient  (Intercept) 1.3567e-03 0.0368331
>           Trial       2.1394e-05 0.0046254 -1.000
>   Residual             1.9524e-01 0.4418560
>  Number of obs: 502, groups: Patient, 10
>
>  Fixed effects:
>                     Estimate Std. Error t value
>  (Intercept)         2.995689   0.060844   49.24
>  Trial              -0.001977   0.002406   -0.82
>  StatusNormal        0.245692   0.084061    2.92
>  Trial:StatusNormal  0.008709   0.003382    2.58
>
>  Correlation of Fixed Effects:
>             (Intr) Trial  SttsNr
>  Trial       -0.661
>  StatusNorml -0.724  0.479
>  Trl:SttsNrm  0.470 -0.711 -0.660
>
>   > a = mcmcsamp( path.lmer, 10000 )
>  Error in .local(object, n, verbose, ...) :
>   crossproduct matrix 1 is not positive definite
>
>     Have I broken something?  Is my model revealed to be invalid?
>  Thanks for any thoughts,
>
>  Dave
>
>  _______________________________________________
>  R-sig-mixed-models at r-project.org mailing list
>  https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>
>
>         [[alternative HTML version deleted]]
>
>
>
>  _______________________________________________
>  R-sig-mixed-models at r-project.org mailing list
>  https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>