Menelaos Stavrinides
2011-Jul-25 12:39 UTC
[R] Wide confidence intervals or Error message in a mixed effects model (nlme)
I am analyzing a dataset on the effects of six pesticides on population growth rate of a predatory mite. The response variable is the population growth rate of the mite (ranges from negative to positive) and the exploratory variable is a categorical variable (treatment). The experiment was blocked in time (3 blocks / replicates per block) and it is unbalanced - at least 1 replicate per block. I am analyzing the data in nlme using model<-lme(growth.rate~treatment,random=~1|block). When I ran intervals(model), the confidence intervals of the variance of the random factor range from 0 to inf. Any comments as to why I get unrealistic confidence intervals for the random factor? In another study, I am investigating the interactions between pesticides in a two-way design: (pesticideA x no pesticide A) crossed with (pesticideB x no pesticide B). The blocking is as above, and the data are unbalanced again. The model is defined as model<-lme(growth.rate~pestA*pestB,random=~1|block). When I run intervals (model), I usually get the following error message: "Error in intervals.lme(model) : Cannot get confidence intervals on var-cov components: Non-positive definite approximate variance-covariance". I have read on the mailing list that the error message appears when the model is not well-specified, but I do not see an alternative way of specifying the model. Any ideas as to why I get wide confidence intervals or the error message? Any recommendations on other possibilities for analyzing the data are greatly appreciated. Thank you very much, Menelaos Menelaos Stavrinides Lecturer, Dept. of Agricultural Sciences Cyprus University of Technology P.O. Box 50329, 3603 Limassol, Cyprus Tel.: + 357 25002186 Fax: + 357 25002767 Email: m.stavrinides@cut.ac.cy <html> <body> <img src="http://www.cut.ac.cy/images/environmentalSign.gif"/> </body> </html> [[alternative HTML version deleted]]
Dieter Menne
2011-Jul-25 14:42 UTC
[R] Wide confidence intervals or Error message in a mixed effects model (nlme)
Menelaos Stavrinides-2 wrote:> > I am analyzing a dataset on the effects of six pesticides on population > growth rate of a predatory mite. The response variable is the population > growth rate of the mite (ranges from negative to positive) and the > exploratory variable is a categorical variable (treatment). The > experiment was blocked in time (3 blocks / replicates per block) and it > is unbalanced - at least 1 replicate per block. I am analyzing the data > in nlme using model<-lme(growth.rate~treatment,random=~1|block). > > > In another study, I am investigating the interactions between pesticides > in a two-way design: (pesticideA x no pesticide A) crossed with > (pesticideB x no pesticide B). The blocking is as above, and the data > are unbalanced again. The model is defined as > model<-lme(growth.rate~pestA*pestB,random=~1|block). When I run > intervals (model), I usually get the following error message: "Error in > intervals.lme(model) : Cannot get confidence intervals on var-cov > components: Non-positive definite approximate variance-covariance". > >It depends on your data which you did not show. At least show str() or summary(), it could also be a forgotten factor(). With simulated data, it is quite easy to get reasonable-looking cases where var-cov is degenerate. Dieter library(nlme) set.seed(47) d = expand.grid(block=LETTERS[1:3],treatment=letters[1:6]) d$growth.rate = as.integer(d$treatment)*0.2+rnorm(nrow(d)) #d.lme = lme(growth.rate~treatment, random=~1|block,data=d) for (i in 1:100){ d1= d[sample(nrow(d) ,nrow(d)-3) ,] d.lme = lme(growth.rate~treatment, random=~1|block,data=d1) iv = try(intervals(d.lme), TRUE) if (inherits(iv, "try-error")){ print(iv) print(table(d1[,c("block","treatment")])) } } -- View this message in context: http://r.789695.n4.nabble.com/Wide-confidence-intervals-or-Error-message-in-a-mixed-effects-model-nlme-tp3692637p3692986.html Sent from the R help mailing list archive at Nabble.com.
Bert Gunter
2011-Jul-25 15:29 UTC
[R] Wide confidence intervals or Error message in a mixed effects model (nlme)
Inline below. -- Bert On Mon, Jul 25, 2011 at 7:42 AM, Dieter Menne <dieter.menne at menne-biomed.de> wrote:> > Menelaos Stavrinides-2 wrote: >> >> I am analyzing a dataset on the effects of six pesticides on population >> growth rate of a predatory mite. The response variable is the population >> growth rate of the mite (ranges from negative to positive) and the >> exploratory variable is a categorical variable (treatment). The >> experiment was blocked in time (3 blocks / replicates per block) and it >> is unbalanced - at least 1 replicate per block. I am analyzing the data >> in nlme using model<-lme(growth.rate~treatment,random=~1|block). >> >> >> In another study, I am investigating the interactions between pesticides >> in a two-way design: (pesticideA x no pesticide A) crossed with >> (pesticideB x no pesticide B). The blocking is as above, and the data >> are unbalanced again. The model is defined as >> model<-lme(growth.rate~pestA*pestB,random=~1|block). When I run >> intervals (model), I usually get the following error message: "Error in >> intervals.lme(model) : Cannot get confidence intervals on var-cov >> components: Non-positive definite approximate variance-covariance". >> >> > > It depends on your data which you did not show. At least show str() or > summary(), it could also be a forgotten factor(). With simulated data, it is > quite easy to get reasonable-looking cases where var-cov is degenerate.Yes, especially when you must estimate the the b(b+1)/2 = 6 covariance parameters for b=3 blocks of the unstructured full variance covariance matrix. You might try using a structured matrix (e.g. ?pdDiag ) to reduce the number of random effects parameters you are estimating. Remember, this is a nonlinear estimation problem, and you need a "lot" of data to get good estimates in nonlinear estimation. I would also address future questions of this sort to the R-sig-mixed-models list, as this is not really an r-help kind of question. -- Bert> > Dieter > > library(nlme) > set.seed(47) > d = expand.grid(block=LETTERS[1:3],treatment=letters[1:6]) > d$growth.rate = as.integer(d$treatment)*0.2+rnorm(nrow(d)) > #d.lme = lme(growth.rate~treatment, random=~1|block,data=d) > > for (i in 1:100){ > ?d1= d[sample(nrow(d) ,nrow(d)-3) ?,] > ?d.lme = lme(growth.rate~treatment, random=~1|block,data=d1) > ?iv = try(intervals(d.lme), TRUE) > ?if (inherits(iv, "try-error")){ > ? ?print(iv) > ? ?print(table(d1[,c("block","treatment")])) > ?} > } > > > > > > > -- > View this message in context: http://r.789695.n4.nabble.com/Wide-confidence-intervals-or-Error-message-in-a-mixed-effects-model-nlme-tp3692637p3692986.html > Sent from the R help mailing list archive at Nabble.com. > > ______________________________________________ > R-help at 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. >-- "Men by nature long to get on to the ultimate truths, and will often be impatient with elementary studies or fight shy of them. If it were possible to reach the ultimate truths without the elementary studies usually prefixed to them, these would not be preparatory studies but superfluous diversions." -- Maimonides (1135-1204) Bert Gunter Genentech Nonclinical Biostatistics