[R] convergence error code in mixed effects models
Ilona Leyer
ileyer at yahoo.de
Mon Dec 17 09:31:01 CET 2007
Thank you for your helpful answer. The only thing
obout what I still wonder is that some of the analyses
work with R some time and some time not without any
changes in the data or commands.
Your text (from both emails) is cut at the end, so
that your information about the convergence error
problem is not complete. Would you send it again?
Than you for your help!
Ilona
--- Douglas Bates <bates at stat.wisc.edu> schrieb:
> Thank you for sending the data. It is very helpful
> in understanding
> the nature of the problem.
>
> For example, in your original description of your
> study you referred
> to week as a factor, which is a completely
> reasonable term, but I
> mistakenly thought that you meant an object of class
> "factor" and that
> was why I replied that you would be estimating too
> many variances and
> covariances.
>
> I can tell you why you are having problems fitting a
> mixed-effects
> model. Strangely it is because there is too little
> variability in the
> patterns across the replicates, especially at the
> early times. The
> leaf number is discrete with a small range (all the
> leaffra
> observations in this example are 4, 5, 6, 7 or 8)
> and non-decreasing
> over time. (I assume the nature of the experiment
> is such that the
> leaf number is necessarily non-decreasing.) That
> doesn't allow for
> many patterns. I'm sure some clever person reading
> this will be able
> to tell us exactly how many different such patterns
> you could get but
> I will simply say "not many".
>
> Notice that the first 10 lines show that the leaffra
> is 4 at week 4 in
> *every* replicate. There isn't a whole lot of
> variation here for the
> random effects to model.
>
> The best place to start is with a plot of the data.
> I changed the
> levels of the rep factor to "f1"-"f5" and "s1"-"s5"
> to indicate that
> each rep is at one level of the treat. (Those who
> are playing along at
> home should be careful of the ordering of the
> original levels because
> ID10, which I now call s5, occurs between ID1 and
> ID2.)
>
> With this set of labels the cross-tabulation and
> treat should be
>
> > xtabs(~ rep + treat, leaf)
> treat
> rep pHf pHs
> f1 4 0
> f2 4 0
> f3 4 0
> f4 4 0
> f5 4 0
> s1 0 4
> s2 0 4
> s3 0 4
> s4 0 4
> s5 0 4
>
> Now look at lattice plots such as
>
> library(lattice)
> xyplot(heightfra ~ week | rep, leaf, type = c("g",
> "p", "r"), layout =
> c(5,2), aspect = 'xy', groups = treat)
>
> (I enclose PDF files of these plots for each of the
> three responses.)
> First you can see that there is very little
> variation at the low end.
> Strangely enough, this causes a problem in fitting
> mixed-effects
> models because the mle's of the variances of the
> random effects for
> the intercept will be zero. The lme function does
> not handle this
> gracefully. The lmer function from the lme4 package
> does a better job
> on this type of model.
>
> Also, note that the pattern of heightfra over time
> is not linear. It
> is consistently concave down. Thuso a mixed-effects
> model that is
> linear in week will miss much of the structure in
> the data.
>
> The point of R is to encourage you to explore your
> data rather than
> subjecting it to a "canned" analysis. You could try
> fitting a
> mixed-effects model to these data in SAS PROC MIXED
> or SPSS MIXED and
> I have no doubt that those packages would give you
> estimates (not to
> mention p-values, something that the author of lmer
> has been woefully
> negligent in not providing :-) but you probably
> won't get much of a
> hint that the model doesn't make sense. I would
> prefer to start with
> the plot and see what the data have to say.
>
>
> The technical problem with convergence in lme is
> that the mle of the
> variance of the intercept term is zero. You can see
> that if you use
> lmer from the lme4 package instead to fit the model.
> the random effect for the intercept is estimate
> On Dec 14, 2007 4:40 AM, Ilona Leyer
> <ileyer at yahoo.de> wrote:
> >
> > Here an simple example:
> >
> > rep treat heightfra leaffra leafvim
> week
> > ID1 pHf 1.54 4 4 4
> > ID2 pHf 1.49 4 4 4
> > ID3 pHf 1.57 4 5 4
> > ID4 pHf 1.48 4 4 4
> > ID5 pHf 1.57 4 4 4
> > ID6 pHs 1.29 4 5 4
> > ID7 pHs 0.97 4 5 4
> > ID8 pHs 2.06 4 4 4
> > ID9 pHs 0.88 4 4 4
> > ID10 pHs 1.47 4 4 4
> > ID1 pHf 3.53 5 6 6
> > ID2 pHf 4.08 6 6 6
> > ID3 pHf 3.89 6 6 6
> > ID4 pHf 3.78 5 6 6
> > ID5 pHf 3.92 6 6 6
> > ID6 pHs 2.76 5 5 6
> > ID7 pHs 3.31 6 7 6
> > ID8 pHs 4.46 6 7 6
> > ID9 pHs 2.19 5 5 6
> > ID10 pHs 3.83 5 5 6
> > ID1 pHf 5.07 7 7 9
> > ID2 pHf 6.42 7 8 9
> > ID3 pHf 5.43 6 8 9
> > ID4 pHf 6.83 6 8 9
> > ID5 pHf 6.26 6 8 9
> > ID6 pHs 4.57 6 9 9
> > ID7 pHs 5.05 6 7 9
> > ID8 pHs 6.27 6 8 9
> > ID9 pHs 3.37 5 7 9
> > ID10 pHs 5.38 6 8 9
> > ID1 pHf 5.58 7 9 12
> > ID2 pHf 7.43 8 9 12
> > ID3 pHf 6.18 8 10 12
> > ID4 pHf 6.91 7 10 12
> > ID5 pHf 6.78 7 10 12
> > ID6 pHs 4.99 6 13 12
> > ID7 pHs 5.50 7 8 12
> > ID8 pHs 6.56 7 10 12
> > ID9 pHs 3.72 6 10 12
> > ID10 pHs 5.94 6 10 12
> >
> >
> > I used the procedure described in Crawley´s new R
> > Book.
> > For two of the tree response variables
> > (heightfra,leaffra) it doesn´t work, while it
> worked
> > with leafvim (but in another R session, yesterday,
> > leaffra worked as well...).
> >
> > Here the commands:
> >
> > > attach(test)
> > > names(test)
> > [1] "week" "rep" "treat"
> "heightfra"
> > "leaffra" "leafvim"
> > > library(nlme)
> > >
> >
>
test<-groupedData(heightfra~week|rep,outer=~treat,test)
> > > model1<-lme(heightfra~treat,random=~week|rep)
> > Error in lme.formula(heightfra ~ treat, random =
> ~week
> > | rep) :
> > nlminb problem, convergence error code =
> 1;
> > message = iteration limit reached without
> convergence
> > (9)
> >
> > >
> >
>
test<-groupedData(leaffra~week|rep,outer=~treat,test)
> > > model2<-lme(leaffra~treat,random=~week|rep)
>
=== message truncated ===
More information about the R-help
mailing list