Re: [R] glm with nesting
"Jeffrey Stratford" <[EMAIL PROTECTED]> writes: > I just had a manuscript returned with the biggest problem being the > analysis. Instead of using principal components in a regression I've > been asked to analyze a few variables separately. So that's what I'm > doing. > I pulled a feather from young birds and we quantified certain aspects of > the color of those feathers. > Since I often have more than one sample > from a nest, I thought I should use a nested design. Notwithstanding comments below, that quote could be aiming for the fortunes package... > Here's the code I've been using and I'd appreciate if someone could look > it over and see if it was correct. > > bb.glm1 <- glm(rtot ~ box/(julian +purbank), data=bbmale, > family="gaussian", na.action=na.omit) > > where rtot = total reflectance, box = nest box (i.e., birdhouse), > julian = day of the year and purbank = the proportion of urban cover in > a 1 km buffer around the nest box. I'm not interested in the box effect > and I've seperated males and female chicks. > > I've asked about nestedness before and I was given code that included > "|" to indicate nestedness but this indicates a grouping does it not? I > suspect that there is something wrong. In the summary I get > > Coefficients: > Estimate Std. Error t value Pr(>|t|) > (Intercept) 2.880e-01 3.224e-03 89.322 <2e-16 *** > box -3.219e-05 6.792e-05 -0.4740.636 > box:julian 7.093e-08 3.971e-07 0.1790.859 > box:purbank -1.735e-05 1.502e-04 -0.1150.908 Several things look wrong here. Most importantly, you appear to have single-degree of freedom effects (t tests) of things that appear not to be linear effects: Certainly, you have more than two nest boxes, but also day of year as a linear term looks suspicious to me. Unless there is something I have missed completely, "box" should be a factor variable, and you might also need trigonometric terms for the julian effect (depending on what sort of time spans we are talking about.) Secondly, notation like box/julian suggests that julian only makes sense within a nest box i.e. 1st of March in one box is completely different from 1st of March in another box (the notation is more commonly used to describe bird number within nests and the like). And with purbank presumably constant for measurements from the same box, the box:purbank term looks strange indeed. If you want to take account of a between-box variation in the effect of covariates, you probably need to add them as variance components, but this requires non-glm software, either lme() or lmer(). However, instructing you on those is outside the scope of this mailing list, and you may need to find a local consultant. > The other question I have is how do I test a null hypothesis - no > explanatory variables? [rtot ~ NULL?] > > 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 > > __ > 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 > and provide commented, minimal, self-contained, reproducible code. > -- O__ Peter Dalgaard Øster Farimagsgade 5, Entr.B c/ /'_ --- Dept. of Biostatistics PO Box 2099, 1014 Cph. K (*) \(*) -- University of Copenhagen Denmark Ph: (+45) 35327918 ~~ - ([EMAIL PROTECTED]) FAX: (+45) 35327907 __ 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 and provide commented, minimal, self-contained, reproducible code.
Re: [R] glm with nesting
Peter and list, Thanks for the response. A did add box as a factor (box <- factor(box)). Julian should be linear - bluebird chicks are bluer as the season progresses from March to August. I did try the following rtot.lme <- lmer(rtot ~ sex +(purban|box:chick) + (purban|box), data=bb, na.action=na.omit) # from H. Doran and rtot.lme2 <- lme(fixed=rtot ~ sex + purban + sexv:purban, data = bb, random = ~1 |box) # from K. Jones <[EMAIL PROTECTED] but these did not work (months ago and I don't remember exactly why) and I have since seperated males and females and added day of the year (julian). But "|" does indicate grouping not nested, correct? Could someone suggest some coding that might work? Thanks again, Jeff >>> Peter Dalgaard <[EMAIL PROTECTED]> 10/05/06 7:14 AM >>> "Jeffrey Stratford" <[EMAIL PROTECTED]> writes: > I just had a manuscript returned with the biggest problem being the > analysis. Instead of using principal components in a regression I've > been asked to analyze a few variables separately. So that's what I'm > doing. > I pulled a feather from young birds and we quantified certain aspects of > the color of those feathers. > Since I often have more than one sample > from a nest, I thought I should use a nested design. Notwithstanding comments below, that quote could be aiming for the fortunes package... > Here's the code I've been using and I'd appreciate if someone could look > it over and see if it was correct. > > bb.glm1 <- glm(rtot ~ box/(julian +purbank), data=bbmale, > family="gaussian", na.action=na.omit) > > where rtot = total reflectance, box = nest box (i.e., birdhouse), > julian = day of the year and purbank = the proportion of urban cover in > a 1 km buffer around the nest box. I'm not interested in the box effect > and I've seperated males and female chicks. > > I've asked about nestedness before and I was given code that included > "|" to indicate nestedness but this indicates a grouping does it not? I > suspect that there is something wrong. In the summary I get > > Coefficients: > Estimate Std. Error t value Pr(>|t|) > (Intercept) 2.880e-01 3.224e-03 89.322 <2e-16 *** > box -3.219e-05 6.792e-05 -0.4740.636 > box:julian 7.093e-08 3.971e-07 0.1790.859 > box:purbank -1.735e-05 1.502e-04 -0.1150.908 Several things look wrong here. Most importantly, you appear to have single-degree of freedom effects (t tests) of things that appear not to be linear effects: Certainly, you have more than two nest boxes, but also day of year as a linear term looks suspicious to me. Unless there is something I have missed completely, "box" should be a factor variable, and you might also need trigonometric terms for the julian effect (depending on what sort of time spans we are talking about.) Secondly, notation like box/julian suggests that julian only makes sense within a nest box i.e. 1st of March in one box is completely different from 1st of March in another box (the notation is more commonly used to describe bird number within nests and the like). And with purbank presumably constant for measurements from the same box, the box:purbank term looks strange indeed. If you want to take account of a between-box variation in the effect of covariates, you probably need to add them as variance components, but this requires non-glm software, either lme() or lmer(). However, instructing you on those is outside the scope of this mailing list, and you may need to find a local consultant. > The other question I have is how do I test a null hypothesis - no > explanatory variables? [rtot ~ NULL?] > > 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 > 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. > -- O__ Peter Dalgaard Øster Farimagsgade 5, Entr.B c/ /'_ --- Dept. of Biostatistics PO Box 2099, 1014 Cph. K (*) \(*) -- University of Copenhagen Denmark Ph: (+45) 35327918 ~~ - ([EMAIL PROTECTED]) FAX: (+45) 35327907 __ 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 and provide commented, minimal, self-contained, reproducible code.
Re: [R] glm with nesting
It's not really possible to help without knowing what errors you received and maybe some reproducible code. I think I remember this, though. From what I recall, there was no distinction between box and chick, so you cannot estimate both variance components. > -Original Message- > From: [EMAIL PROTECTED] > [mailto:[EMAIL PROTECTED] On Behalf Of > Jeffrey Stratford > Sent: Thursday, October 05, 2006 9:27 AM > To: [EMAIL PROTECTED] > Cc: r-help@stat.math.ethz.ch > Subject: Re: [R] glm with nesting > > Peter and list, > > Thanks for the response. A did add box as a factor (box <- > factor(box)). Julian should be linear - bluebird chicks are > bluer as the season progresses from March to August. > > I did try the following > > rtot.lme <- lmer(rtot ~ sex +(purban|box:chick) + > (purban|box), data=bb, > na.action=na.omit) # from H. Doran > > and > > rtot.lme2 <- lme(fixed=rtot ~ sex + purban + sexv:purban, > data = bb, random = ~1 |box) # from K. Jones <[EMAIL PROTECTED] > > but these did not work (months ago and I don't remember > exactly why) and I have since seperated males and females and > added day of the year (julian). But "|" does indicate > grouping not nested, correct? > > Could someone suggest some coding that might work? > > Thanks again, > > Jeff > > > >>> Peter Dalgaard <[EMAIL PROTECTED]> 10/05/06 7:14 AM >>> > "Jeffrey Stratford" <[EMAIL PROTECTED]> writes: > > > I just had a manuscript returned with the biggest problem being the > > analysis. Instead of using principal components in a > regression I've > > been asked to analyze a few variables separately. So that's > what I'm > > doing. > > > I pulled a feather from young birds and we quantified > certain aspects > > of the color of those feathers. > > > Since I often have more than one sample from a nest, I thought I > > should use a nested design. > > Notwithstanding comments below, that quote could be aiming > for the fortunes package... > > > > Here's the code I've been using and I'd appreciate if someone could > > look it over and see if it was correct. > > > > bb.glm1 <- glm(rtot ~ box/(julian +purbank), data=bbmale, > > family="gaussian", na.action=na.omit) > > > > where rtot = total reflectance, box = nest box (i.e., birdhouse), > > julian = day of the year and purbank = the proportion of > urban cover > > in a 1 km buffer around the nest box. I'm not interested > in the box effect > > and I've seperated males and female chicks. > > > > I've asked about nestedness before and I was given code > that included > > "|" to indicate nestedness but this indicates a grouping > does it not? > > I suspect that there is something wrong. In the summary I get > > > > Coefficients: > > Estimate Std. Error t value Pr(>|t|) > > (Intercept) 2.880e-01 3.224e-03 89.322 <2e-16 *** > > box -3.219e-05 6.792e-05 -0.4740.636 > > box:julian 7.093e-08 3.971e-07 0.1790.859 > > box:purbank -1.735e-05 1.502e-04 -0.1150.908 > > Several things look wrong here. > > Most importantly, you appear to have single-degree of freedom > effects (t tests) of things that appear not to be linear > effects: Certainly, you have more than two nest boxes, but > also day of year as a linear term looks suspicious to me. > Unless there is something I have missed completely, "box" > should be a factor variable, and you might also need > trigonometric terms for the julian effect (depending on what > sort of time spans we are talking about.) > > Secondly, notation like box/julian suggests that julian only > makes sense within a nest box i.e. 1st of March in one box is > completely different from 1st of March in another box (the > notation is more commonly used to describe bird number within > nests and the like). And with purbank presumably constant for > measurements from the same box, the box:purbank term looks > strange indeed. > > If you want to take account of a between-box variation in the > effect of covariates, you probably need to add them as > variance components, but this requires non-glm software, > either lme() or lmer(). However, instructing you on those is > outside the scope of this mailing list, and you may need to > find a local consultant. > > > The other question I have is how do I test a null hypothesis - no > > explanatory v
Re: [R] glm with nesting
Harold and list,
I've changed a few things since the last time so I'm really starting
from scratch.
I start with
bbmale <- read.csv("c:\\eabl\\2004\\feathers\\male_feathers2.csv",
header=TRUE)
box <-factor(box)
chick <- factor(chick)
Here's a sample of the data
box,chick,julian,cltchsz,mrtot,cuv,cblue,purbank,purban2,purban1,pgrassk,pgrass2,pgrass1,grassdist,grasspatchk
1,2,141,2,21.72290152,0.305723811,0.327178868,0.003813435,0.02684564,0.06896552,0.3282487,0.6845638,0.7586207,0,3.73
4,1,164,4,18.87699007,0.281863299,0.310935559,0.06072162,0.2080537,0.06896552,0.01936052,0,0,323.1099,0.2284615
4,2,164,4,19.64359348,0.294117388,0.316049817,0.06072162,0.2080537,0.06896552,0.01936052,0,0,323.1099,0.2284615
7,1,118,4,13.48699876,0.303649408,0.31765218,0.3807568,0.4362416,0.6896552,0.06864183,0.03355705,0,94.86833,0.468
12,1,180,4,21.42196378,0.289731361,0.317562294,0.09238011,0.1342282,0,0.2430127,0.8322148,1,0,1.199032
12,2,180,4,18.79487905,0.286052077,0.316367349,0.09238011,0.1342282,0,0.2430127,0.8322148,1,0,1.199032
12,3,180,4,12.83127682,0.260197475,0.292636914,0.09238011,0.1342282,0,0.2430127,0.8322148,1,0,1.199032
15,1,138,4,20.07161467,0.287632782,0.318671887,0.07046477,0.03355705,0.03448276,0.2755622,0.6577181,0.8275862,0,1.503818
15,2,138,4,17.61146256,0.305581768,0.315848051,0.07046477,0.03355705,0.03448276,0.2755622,0.6577181,0.8275862,0,1.503818
15,3,138,4,20.36397134,0.271795667,0.30539683,0.07046477,0.03355705,0.03448276,0.2755622,0.6577181,0.8275862,0,1.503818
15,4,138,4,20.81940158,0.269468041,0.304160648,0.07046477,0.03355705,0.03448276,0.2755622,0.6577181,0.8275862,0,1.503818
As you can see I have multiple boxes (> 70). Sometimes I have multiple
chicks per box each having their own response to mrtot, cuv, and cblue
but the same landscape variables for that box. Chick number is randomly
assigned and is not an effect I'm interested in. I'm really not
interested in the box effect either. I would like to know if landscape
affects the color of chicks (which may be tied into chick
health/physiology). We also know that chicks get bluer as the season
progresses and that clutch size (cltchsz) has an effect so I'm including
that as covariates.
Hopefully, this clears things up a bit.
I do have the MASS and MEMS (Pineiro's) texts in hand.
Many thanks,
Jeff
__
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
and provide commented, minimal, self-contained, reproducible code.
Re: [R] glm with nesting
Jeffrey:
Please... May I repeat what Peter Dalgaard already said: consult a local
statistician. The structure of your study is sufficiently complicated that
your stat 101 training is inadequate. Get professional help, which this list
is not set up to provide (though it often does, through the good offices and
patience of many wise contributors).
Bert Gunter
Genentech Nonclinical Statistics
South San Francisco, CA 94404
> -Original Message-
> From: [EMAIL PROTECTED]
> [mailto:[EMAIL PROTECTED] On Behalf Of
> Jeffrey Stratford
> Sent: Thursday, October 05, 2006 7:46 AM
> To: [EMAIL PROTECTED]; r-help@stat.math.ethz.ch
> Subject: Re: [R] glm with nesting
>
> Harold and list,
>
> I've changed a few things since the last time so I'm really starting
> from scratch.
>
> I start with
>
> bbmale <- read.csv("c:\\eabl\\2004\\feathers\\male_feathers2.csv",
> header=TRUE)
> box <-factor(box)
> chick <- factor(chick)
>
> Here's a sample of the data
>
> box,chick,julian,cltchsz,mrtot,cuv,cblue,purbank,purban2,purba
> n1,pgrassk,pgrass2,pgrass1,grassdist,grasspatchk
> 1,2,141,2,21.72290152,0.305723811,0.327178868,0.003813435,0.02
> 684564,0.06896552,0.3282487,0.6845638,0.7586207,0,3.73
> 4,1,164,4,18.87699007,0.281863299,0.310935559,0.06072162,0.208
> 0537,0.06896552,0.01936052,0,0,323.1099,0.2284615
> 4,2,164,4,19.64359348,0.294117388,0.316049817,0.06072162,0.208
> 0537,0.06896552,0.01936052,0,0,323.1099,0.2284615
> 7,1,118,4,13.48699876,0.303649408,0.31765218,0.3807568,0.43624
> 16,0.6896552,0.06864183,0.03355705,0,94.86833,0.468
> 12,1,180,4,21.42196378,0.289731361,0.317562294,0.09238011,0.13
> 42282,0,0.2430127,0.8322148,1,0,1.199032
> 12,2,180,4,18.79487905,0.286052077,0.316367349,0.09238011,0.13
> 42282,0,0.2430127,0.8322148,1,0,1.199032
> 12,3,180,4,12.83127682,0.260197475,0.292636914,0.09238011,0.13
> 42282,0,0.2430127,0.8322148,1,0,1.199032
> 15,1,138,4,20.07161467,0.287632782,0.318671887,0.07046477,0.03
> 355705,0.03448276,0.2755622,0.6577181,0.8275862,0,1.503818
> 15,2,138,4,17.61146256,0.305581768,0.315848051,0.07046477,0.03
> 355705,0.03448276,0.2755622,0.6577181,0.8275862,0,1.503818
> 15,3,138,4,20.36397134,0.271795667,0.30539683,0.07046477,0.033
> 55705,0.03448276,0.2755622,0.6577181,0.8275862,0,1.503818
> 15,4,138,4,20.81940158,0.269468041,0.304160648,0.07046477,0.03
> 355705,0.03448276,0.2755622,0.6577181,0.8275862,0,1.503818
>
> As you can see I have multiple boxes (> 70). Sometimes I
> have multiple
> chicks per box each having their own response to mrtot, cuv,
> and cblue
> but the same landscape variables for that box. Chick number
> is randomly
> assigned and is not an effect I'm interested in. I'm really not
> interested in the box effect either. I would like to know if
> landscape
> affects the color of chicks (which may be tied into chick
> health/physiology). We also know that chicks get bluer as the season
> progresses and that clutch size (cltchsz) has an effect so
> I'm including
> that as covariates.
>
> Hopefully, this clears things up a bit.
>
> I do have the MASS and MEMS (Pineiro's) texts in hand.
>
> Many thanks,
>
> Jeff
>
> __
> 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
> and provide commented, minimal, self-contained, reproducible code.
>
__
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
and provide commented, minimal, self-contained, reproducible code.
