An embedded and charset-unspecified text was scrubbed... Name: not available URL: <https://stat.ethz.ch/pipermail/r-help/attachments/20090310/13bc5c32/attachment-0002.pl>
help structuring mixed model using lmer()
6 messages · Simon Pickett, Douglas Bates, Mark Difford +1 more
On Tue, Mar 10, 2009 at 10:15 AM, Simon Pickett <simon.pickett at bto.org> wrote:
This is partly a statistical question as well as a question about R, but I am stumped!
I have count data from various sites across years. (Not all of the sites in the study appear in all years). Each site has its own habitat score "habitat" that remains constant across all years.
I want to know if counts declined faster on sites with high "habitat" scores.
I can construct a model that tests for the effect of habitat as a main effect, controlling for year
model1<-lmer(count~habitat+yr+(1|site), family=quasibinomial,data=m) model2<-lmer(count~yr+(1|site), family=quasibinomial,data=m) anova(model1,model2)
I'm curious as to why you use the quasibinomial family for count data. When you say "count data" do you mean just presence/absence or an actual count of the number present? Generally the binomial and quasibinomial families are used when you have a binary response, and the poisson or quasipoisson family are used for responses that are counts.
Cheers, Actually I was using quasipoisson for my models, but for the puposes of my example, it doesnt really matter. I am trying to work out a way of quantifying whether the slopes (for years) are covary with habitat scores. The more I think about it, the more I am convinced that it isnt possible do to that using a glm approach. I think I have to run separate models for each site, calculate the gradient, then do a lm with gradient explained by habitat score.... Thanks, Simon.
On Tue, Mar 10, 2009 at 10:15 AM, Simon Pickett <simon.pickett at bto.org> wrote:
This is partly a statistical question as well as a question about R, but I am stumped!
I have count data from various sites across years. (Not all of the sites in the study appear in all years). Each site has its own habitat score "habitat" that remains constant across all years.
I want to know if counts declined faster on sites with high "habitat" scores.
I can construct a model that tests for the effect of habitat as a main effect, controlling for year
model1<-lmer(count~habitat+yr+(1|site), family=quasibinomial,data=m) model2<-lmer(count~yr+(1|site), family=quasibinomial,data=m) anova(model1,model2)
I'm curious as to why you use the quasibinomial family for count data. When you say "count data" do you mean just presence/absence or an actual count of the number present? Generally the binomial and quasibinomial families are used when you have a binary response, and the poisson or quasipoisson family are used for responses that are counts.
Hi Simon, Have a look at Chap. 11 of "An Introduction to R" (one of R's manuals), which explains the different ways of specifying models using formulae. Briefly, y ~ x1 * x2 expands to y ~ x1 + x2 + x1:x2, where the last term (interaction term) amounts to a test of slope. Normally you would read its significance from F/chisq/p-value. Many practitioners consider the L.Ratio test to be a better option. For the fixed effects part in lmer() do: mod1 <- y ~ x1 + x2 == y ~ x1 + x2 mod2 <- y ~ x1 * x2 == y ~ x1 + x2 + x1:x2 anova(mod1, mod2) This will tell you if you need to worry about interaction or whether slopes are parallel. Regards, Mark.
Simon Pickett-4 wrote:
Cheers, Actually I was using quasipoisson for my models, but for the puposes of my example, it doesnt really matter. I am trying to work out a way of quantifying whether the slopes (for years) are covary with habitat scores. The more I think about it, the more I am convinced that it isnt possible do to that using a glm approach. I think I have to run separate models for each site, calculate the gradient, then do a lm with gradient explained by habitat score.... Thanks, Simon.
On Tue, Mar 10, 2009 at 10:15 AM, Simon Pickett <simon.pickett at bto.org> wrote:
This is partly a statistical question as well as a question about R, but I am stumped!
I have count data from various sites across years. (Not all of the sites in the study appear in all years). Each site has its own habitat score "habitat" that remains constant across all years.
I want to know if counts declined faster on sites with high "habitat" scores.
I can construct a model that tests for the effect of habitat as a main effect, controlling for year
model1<-lmer(count~habitat+yr+(1|site), family=quasibinomial,data=m) model2<-lmer(count~yr+(1|site), family=quasibinomial,data=m) anova(model1,model2)
I'm curious as to why you use the quasibinomial family for count data. When you say "count data" do you mean just presence/absence or an actual count of the number present? Generally the binomial and quasibinomial families are used when you have a binary response, and the poisson or quasipoisson family are used for responses that are counts.
______________________________________________ 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.
View this message in context: http://www.nabble.com/help-structuring-mixed-model-using-lmer%28%29-tp22436596p22441985.html Sent from the R help mailing list archive at Nabble.com.
Hi Simon, Have a look at Chap. 11 of "An Introduction to R" (one of R's manuals), which explains the different ways of specifying models using formulae. Briefly, y ~ x1 * x2 expands to y ~ x1 + x2 + x1:x2, where the last term (interaction term) amounts to a test of slope. Normally you would read its significance from F/chisq/p-value. Many practitioners consider the L.Ratio test to be a better option. For the fixed effects part in lmer() do: mod1 <- y ~ x1 + x2 == y ~ x1 + x2 mod2 <- y ~ x1 * x2 == y ~ x1 + x2 + x1:x2 anova(mod1, mod2) This will tell you if you need to worry about interaction or whether slopes are parallel. Regards, Mark.
Simon Pickett-4 wrote:
Cheers, Actually I was using quasipoisson for my models, but for the puposes of my example, it doesnt really matter. I am trying to work out a way of quantifying whether the slopes (for years) are covary with habitat scores. The more I think about it, the more I am convinced that it isnt possible do to that using a glm approach. I think I have to run separate models for each site, calculate the gradient, then do a lm with gradient explained by habitat score.... Thanks, Simon.
On Tue, Mar 10, 2009 at 10:15 AM, Simon Pickett <simon.pickett at bto.org> wrote:
This is partly a statistical question as well as a question about R, but I am stumped!
I have count data from various sites across years. (Not all of the sites in the study appear in all years). Each site has its own habitat score "habitat" that remains constant across all years.
I want to know if counts declined faster on sites with high "habitat" scores.
I can construct a model that tests for the effect of habitat as a main effect, controlling for year
model1<-lmer(count~habitat+yr+(1|site), family=quasibinomial,data=m) model2<-lmer(count~yr+(1|site), family=quasibinomial,data=m) anova(model1,model2)
I'm curious as to why you use the quasibinomial family for count data. When you say "count data" do you mean just presence/absence or an actual count of the number present? Generally the binomial and quasibinomial families are used when you have a binary response, and the poisson or quasipoisson family are used for responses that are counts.
______________________________________________ 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.
View this message in context: http://www.nabble.com/help-structuring-mixed-model-using-lmer%28%29-tp22436596p22442012.html Sent from the R help mailing list archive at Nabble.com.
Mark Difford <mark_difford <at> yahoo.co.uk> writes:
Briefly, y ~ x1 * x2 expands to y ~ x1 + x2 + x1:x2, where the last term (interaction term) amounts to a test of slope. Normally you would read its significance from F/chisq/p-value. Many practitioners consider the L.Ratio test to be a better option. For the fixed effects part in lmer() do: mod1 <- y ~ x1 + x2 == y ~ x1 + x2 mod2 <- y ~ x1 * x2 == y ~ x1 + x2 + x1:x2 anova(mod1, mod2) This will tell you if you need to worry about interaction or whether slopes are parallel.
... except that you'd better hope that you have a large number
of random units (years? I forget now), or the LR test will be
unreliable -- see Pinheiro and Bates 2000, and refer further questions
to the r-sig-mixed-models list ...
cheers
Ben Bolker