[R-sig-ME] differing variances within different random effects levels

Jarrod Hadfield j.hadfield at ed.ac.uk
Mon Feb 14 23:18:21 CET 2011


Hi Ben,

In MCMCglmm the appropriate model (without random slopes) would look like:

random=~idh(species):indiv

which fits species specific individual variances. The use of idh  
(rather than us) forces the covariances between individuals in  
different species to zero, which is appropriate given the nested design.

A clumsier way to specify the model would be:

random=~us(at.level(species, 1)):indiv+us(at.level(species, 2)):indiv +  ....

The "at.level" function used to work with lmer (at least in the past):

(at.level(species, 1)|indiv)+(at.level(species, 2)|indiv) + ...

For random interaction models I think the at.level has to be used in MCMCglmm:

random=~us(at.level(species, 1):treatment):indiv+us(at.level(species,  
2):treatment):indiv +  ....

not elegant, and probably too many parameters! Perhaps  
double-hierarchical models would be a better bet if S is large?

Cheers,

Jarrod


Quoting Ben Bolker <bbolker at gmail.com>:

>
>   This is a question that's been bothering me for a while.  The answer
> might be "read sections 4.2.2 and 4.2.3 of Pinheiro and Bates 2000
> again, more carefully", but if anyone has done this before I wouldn't
> mind some hints ...
>
>   Consider a simplified version of the previous poster's example.
> Suppose I have a randomized block design where T treatments are done for
> each of N individuals within each of S species. My basic model for this
> would be
>
>   response~treatment+(treatment|species/indiv)
>
> (or fixed=response~treatment, random=~treatment|species/indiv in lme)
>
> which would give me random effects
>
> 1|species
> 1|indiv:species
> treatment|species
> treatment|indiv:species
>
> (and correlations between intercept and treatment effects at species and
> indiv:species level).
>
>   Suppose now that I think the *variance* among individuals varies among
> species.  As far as I can tell, lme/lmer assume that this variance is
> constant across species.  If I wanted variance across *observations*
> (i.e. residual variance) to vary across species I could do it via
> something like weights=varIdent(~species), but this only works for
> residual variance.
>
>   I suspect I could work out some way to do this with MCMCglmm, and I
> don't have an immediate need for it, but I thought I would ask ...
>
>  thanks
>     Ben Bolker
>
> _______________________________________________
> R-sig-mixed-models at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>
>



-- 
The University of Edinburgh is a charitable body, registered in
Scotland, with registration number SC005336.




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