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)
Error in lme.formula(leaffra ~ treat, random = ~week |
rep) :
nlminb problem, convergence error code = 1;
message = iteration limit reached without convergence
(9)
>
test<-groupedData(leafvim~week|rep,outer=~treat,test)
> model3<-lme(leafvim~treat,random=~week|rep)
> summary(model)
Error in summary(model) : object "model" not found
> summary(model3)
Linear mixed-effects model fit by REML
Data: NULL
AIC BIC logLik
129.6743 139.4999 -58.83717
Random effects:
Formula: ~week | rep
Structure: General positive-definite, Log-Cholesky
parametrization
StdDev Corr
(Intercept) 4.4110478 (Intr)
week 0.7057311 -0.999
Residual 0.5976143
Fixed effects: leafvim ~ treat
Value Std.Error DF t-value p-value
(Intercept) 5.924659 0.1653596 30 35.82893 0.0000
treatpHs 0.063704 0.2338538 8 0.27241 0.7922
Correlation:
(Intr)
treatpHs -0.707
Standardized Within-Group Residuals:
Min Q1 Med Q3
Max
-1.34714254 -0.53042878 -0.01769195 0.40644540
2.29301560
Number of Observations: 40
Number of Groups: 10
Is there a solution for this problem?
Thanks!!!
Ilona
--- Douglas Bates <[EMAIL PROTECTED]> schrieb:
> On Dec 13, 2007 4:15 PM, Ilona Leyer
> <[EMAIL PROTECTED]> wrote:
> > Dear All,
> > I want to analyse treatment effects with time
> series
> > data: I measured e.g. leaf number (five replicate
> > plants) in relation to two soil pH - after 2,4,6,8
> > weeks. I used mixed effects models, but some
> analyses
> > didn´t work. It seems for me as if this is a
> randomly
> > occurring problem since sometimes the same model
> works
> > sometimes not.
> >
> > An example:
> > > names(test)
> > [1] "rep" "treat" "leaf" "week"
> > > library (lattice)
> > > library (nlme)
> > >
> test<-groupedData(leaf~week|rep,outer=~treat,test)
> > > model<-lme(leaf~treat,random=~leaf|rep)
> > Error in lme.formula(leaf~ treat, random =
> ~week|rep)
>
> Really!? You gave lme a model with random = ~ leaf |
> rep (and no data
> specification) and it tried to fit a model with
> random = ~ week | rep?
> Are you sure that is an exact transcript?
>
> > :
> > nlminb problem, convergence error code =
> 1;
> > message = iteration limit reached without
> convergence
> > (9)
>
> > Has anybody an idea to solve this problem?
>
> Oh, I have lots of ideas but without a reproducible
> example I can't
> hope to decide what might be the problem.
>
> It appears that the model may be over-parameterized.
> Assuming that
> there are 4 different values of week then ~ week |
> rep requires
> fitting 10 variance-covariance parameters. That's a
> lot.
> The error code indicates that the optimizer is
> taking
>
______________________________________________
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.
