Dear R-users,
I did some more research and I'm still not sure how to set up an ANCOVA
with nestedness. Specifically I'm not sure how to express chicks nested
within boxes. I will be getting Pinheiro & Bates (Mixed Effects Models
in S and S-Plus) but it will not arrive for another two weeks from our
interlibrary loan.
The goal is to determine if there are urbanization (purban) effects on
chick health (rtot) and if there are differences between sexes (sex) and
the effect of being in the same clutch (box).
The model is rtot = sex + purban + (chick)box.
I've loaded the package lme4. And the code I have so far is
bb <- read.csv("C:\\eabl\\eabl_feather04.csv", header=TRUE)
bb$sex <- factor(bb$sex)
rtot.lme <- lme(bb$rtot~bb$sex, bb$purban|bb$chick/bb$box,
na.action=na.omit)
but this is not working.
Any suggestions would be greatly appreciated.
Thanks,
Jeff
****************************************
Jeffrey A. Stratford, Ph.D.
Postdoctoral Associate
331 Funchess Hall
Department of Biological Sciences
Auburn University
Auburn, AL 36849
334-329-9198
FAX 334-844-9234
http://www.auburn.edu/~stratja
Dear Jeff:
I see the issues in your code and have provided what I think will solve
your problem. It is often much easier to get help on this list when you
provide a small bit of data that can be replicated and you state what
the error messages are that you are receiving. OK, with that said, here
is what I see. First, you do not need to use the syntax bb$sex in your
model, this can be significantly simplified. Second, you do not have a
random statement in your model.
Here is your original model:
lme(bb$rtot~bb$sex, bb$purban|bb$chick/bb$box, na.action=na.omit)
Here is what it should be:
lme(fixed = rtot~sex, random=~purban|chick/box, na.action=na.omit,
data=bb)
Notice there is a fixed and random call. You can simplify this as
lme(rtot~sex, random=~purban|chick/box, na.action=na.omit, bb)
Note, you can eliminate the "fixed=" portion but not the random
statement.
Last, if you want to do this in lmer, the newer function for mixed
models in the Matrix package, you would do
lmer(rtot~sex + (purban|box:chick) + (purban|box), na.action=na.omit,
data=bb)
Hope this helps.
Harold
-----Original Message-----
From: r-help-bounces at stat.math.ethz.ch
[mailto:r-help-bounces at stat.math.ethz.ch] On Behalf Of Jeffrey Stratford
Sent: Tuesday, January 24, 2006 11:34 AM
To: r-help at stat.math.ethz.ch
Subject: [R] nested ANCOVA: still confused
Dear R-users,
I did some more research and I'm still not sure how to set up an ANCOVA
with nestedness. Specifically I'm not sure how to express chicks nested
within boxes. I will be getting Pinheiro & Bates (Mixed Effects Models
in S and S-Plus) but it will not arrive for another two weeks from our
interlibrary loan.
The goal is to determine if there are urbanization (purban) effects on
chick health (rtot) and if there are differences between sexes (sex) and
the effect of being in the same clutch (box).
The model is rtot = sex + purban + (chick)box.
I've loaded the package lme4. And the code I have so far is
bb <- read.csv("C:\\eabl\\eabl_feather04.csv", header=TRUE) bb$sex
<-
factor(bb$sex) rtot.lme <- lme(bb$rtot~bb$sex,
bb$purban|bb$chick/bb$box,
na.action=na.omit)
but this is not working.
Any suggestions would be greatly appreciated.
Thanks,
Jeff
****************************************
Jeffrey A. Stratford, Ph.D.
Postdoctoral Associate
331 Funchess Hall
Department of Biological Sciences
Auburn University
Auburn, AL 36849
334-329-9198
FAX 334-844-9234
http://www.auburn.edu/~stratja
______________________________________________
R-help at stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide!
http://www.R-project.org/posting-guide.html
R-users and Harold.
First, thanks for the advice; I'm almost there.
The code I'm using now is
library(nlme)
bb <- read.csv("E:\\eabl_feather04.csv", header=TRUE)
bb$sexv <- factor(bb$sexv)
rtot.lme <- lme(fixed=rtot~sexv, random=~purban2|chick/box,
na.action=na.omit, data=bb)
A sample of the data looks like this
box chick rtot purban2 sexv
1 1 6333.51 0.026846 f
1 2 8710.884 0.026846 m
2 1 5810.007 0.161074 f
2 2 5524.33 0.161074 f
2 3 4824.474 0.161074 f
2 4 5617.641 0.161074 f
2 5 6761.724 0.161074 f
4 1 7569.673 0.208054 m
4 2 7877.081 0.208054 m
4 4 7455.55 0.208054 f
7 1 5408.287 0.436242 m
10 1 6991.727 0.14094 f
12 1 8590.207 0.134228 f
12 2 7536.747 0.134228 m
12 3 5145.342 0.134228 m
12 4 6853.628 0.134228 f
15 1 8048.717 0.033557 m
15 2 7062.196 0.033557 m
15 3 8165.953 0.033557 m
15 4 8348.58 0.033557 m
16 2 6534.775 0.751678 m
16 3 7468.827 0.751678 m
16 4 5907.338 0.751678 f
21 1 7761.983 0.221477 m
21 2 6634.115 0.221477 m
21 3 6982.923 0.221477 m
21 4 7464.075 0.221477 m
22 1 6756.733 0.281879 f
23 2 8231.496 0.134228 m
The error I'm getting is
Error in logLik.lmeStructInt(lmeSt, lmePars) :
Calloc could not allocate (590465568 of 8) memory
In addition: Warning messages:
1: Fewer observations than random effects in all level 2 groups in:
lme.formula(fixed = rtot ~ sexv, random = ~purban2 | chick/box,
2: Reached total allocation of 382Mb: see help(memory.size)
There's nothing "special" about chick 1, 2, etc. These were
simply the
order of the birds measured in each box so chick 1 in box 1 has nothing
to do with chick 1 in box 2.
Many thanks,
Jeff
****************************************
Jeffrey A. Stratford, Ph.D.
Postdoctoral Associate
331 Funchess Hall
Department of Biological Sciences
Auburn University
Auburn, AL 36849
334-329-9198
FAX 334-844-9234
http://www.auburn.edu/~stratja
****************************************>>> "Doran, Harold" <HDoran at air.org> 01/24/06 2:04
PM >>>
Dear Jeff:
I see the issues in your code and have provided what I think will solve
your problem. It is often much easier to get help on this list when you
provide a small bit of data that can be replicated and you state what
the error messages are that you are receiving. OK, with that said, here
is what I see. First, you do not need to use the syntax bb$sex in your
model, this can be significantly simplified. Second, you do not have a
random statement in your model.
Here is your original model:
lme(bb$rtot~bb$sex, bb$purban|bb$chick/bb$box, na.action=na.omit)
Here is what it should be:
lme(fixed = rtot~sex, random=~purban|chick/box, na.action=na.omit,
data=bb)
Notice there is a fixed and random call. You can simplify this as
lme(rtot~sex, random=~purban|chick/box, na.action=na.omit, bb)
Note, you can eliminate the "fixed=" portion but not the random
statement.
Last, if you want to do this in lmer, the newer function for mixed
models in the Matrix package, you would do
lmer(rtot~sex + (purban|box:chick) + (purban|box), na.action=na.omit,
data=bb)
Hope this helps.
Harold
-----Original Message-----
From: r-help-bounces at stat.math.ethz.ch
[mailto:r-help-bounces at stat.math.ethz.ch] On Behalf Of Jeffrey Stratford
Sent: Tuesday, January 24, 2006 11:34 AM
To: r-help at stat.math.ethz.ch
Subject: [R] nested ANCOVA: still confused
Dear R-users,
I did some more research and I'm still not sure how to set up an ANCOVA
with nestedness. Specifically I'm not sure how to express chicks nested
within boxes. I will be getting Pinheiro & Bates (Mixed Effects Models
in S and S-Plus) but it will not arrive for another two weeks from our
interlibrary loan.
The goal is to determine if there are urbanization (purban) effects on
chick health (rtot) and if there are differences between sexes (sex) and
the effect of being in the same clutch (box).
The model is rtot = sex + purban + (chick)box.
I've loaded the package lme4. And the code I have so far is
bb <- read.csv("C:\\eabl\\eabl_feather04.csv", header=TRUE) bb$sex
<-
factor(bb$sex) rtot.lme <- lme(bb$rtot~bb$sex,
bb$purban|bb$chick/bb$box,
na.action=na.omit)
but this is not working.
Any suggestions would be greatly appreciated.
Thanks,
Jeff
****************************************
Jeffrey A. Stratford, Ph.D.
Postdoctoral Associate
331 Funchess Hall
Department of Biological Sciences
Auburn University
Auburn, AL 36849
334-329-9198
FAX 334-844-9234
http://www.auburn.edu/~stratja
______________________________________________
R-help at stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide!
http://www.R-project.org/posting-guide.html
OK, we're getting somewhere. First, it looks as though (by the error
message) that you have a big dataset. My first recommendation is to use lmer
instead of lme, you will see a significant benefit in terms of computional
speed.
For the model this would be
lmer(rtot ~ sexv +(purban|box:chick) + (purban|box), bb, na.action=na.omit)
Now, you have run out of memory. I don't know what operating system you are
using, so go and see the appropriate FAQ for increasing memory for your OS.
Second, I made a mistake in my reply. Your random statement should be
random=~purban|box/chick denoting that chicks are nested in boxes, not boxes
nested in chicks, sorry about that.
Now, why is it that each chick within box has the same value for purban? If this
is so, why are you fitting that as a random effect? It seems not to vary across
individual chicks, right? It seems there is only an effect of box and not an
effect for chicks. Why not just fit a random effect only for box such as:
rtot.lme <- lme(fixed=rtot~sexv, random=~purban2|box, na.action=na.omit,bb)
or in lmer
lmer(rtot ~ sexv + (purban|box), bb, na.action=na.omit)
Harold
-----Original Message-----
From: Jeffrey Stratford [mailto:stratja@auburn.edu]
Sent: Tue 1/24/2006 8:57 PM
To: Doran, Harold; r-help@stat.math.ethz.ch
Cc:
Subject: RE: [R] nested ANCOVA: still confused
R-users and Harold.
First, thanks for the advice; I'm almost there.
The code I'm using now is
library(nlme)
bb <- read.csv("E:\\eabl_feather04.csv", header=TRUE)
bb$sexv <- factor(bb$sexv)
rtot.lme <- lme(fixed=rtot~sexv, random=~purban2|chick/box,
na.action=na.omit, data=bb)
A sample of the data looks like this
box chick rtot purban2 sexv
1 1 6333.51 0.026846 f
1 2 8710.884 0.026846 m
2 1 5810.007 0.161074 f
2 2 5524.33 0.161074 f
2 3 4824.474 0.161074 f
2 4 5617.641 0.161074 f
2 5 6761.724 0.161074 f
4 1 7569.673 0.208054 m
4 2 7877.081 0.208054 m
4 4 7455.55 0.208054 f
7 1 5408.287 0.436242 m
10 1 6991.727 0.14094 f
12 1 8590.207 0.134228 f
12 2 7536.747 0.134228 m
12 3 5145.342 0.134228 m
12 4 6853.628 0.134228 f
15 1 8048.717 0.033557 m
15 2 7062.196 0.033557 m
15 3 8165.953 0.033557 m
15 4 8348.58 0.033557 m
16 2 6534.775 0.751678 m
16 3 7468.827 0.751678 m
16 4 5907.338 0.751678 f
21 1 7761.983 0.221477 m
21 2 6634.115 0.221477 m
21 3 6982.923 0.221477 m
21 4 7464.075 0.221477 m
22 1 6756.733 0.281879 f
23 2 8231.496 0.134228 m
The error I'm getting is
Error in logLik.lmeStructInt(lmeSt, lmePars) :
Calloc could not allocate (590465568 of 8) memory
In addition: Warning messages:
1: Fewer observations than random effects in all level 2 groups in:
lme.formula(fixed = rtot ~ sexv, random = ~purban2 | chick/box,
2: Reached total allocation of 382Mb: see help(memory.size)
There's nothing "special" about chick 1, 2, etc. These were
simply the
order of the birds measured in each box so chick 1 in box 1 has nothing
to do with chick 1 in box 2.
Many thanks,
Jeff
****************************************
Jeffrey A. Stratford, Ph.D.
Postdoctoral Associate
331 Funchess Hall
Department of Biological Sciences
Auburn University
Auburn, AL 36849
334-329-9198
FAX 334-844-9234
http://www.auburn.edu/~stratja
****************************************>>> "Doran, Harold" <HDoran@air.org> 01/24/06 2:04 PM
>>>
Dear Jeff:
I see the issues in your code and have provided what I think will solve
your problem. It is often much easier to get help on this list when you
provide a small bit of data that can be replicated and you state what
the error messages are that you are receiving. OK, with that said, here
is what I see. First, you do not need to use the syntax bb$sex in your
model, this can be significantly simplified. Second, you do not have a
random statement in your model.
Here is your original model:
lme(bb$rtot~bb$sex, bb$purban|bb$chick/bb$box, na.action=na.omit)
Here is what it should be:
lme(fixed = rtot~sex, random=~purban|chick/box, na.action=na.omit,
data=bb)
Notice there is a fixed and random call. You can simplify this as
lme(rtot~sex, random=~purban|chick/box, na.action=na.omit, bb)
Note, you can eliminate the "fixed=" portion but not the random
statement.
Last, if you want to do this in lmer, the newer function for mixed
models in the Matrix package, you would do
lmer(rtot~sex + (purban|box:chick) + (purban|box), na.action=na.omit,
data=bb)
Hope this helps.
Harold
-----Original Message-----
From: r-help-bounces@stat.math.ethz.ch
[mailto:r-help-bounces@stat.math.ethz.ch] On Behalf Of Jeffrey Stratford
Sent: Tuesday, January 24, 2006 11:34 AM
To: r-help@stat.math.ethz.ch
Subject: [R] nested ANCOVA: still confused
Dear R-users,
I did some more research and I'm still not sure how to set up an ANCOVA
with nestedness. Specifically I'm not sure how to express chicks nested
within boxes. I will be getting Pinheiro & Bates (Mixed Effects Models
in S and S-Plus) but it will not arrive for another two weeks from our
interlibrary loan.
The goal is to determine if there are urbanization (purban) effects on
chick health (rtot) and if there are differences between sexes (sex) and
the effect of being in the same clutch (box).
The model is rtot = sex + purban + (chick)box.
I've loaded the package lme4. And the code I have so far is
bb <- read.csv("C:\\eabl\\eabl_feather04.csv", header=TRUE) bb$sex
<-
factor(bb$sex) rtot.lme <- lme(bb$rtot~bb$sex,
bb$purban|bb$chick/bb$box,
na.action=na.omit)
but this is not working.
Any suggestions would be greatly appreciated.
Thanks,
Jeff
****************************************
Jeffrey A. Stratford, Ph.D.
Postdoctoral Associate
331 Funchess Hall
Department of Biological Sciences
Auburn University
Auburn, AL 36849
334-329-9198
FAX 334-844-9234
http://www.auburn.edu/~stratja
______________________________________________
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide!
http://www.R-project.org/posting-guide.html
[[alternative HTML version deleted]]
Harold, Kingsford, and R-users, I settled on using the lmer function. I think the memory issue was more a function of my poor coding than an actual memory problem. I also switched the label from "box" to "clutch" to avoid any potential confusion with other functions. This coding seems to have worked:> eabl <- lmer(rtot ~ sexv + (purban2|clutch), maxIter=1000, data=bb,na.action=na.omit) However, I have two remaining questions: (1)how concerned should I be with the warning message below and (2) is there a way to invoke output to get an estimate of the effect of purban2 (the proportion of urban cover 200 m around a box) on feather color (rtot) and if there is a difference between the sexes? I used the summary function and it doesn't tell me much (see output below). I'll read up mixed models when Pinheiro arrives but any suggestions for diagnostics? I'm going to repeat this study and expand it by doubling or tripling the number of birds. Warning message: nlminb returned message false convergence (8) in: "LMEoptimize<-"(`*tmp*`, value = list(maxIter = 200, tolerance 1.49011611938477e-08,> summary(eabl)Linear mixed-effects model fit by REML Formula: rtot ~ sexv + (purban2 | clutch) Data: bb AIC BIC logLik MLdeviance REMLdeviance 5164.284 6997.864 -2052.142 4128.792 4104.284 Random effects: Groups Name Variance Std.Dev. Corr clutch (Intercept) 502829 709.10 purban20 1341990 1158.44 -0.477 purban20.006711409 5683957 2384.10 -0.226 0.082 purban20.01342282 1772922 1331.51 -0.386 0.176 0.067 . . . . . # of obs: 235, groups: clutch, 74 Fixed effects: Estimate Std. Error t value (Intercept) 5950.01 241.59 24.628 sexvm 1509.07 145.73 10.355 Correlation of Fixed Effects: (Intr) sexvm -0.304 Thanks many time over, Jeff **************************************** Jeffrey A. Stratford, Ph.D. Postdoctoral Associate 331 Funchess Hall Department of Biological Sciences Auburn University Auburn, AL 36849 334-329-9198 FAX 334-844-9234 http://www.auburn.edu/~stratja ****************************************>>> "Doran, Harold" <HDoran at air.org> 01/25/06 6:37 AM >>>OK, we're getting somewhere. First, it looks as though (by the error message) that you have a big dataset. My first recommendation is to use lmer instead of lme, you will see a significant benefit in terms of computional speed. For the model this would be lmer(rtot ~ sexv +(purban|box:chick) + (purban|box), bb, na.action=na.omit) Now, you have run out of memory. I don't know what operating system you are using, so go and see the appropriate FAQ for increasing memory for your OS. Second, I made a mistake in my reply. Your random statement should be random=~purban|box/chick denoting that chicks are nested in boxes, not boxes nested in chicks, sorry about that. Now, why is it that each chick within box has the same value for purban? If this is so, why are you fitting that as a random effect? It seems not to vary across individual chicks, right? It seems there is only an effect of box and not an effect for chicks. Why not just fit a random effect only for box such as: rtot.lme <- lme(fixed=rtot~sexv, random=~purban2|box, na.action=na.omit,bb) or in lmer lmer(rtot ~ sexv + (purban|box), bb, na.action=na.omit) Harold -----Original Message----- From: Jeffrey Stratford [mailto:stratja at auburn.edu] Sent: Tue 1/24/2006 8:57 PM To: Doran, Harold; r-help at stat.math.ethz.ch Cc: Subject: RE: [R] nested ANCOVA: still confused R-users and Harold. First, thanks for the advice; I'm almost there. The code I'm using now is library(nlme) bb <- read.csv("E:\\eabl_feather04.csv", header=TRUE) bb$sexv <- factor(bb$sexv) rtot.lme <- lme(fixed=rtot~sexv, random=~purban2|chick/box, na.action=na.omit, data=bb) A sample of the data looks like this box chick rtot purban2 sexv 1 1 6333.51 0.026846 f 1 2 8710.884 0.026846 m 2 1 5810.007 0.161074 f 2 2 5524.33 0.161074 f 2 3 4824.474 0.161074 f 2 4 5617.641 0.161074 f 2 5 6761.724 0.161074 f 4 1 7569.673 0.208054 m 4 2 7877.081 0.208054 m 4 4 7455.55 0.208054 f 7 1 5408.287 0.436242 m 10 1 6991.727 0.14094 f 12 1 8590.207 0.134228 f 12 2 7536.747 0.134228 m 12 3 5145.342 0.134228 m 12 4 6853.628 0.134228 f 15 1 8048.717 0.033557 m 15 2 7062.196 0.033557 m 15 3 8165.953 0.033557 m 15 4 8348.58 0.033557 m 16 2 6534.775 0.751678 m 16 3 7468.827 0.751678 m 16 4 5907.338 0.751678 f 21 1 7761.983 0.221477 m 21 2 6634.115 0.221477 m 21 3 6982.923 0.221477 m 21 4 7464.075 0.221477 m 22 1 6756.733 0.281879 f 23 2 8231.496 0.134228 m The error I'm getting is Error in logLik.lmeStructInt(lmeSt, lmePars) : Calloc could not allocate (590465568 of 8) memory In addition: Warning messages: 1: Fewer observations than random effects in all level 2 groups in: lme.formula(fixed = rtot ~ sexv, random = ~purban2 | chick/box, 2: Reached total allocation of 382Mb: see help(memory.size) There's nothing "special" about chick 1, 2, etc. These were simply the order of the birds measured in each box so chick 1 in box 1 has nothing to do with chick 1 in box 2. Many thanks, Jeff **************************************** Jeffrey A. Stratford, Ph.D. Postdoctoral Associate 331 Funchess Hall Department of Biological Sciences Auburn University Auburn, AL 36849 334-329-9198 FAX 334-844-9234 http://www.auburn.edu/~stratja ****************************************>>> "Doran, Harold" <HDoran at air.org> 01/24/06 2:04 PM >>>Dear Jeff: I see the issues in your code and have provided what I think will solve your problem. It is often much easier to get help on this list when you provide a small bit of data that can be replicated and you state what the error messages are that you are receiving. OK, with that said, here is what I see. First, you do not need to use the syntax bb$sex in your model, this can be significantly simplified. Second, you do not have a random statement in your model. Here is your original model: lme(bb$rtot~bb$sex, bb$purban|bb$chick/bb$box, na.action=na.omit) Here is what it should be: lme(fixed = rtot~sex, random=~purban|chick/box, na.action=na.omit, data=bb) Notice there is a fixed and random call. You can simplify this as lme(rtot~sex, random=~purban|chick/box, na.action=na.omit, bb) Note, you can eliminate the "fixed=" portion but not the random statement. Last, if you want to do this in lmer, the newer function for mixed models in the Matrix package, you would do lmer(rtot~sex + (purban|box:chick) + (purban|box), na.action=na.omit, data=bb) Hope this helps. Harold -----Original Message----- From: r-help-bounces at stat.math.ethz.ch [mailto:r-help-bounces at stat.math.ethz.ch] On Behalf Of Jeffrey Stratford Sent: Tuesday, January 24, 2006 11:34 AM To: r-help at stat.math.ethz.ch Subject: [R] nested ANCOVA: still confused Dear R-users, I did some more research and I'm still not sure how to set up an ANCOVA with nestedness. Specifically I'm not sure how to express chicks nested within boxes. I will be getting Pinheiro & Bates (Mixed Effects Models in S and S-Plus) but it will not arrive for another two weeks from our interlibrary loan. The goal is to determine if there are urbanization (purban) effects on chick health (rtot) and if there are differences between sexes (sex) and the effect of being in the same clutch (box). The model is rtot = sex + purban + (chick)box. I've loaded the package lme4. And the code I have so far is bb <- read.csv("C:\\eabl\\eabl_feather04.csv", header=TRUE) bb$sex <- factor(bb$sex) rtot.lme <- lme(bb$rtot~bb$sex, bb$purban|bb$chick/bb$box, na.action=na.omit) but this is not working. Any suggestions would be greatly appreciated. Thanks, Jeff **************************************** Jeffrey A. Stratford, Ph.D. Postdoctoral Associate 331 Funchess Hall Department of Biological Sciences Auburn University Auburn, AL 36849 334-329-9198 FAX 334-844-9234 http://www.auburn.edu/~stratja ______________________________________________ R-help at stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html