[R-sig-ME] Multiresponse MCMCglmm model specification and questions

Rafael Maia queirozrafaelmv at yahoo.com.br
Fri Jan 17 18:05:16 CET 2014


Dear list members,

I am attempting to run a multiresponse model to investigate factors influencing the evolution of a trait measured in both males and females of multiple species. Therefore I chose a multiresponse approach using MCMCglmm, but my model has a series of little caveats (don’t they all?) that I’d like to make sure I’m addressing correctly. Here is a brief description of the model I’m trying to fit:

Ymale, Yfemale: the same response variable, measured in males and females

X1: a variable (factor, 2 levels) that I expect to affect both sexes differently
X2: a variable (factor, 2 levels) that I expect to affect both sexes equally (or estimate the joint effect across sexes)
X3: a continuous variable measured only for males (which I therefore expect to only have an effect in male measurements)

so I am currently specifying the model as below:

mod <- cbind(Ymale, Yfemale) ~ trait*X1 + X2 + at.level(trait, 1):X3

ainv <- inverseA(phylo)$Ainv

Prior <- list(R=list(V=diag(2), nu=0.002), G=list(G1=list(V=diag(2), nu=0.002)))

result <- MCMCglmm(mod, 
    random = ~us(trait):spp, rcov = ~idh(trait):units, 
    data=dat, prior = Prior, family = rep(‘gaussian’, 2),
    nitt=2100000, burnin=100000, thin=1000)

with that, I have a couple questions (besides: does this model specification look reasonable?):

1. Prior: I wanna make sure I’m specifying an inverse-gamma with scale=shape=0.001 for the variance components. Based on footnote 1 of the Course Notes (p.102) I think I am based on the changes, but I’ve seen several times in mailing lists responses and other papers nu=1.002 being used for a matrix of same dimensions. My results are robust to either specification but I’d like to make sure I’m doing what I think I’m doing, given that there have been changes on this regard.

2. I am calculating the correlation between ymale and female as:
    result$VCV[,2]/sqrt(result$VCV[,1]*result$VCV[,4])
as described in the tutorial vignette. My values, however, are extremely high. The two traits are indeed strongly correlated among species, but I’m getting credible intervals of 0.96 - 0.99 (for comparison, a raw non-phylogenetic correlation of the traits is 0.75). Am I forgetting to include any variance components to calculate the intersexual correlation? Is this increase in the correlation expected given the inclusion of the fixed and random effects in the model?

3. I have a hypothesis that the correlation between male and female traits should be stronger for one level of X1 than for the other. Is it possible to specify random terms in order to calculate (co)variances of the traits conditional on X1?
    
I deeply appreciate any help you can provide. Cheers!


Abraços,
Rafael Maia
---
http://www.rafaelmaia.net/
PhD Candidate, Integrated Bioscience
University of Akron
"A little learning is a dangerous thing; drink deep, or taste not the Pierian spring." (A. Pope)



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