[R-sig-ME] MCMCglmm ignores custom contrasts in multi-response models

[rafter sass ferguson]

MCMCglmm appears to consistently ignore set contrasts and revert to default
contrasts when fitting a multivariate model (and not when fitting a
univariate model).

Reproducible example:

### recreating the data from the Course Notes multi-response vignette
require(MCMCglmm)

set.seed(102)
n_indiv <- 500
n_per_indiv <- 4
n_tot <- n_indiv * n_per_indiv
id <- gl(n_indiv, n_per_indiv)
av_wealth <- rlnorm(n_indiv, 0, 1)
ac_wealth <- av_wealth[id] + rlnorm(n_tot, 0, 1)
av_ratio <- rbeta(n_indiv, 10, 10)
ac_ratio <- rbeta(n_tot, 2 * av_ratio[id], 2 * (1 - av_ratio[id]))
y.car <- (ac_wealth * ac_ratio)^0.25
y.hol <- (ac_wealth * (1 - ac_ratio))^0.25

Spending <- data.frame(y.hol, y.car, id)

### adding an ordered factor variable 'education'
Spending$edu <- factor(sample(c("1_HS", "2_2YrCol", "3_4YrCol", "4_Grad"),
size=n_indiv, replace=TRUE), levels=c("HS", "2YrCol", "4YrCol", "Grad"),
ordered=T)
### setting custom contrast
contrasts(Spending$edu) <- contr.helmert(4)

### fit samemodels to vignette, but adding education as predictor
m1_id <- MCMCglmm(y.car ~ y.hol + edu, random = ~id, data = Spending,
verbose = FALSE)
m1_idyr <- MCMCglmm(cbind(y.hol, y.car) ~ trait - 1 + trait:edu, random =
~us(trait):id,
                    rcov = ~us(trait):units, data = Spending, family =
c("gaussian", "gaussian"),
                    verbose = T)

summary(m1_id) ### Helmert contrasts retained for univariate
summary(m1_idyr) ### reverts to default polynomial contrasts for
multivariate

### making sure it's not specific to ordered factors - change edu to an
unordered
### factor with custom contrasts
Spending$edu <- factor(Spending$edu, ordered=F)
contrasts(Spending$edu) <- contr.treatment(4, base=3)

### re-fit models
m1_id2 <- MCMCglmm(y.car ~ y.hol + edu, random = ~id, data = Spending,
verbose = FALSE)
m1_idyr2 <- MCMCglmm(cbind(y.hol, y.car) ~ trait - 1 + trait:edu, random =
~us(trait):id,
                    rcov = ~us(trait):units, data = Spending, family =
c("gaussian", "gaussian"),
                    verbose = F)

### same results:
summary(m1_id2) ### custom contrasts (base=3) retained for univariate
summary(m1_idyr2) ### reverts to default (base=1) for multivariate






Rafter Sass Ferguson, MS
PhD Candidate | Crop Sciences Department
University of Illinois in Urbana-Champaign
liberationecology.org
518 567 7407

	[[alternative HTML version deleted]]