[R-sig-ME] MCMCglmm_2.13

[Malcolm Fairbrother]

Dear Jarrod,

Can I please just confirm what you meant here by "covariance structures for the random effects"? If I understand correctly, this technique might be very useful for quite a lot of things--though I'm also trying to think it through.

Specifically, the nlme package for example allows for error structures like AR1, etc., in the *residuals* of a mixed model. Are you saying that MCMCglmm's "ginverse" function allows for similar sorts of structure to be modelled in the *random effects*? Consider:

library(nlme); library(MCMCglmm)

# models with uncorrelated random effects (16 Rats observed on 11 occasions):
m1 <- lme(weight ~ Time, random = ~ 1 | Rat, BodyWeight)
m2 <- MCMCglmm(weight ~ Time, random = ~ Rat, data=BodyWeight, verbose=F) # m1 and m2 are the same

# same model as the above, but allowing for AR1 residuals errors:
cs1AR1 <- corAR1(0.8, form = ~ Time | Rat)
cs1AR1 <- Initialize(cs1AR1, BodyWeight)
corMatrix(cs1AR1)[[1]] # shows the correlation structure for each Rat (all 16 are the same, since each Rat is measured on the same occasions)
m3 <- lme(weight ~ Time, random = ~ Time | Rat, correlation = cs1AR1, BodyWeight)

# now, for an illustrative experiment, treat *Time* as the random grouping factor, instead of Rat:
m4 <- lme(weight ~ Diet, random = ~ 1 | Time, BodyWeight)
m5 <- MCMCglmm(weight ~ Diet, random=~Time, data=BodyWeight, verbose=F, nitt=53000, pr=T)
# estimates for the variance of the random intercepts is somewhat different, but otherwise m4 and m5 yield similar results

# now, create a "ginverse" matrix, and incorporate it in a call to MCMCglmm:
mat <- as(corMatrix(cs1AR1)[[1]], "sparseMatrix")
rownames(mat) <- unique(BodyWeight$Time)
m6 <- MCMCglmm(weight ~ Diet, random=~Time, ginverse=list(Time=mat), data=BodyWeight, verbose=F, nitt=53000, pr=T)

The variance for the random effects for Time is now substantially reduced. Does this mean that the "estimates" for the U_j terms has taken into account their correlation, whereas m5 treated them as independent? (Maybe I should have inverted "mat" before including it in the call?)

Do you see what I'm getting at? Maybe I've not understood, but this seemed to raise the possibility of allowing for random effects to be structured in many interesting ways (serially in time, spatially, etc.), which I don't think other packages can handle…

Many thanks for any clarification, as always,
Malcolm




> Date: Mon, 25 Jul 2011 12:41:43 +0100
> From: Jarrod Hadfield <j.hadfield at ed.ac.uk>
> To: r-sig-mixed-models at r-project.org
> Subject: [R-sig-ME] MCMCglmm_2.13
> 
> Hi,
> 
> A new version of MCMCglmm (MCMCglmm_2.13) has been submitted to CRAN  
> and should be on-line soon.
> 
> NEW FEATURES:
> 
> a) Covariance structures for the random effects of the form  
> kronecker(V,A) where A is some arbitrary matrix can now be fitted.  
> Previous versions only allowed A matrices defined by a pedigree or  
> phylogeny, and even then only single phylogenies/pedigrees could be  
> passed. The A matrices (actually their sparse inverses) are passed via  
> the argument ginverse.
> 
> b) c. 10% increase in speed for phylogenetic and pedigree models.
> 
> BUG FIXES:
> 
> i) If a column sum of Z was zero the associated term could be  
> incorrectly omitted from the mixed model equations. This is unlikely  
> to affect the vast majority of users.
> 
> ii) If a prior was specified with too many structures then the  
> redundant prior structures were ignored. An error message is now issued.
> 
> iii) MCMCglmm failed with a cryptic error message if the internal  
> nodes of a phylogeny were labelled but not uniquely. A more helpful  
> error message is now given.
> 
> Cheers,
> 
> Jarrod