[R-sig-ME] mixed models with two random variables?

[Malcolm Fairbrother]

Hi Barbara,

The lme4 package may be a bit easier for you to use, and is more current than nlme in various ways (though nlme does do a few things lme4 can't do). I think the problem is the right hand side of the random part of your model: I don't believe you can have a "+" there, as you have in your first two models. Your third model definitely doesn't make sense, because you're treating fReserve as both a covariate and a random classification.

It's a bit hard to understand your query generally. Variables treated as fixed effects can be either categorical or continuous--not a problem. But I don't understand the random classification in your case, which is the part for which you probably want to calculate a variance term (and why you need a mixed model). A random classification is usually coded as 1 to 20, or A to Z, or something like that, whereas your Roughness and DivBoulders look like continuous variables.

However, if you are trying to model observations cross-classified in both Roughness and DivBoulders, and treating fReserve as the only covariate/predictor, then try lme4's lmer, along the lines of:

lmer(Biomass ~ fReserve + (1 | DivBoulders) + (1 | Roughness), myData)

Hope that helps. Follow up with further clarification if not.

Cheers,
Malcolm



> Date: Fri, 14 Sep 2012 18:34:28 +0100
> From: barbara costa <rbarbarahc at gmail.com>
> To: r-sig-mixed-models at r-project.org
> Subject: [R-sig-ME] mixed models with two random variables?
> 
> Dear R users,
> Does anyone knows how to run a glmm with one fixed factor and 2 random
> numeric variables (indices)? Is there any way to force in the model a
> separate interaction of those random variables with the fixed one?
> I hope you can help me.
> 
> #eg.
> Reserve <- rep(c("In","Out"), 100)
> fReserve <- factor(Reserve)
> DivBoulders <- rep (c(1.23,2.4,1.26,1.78,1.97,1.35,1.23,2.4,1.26,1.78), 20)
> Roughness <- rep(c(3.45,2.56,1.32,5.67,3.73,3.57,2.66,1.52,7.67,2.73),20)
> Biomass <- rep(c(8,5.3,3.5,12,25.4,10.1,9.8,2.4,5.6,5.3),20)
> 
> myData <- data.frame (fReserve ,DivBoulders ,Roughness ,Biomass )
> 
> #glm
> glm1 <- glm (Biomass ~ fReserve * Roughness + fReserve *  DivBoulders  ,
> family= Gamma, data= myData)
> 
> #glmm:
> library (nlme)
> 
> lme1 <- lme (Biomass  ~  fReserve , random= ~1 + fReserve | Roughness
> +  DivBoulders  , data=  myData ) # random intercept and slope  - I suspect
> it's my case
> #Error in getGroups.data.frame(dataMix, groups) :
> # Invalid formula for groups
> # if I only use one random variable I have:
> #Error in chol.default((value + t(value))/2) :
>  #the leading minor of order 2 is not positive definite
> 
> 
> lme2 <- lme( Biomass  ~  fReserve , random = ~1 | Roughness
> +  DivBoulders  ,data=myData) #random intercept
> #Error in getGroups.data.frame(dataMix, groups) :
> # Invalid formula for groups
> # if I only use one random variable my result is fine!
> 
> 
> lme3 <- lme (Biomass   ~  fReserve  , random=  ~ Roughness +   DivBoulders
> | fReserve   , data= myData ) # from help (lme)
> summary (lme3)
> # I have a result. Is the model correct for what I want?
> # BUT my fixed effect (reserve) does not have a p-value due to zero degrees
> of freedom. However in glm1 or lm2 with one random variable it has.
> 
> Can you help me please?
> 
> Thanks a lot in advance,
> Barbara