[R-sig-ME] MCMCglmm Contrasting h2 in univariate and bivariate model

[Cameron So]

Hello MCMCglmm users,

I am estimating the Va and heritability of a trait in separate univariate and bivariate models.
Some background beforehand: In my experiment, individuals are exposed to separate treatments.
The pedigree for each treatment is quite similar given I am using plants + multiple replicates.
Therefore, I have two univariate models: a subset of the dataset for each treatment. This has allowed
me to compare the Va and heritability estimates between treatments using the univariate models.
In the multivariate model, I include the entire dataset so I can evaluate cross-environment COVa.

I am noticing that my heritability estimates between the univariate and bivariate models differ
substantially. Could anyone offer any advice on what I am observing?

NOTE: I am not having issues of autocorrelation, nor convergence. Additionally, the univariate
Models include a second random effect term (matID) but this shouldn�t affect the additive (animal)
term that drastically. The priors obviously do differ between the univariate and bivariate models
because of the binary character of �treatment�.

Below are my models:

## UNIVARITE MODEL EXAMPLE ##

prior8.3b <- list(R = list(V = 1, nu = 0.002),
                 G = list(G1 = list(V = 1, nu = 1, alpha.mu = 0, alpha.V = 1000),
                          G2 = list(V = 1, nu = 1, alpha.mu = 0, alpha.V = 1000)))

HW_model8.3b <- MCMCglmm(leaf ~ plot, random = ~animal + matID,
                     ginverse = list(animal = Ainv),
                     family = "gaussian", data = heated.leaf, prior = prior8.3b,
                     nitt = 2100000, thin = 1000, burnin = 100000, verbose = T, pr = TRUE)

HW_herit8.3b <- HW_model8.3b$VCV[, "animal"]/(HW_model8.3b$VCV[, "animal"] + HW_model8.3b$VCV[, "matID"] + HW_model8.3b$VCV[, "units"])
mean(HW_herit8.3b)
#h2 = ~ 0.176

## BIVARIATE MODEL EXAMPLE ##


prior6.2 <- list(R = list(V = diag(1), fix = 1),
                 G = list(G1 = list(V = diag(2), nu = 2, alpha.mu = c(0,0), alpha.V = diag(c(1,1)))))

PL_model6.2 <- MCMCglmm(leaf ~ plot:treatment + treatment, random = ~us(treatment):animal,
                        ginverse = list(animal = Ainv), rcov=~units,
                        family = "gaussian", data = plastic.leaf, prior = prior6.2,
                        nitt = 4100000, thin = 2000, burnin = 100000, verbose = T, pr = TRUE, trunc = TRUE)

herit_PL9.2_A <-PL_model9.2$VCV[,'treatmentA:treatmentA.animal']/  (PL_model9.2$VCV[,'treatmentA:treatmentA.animal'] + 1)
mean(herit_PL9.2_A)
#h2 = ~0.73


Thank you for any help in advance!

Cameron



	[[alternative HTML version deleted]]