[R-sig-ME] Nested longitudinal

Ben Bolker bbolker at gmail.com
Wed Apr 25 20:13:06 CEST 2012


  a. these weren't my comments (mostly), they were Rob Kushler's

  b. in general it's best to keep these conversations on-list, in case
someone else wants to chime in and to keep the information public and
archivable for future reference.  I've cc'd back to r-sig-mixed-models.

  Ben Bolker


On 12-04-25 02:06 PM, arun wrote:
> 
> 
> Hi Ben, Thanks for the suggestions.  I did the analysis again by
> collapsing the data by time, but I made a small change.  Instead of
> collapsing the whole data, I made it collapsed in 5 minute intervals
> (1-5), (6-10), (11-15).  So, I have 3 levels -5, 10, 15 for time.
> The new dataset is run using the model after fitting the observation
> level random effects:
> 
> 
> Behavdat2$resid<-as.factor(1:dim(Behavdat2)[1])
> 
> lmer gives error message when I used family=quasibinomial.

  Yes.  See the last paragraph of my comments below.
> 
> 
> Response2 <- 5-Behavdat2$Response 
> Response3 <-cbind(Response,Response2)
> 
> My model:
> 
> (fm1<-lmer(Response3~Wavelength*(Start_Resp/Time)+(1|Subject)+(1|resid),
> family=binomial,data=Behavdat2))
> 
> I didn't use the Strain as only a subset of the data with 1 level of
> strain was analyzed.
> 
> In the model, I specified Time as being nested within Start_Resp.  I
> think Time could be nested under Subject, but I would like to run
> Time as a fixed effect with 4 levels and my aim will  be to look for
> whether there is any change in light/dark response to various
> wavelengths (3 levels used - Red, Blue, Green) as time changes.
> 
> Is it a correct model to find the changes in light response in
> different time levels for different wavelengths?  Do I need to have
> an interaction model with Wavelength*Start_Resp*Time, given that the
> correlations within subject will be taken care by 1|resid statement
> in the model?
> 
> I appreciate your comments. Thanks AK
> 
> 
> My  dataset: Number Wavelength Strain Subject Start_Resp time
> Response 1 Red G 1 L 5 3 2 Red G 1 L 10 0 3 Red G 1 L 15 0
> 
> 
> 
> 
> 
> 
> 
> 5 Red G 2 L 5 3 6 Red G 2 L 10 3 7 Red G 2 L 15 3
> 
> ===== ======
> 
> ===============
> 
> ====-------------------------------------------
> 
> 
> 
> 
> 
> From: Ben Bolker <bbolker at gmail.com> To:
> r-sig-mixed-models at r-project.org Sent: Tuesday, April 24, 2012 6:14
> PM Subject: Re: [R-sig-ME] Nested longitudinal
> 
> Robert Kushler <kushler at ...> writes:
> 
>> 
>> 
>> I'm breaking my own rules by offering answers before getting
>> answers to my own questions.
> 
> Hmmm.  I didn't see any outstanding unanswered questions from you on
> the list -- did I miss some?  Or have you not asked them yet?
> 
>> 1) The simple general answer to your basic question is that the
>> "/" character is used to represent nested factors, while a "*" is
>> used to indicate "crossed" factors that might interact.  Your use
>> of "+" in the model formula below means you are *assuming* that the
>> factors do not interact - and assuming something doesn't make it
>> true.
> 
> [snip]
> 
>> 4) It seems to me that there will be very strong serial
>> correlation in the 15 measurements on an individual subject.
>> Unfortunately lmer doesn't include this as a modeling option.  Your
>> current syntax does one of the following: (a) fits a linear "time
>> effect" with a random slope and intercept or (b) if time is a
>> factor you are trying to estimate an "unstructured" (sorry, Doug
>> Bates) 15 by 15 matrix. Both approaches are problematic here.  I
>> suggest you collapse the data by time and record y = number of
>> minutes (out of 15) spent in light (or dark - doesn't matter)
>> areas, and then use "cbind(y,15-y)" as the response.  You probably
>> should also try some alternatives to the binomial family (e.g.,
>> "quasibinomial").
> 
> lmer doesn't include it at least in part because it's difficult to 
> fit into a conditional GLMM framework -- the most sensible way to 
> define this would be to allow a per-observation random effect (which 
> would then also make the binomial overdispersed by definition), and 
> then specify that the individual-level random effects were
> themselves temporally correlated.
> 
> Usually when you find packages that can incorporate temporal 
> correlation in a GLMM they either require that you specify the full 
> statistical model yourself (WinBUGS/JAGS, AD Model Builder), *or*
> they are in some sense marginal models (glmmPQL in the MASS package
> [I know] and ASREML [I think] allow 'R-side' structures such as
> temporal correlation in GLMMs, but they use penalized
> quasi-likelihood for estimation, which may under some circumstances
> be problematic).
> 
> For what it's worth you can't use quasi- families in glmer(); you can
> either add an individual-level random effect, or use MASS::glmmPQL if
> you want quasi- (see http://glmm.wikidot.com/faq for more discussion
> of overdispersion in GLMMs).
> 
> 
>> Regards,   Rob Kushler
>> 
>> On 4/24/2012 12:40 AM, arun wrote:
> 
> [snip]
> 
>>> A brief introduction about the work: It is a light/dark choice 
>>> test conducted in insect larvae.  The response is binary (0- 
>>> present in dark area, 1-present in light area) and the
>>> experiment is run for 15 min, so there are 15 measurements per
>>> individual larva at 1 min intervals.  The factors which affect
>>> this study are Strain (2 levels-G and S), wavelength of light (4
>>> levels-blue, green, UV, red), and starting response at 0 min (two
>>> levels- animal present in dark-D or light-L).  This is how I
>>> think it is nested.  Strain nested inside Wavelength, Subject
>>> (individual) nested within strain, Starting response within
>>> subject, and time > within Starting response.  The data looks
>>> like this:
> 
> [snip]
> 
>>> 
>>> 
>>> (fm2<-lmer(Response~Wavelength+Startingresponse+Strain+ time +
>> (time|Subject),family=binomial, data=Behavdat))
>>> 
>>> I am not sure how to specify the nested structures within the 
>>> model.
> 
> _______________________________________________ 
> R-sig-mixed-models at r-project.org mailing list 
> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models



More information about the R-sig-mixed-models mailing list