For my current postdoc I have the opportunity to do some quantitative genetics with pedigree data in the animal model, something I’ve never done (having only used data from hemiclones/chromosome substitution line experiments before). I spent a couple of days learning how, and I turned to “An Ecologist’s Guide to the Animal Model” (Wilson et al 2009) for assistance. The tutorial included in the supplementary material (see revised version here, not the original supplementary material) proved to be an excellent and easy introduction. Being familiar with MCMCglmm, I used that as the environment for me to learn in, and I highly recommend this tutorial for anyone learning about the animal model.

For myself, I’m interested in the evolution of sex differences, which naturally calls for estimates of intersexual genetic correlations and G matrices. I decided to go beyond the tutorial and construct G matrices for the *Gryphon* data, which is included with the guide. To do this I had to generate sex-specific trait columns, and use a model with four response variables (female body weight [BWT_F], male body weight [BWT_M], female tarsus length [TAR_F], male tarsus length [TAR_M]). From this model I am able to estimate additive genetic variance, within-sex cross-trait covariance, cross-sex within-trait covariance, and cross-sex cross-trait covariance; these are the components of a G matrix. I can also derive correlations (e.g. intersexual genetic correlations, *r*_{MF}), and within- and cross-sex heritability (see Lehtovaara et al 2013 Am Nat). The following annotated R script does all of this, and run for 1.3M iterations it produces the G matrix above (note that the numbers above the diagonal of G, the shaded cells, are the genetic correlations; the sex-specific additive genetic variances are in bold). I could also derive an M matrix; one with the variance and covariance explained by the maternal effects. Note, I haven’t worried about priors here (see the MCMCglmm documentation for advice), this was an exercise in extending the model and extracting intersexual genetic correlation /G matrices. It’s a work in progress – I would greatly appreciate feedback/corrections!

# Version 12-01-2017
# Load Package
library("MCMCglmm")
# Load Data (download the data and change the filepaths to suit)
Data = read.table("\\gryphon.txt", header = T)
Ped = read.table("\\gryphonped.txt", header = T)
# Prep Data
colnames(Data)[1] = "animal"
Data$animal = as.factor(Data$animal)
Data$MOTHER = as.factor(Data$MOTHER)
Data$BYEAR = as.factor(Data$BYEAR)
Data$SEX = as.factor(Data$SEX)
Data$BWT = as.numeric(Data$BWT)
Data$TARSUS = as.numeric(Data$TARSUS)
for(x in 1:3) Ped[,x] = as.factor(Ped[,x])
# Make Sex-Specific Traits
Data$BWT_M = ifelse(Data$SEX == 2, Data$BWT, NA)
Data$BWT_F = ifelse(Data$SEX == 1, Data$BWT, NA)
Data$TAR_M = ifelse(Data$SEX == 2, Data$TARSUS, NA)
Data$TAR_F = ifelse(Data$SEX == 1, Data$TARSUS, NA)
# Construct Priors (see Wilson et al 2009 and MCMCglmm documentation for guidance on priors)
phen.var = matrix(c(var(Data$BWT_F, na.rm = TRUE), 0, 0, 0,
0, var(Data$BWT_M, na.rm = TRUE), 0, 0,
0, 0, var(Data$TAR_F, na.rm = TRUE), 0,
0, 0, 0, var(Data$TAR_M, na.rm = TRUE)
),4,4)
prior4.1= list(G = list( G1 = list(V = phen.var/4, n = 2),
G2 = list(V = phen.var/4, n = 2),
G3 = list(V = phen.var/4, n = 2)),
R = list(V = phen.var/4, n = 2))
# Model
model4.1 = MCMCglmm(cbind(BWT_F, BWT_M, TAR_F, TAR_M) ~ trait-1,
random = ~us(trait):animal + us(trait):BYEAR + us(trait):MOTHER,
rcov = ~us(trait):units,
family = rep("gaussian", times = 4),
pedigree = Ped,
data = Data,
prior = prior4.1)
# Additive Genetic Variance (Diagonals of G)
posterior.mode(VA_BWT_F <- model4.1$VCV[,"traitBWT_F:traitBWT_F.animal"])
posterior.mode(VA_BWT_M <- model4.1$VCV[,"traitBWT_M:traitBWT_M.animal"])
posterior.mode(VA_TAR_F <- model4.1$VCV[,"traitTAR_F:traitTAR_F.animal"])
posterior.mode(VA_TAR_M <- model4.1$VCV[,"traitTAR_M:traitTAR_M.animal"])
# Cross Trait Covariance (Off-diagonals of GM and GF)
posterior.mode(COV_BWT_F.TAR_F <- model4.1$VCV[,"traitBWT_F:traitTAR_F.animal"])
posterior.mode(COV_BWT_M.TAR_M <- model4.1$VCV[,"traitBWT_M:traitTAR_M.animal"])
# Cross Sex Covariance (Diagonals of B)
posterior.mode(COV_BWT_F.BWT_M <- model4.1$VCV[,"traitBWT_F:traitBWT_M.animal"])
posterior.mode(COV_TAR_F.TAR_M <- model4.1$VCV[,"traitTAR_F:traitTAR_M.animal"])
# Cross Trait Cross Sex Covariance (Off-diagonals of B)
posterior.mode(COV_BWT_F.TAR_M <- model4.1$VCV[,"traitBWT_F:traitTAR_M.animal"])
posterior.mode(COV_BWT_M.TAR_F <- model4.1$VCV[,"traitBWT_M:traitTAR_F.animal"])
# Cross Trait Genetic Correlation
posterior.mode(COR_BWT_F.TAR_F <- COV_BWT_F.TAR_F/sqrt(VA_BWT_F * VA_TAR_F))
posterior.mode(COR_BWT_M.TAR_M <- COV_BWT_M.TAR_M/sqrt(VA_BWT_M * VA_TAR_M))
# Cross Sex Genetic Correlation (AKA Intersexual Genetic Correlation, or rMF)
posterior.mode(COR_BWT_F.BWT_M <- COV_BWT_F.BWT_M/sqrt(VA_BWT_F * VA_BWT_M))
posterior.mode(COR_TAR_F.TAR_M <- COV_TAR_F.TAR_M/sqrt(VA_TAR_F * VA_TAR_M))
# Cross Trait Cross Sex Genetic Correlation
posterior.mode(COR_BWT_F.TAR_M <- COV_BWT_F.TAR_M/sqrt(VA_BWT_F * VA_TAR_M))
posterior.mode(COR_BWT_M.TAR_F <- COV_BWT_M.TAR_F/sqrt(VA_BWT_F * VA_TAR_M))
# Within-Sex Heritability (D or S = Dam or Sire, D or S = Daughter or Son)
posterior.mode(DDH2_BWT_F <- model4.1$VCV[,"traitBWT_F:traitBWT_F.animal"] / (model4.1$VCV[,"traitBWT_F:traitBWT_F.animal"] + model4.1$VCV[,"traitBWT_F:traitBWT_F.MOTHER"] + model4.1$VCV[,"traitBWT_F:traitBWT_F.BYEAR"] + model4.1$VCV[,"traitBWT_F:traitBWT_F.units"]))
posterior.mode(SSH2_BWT_M <- model4.1$VCV[,"traitBWT_M:traitBWT_M.animal"] / (model4.1$VCV[,"traitBWT_M:traitBWT_M.animal"] + model4.1$VCV[,"traitBWT_M:traitBWT_M.MOTHER"] + model4.1$VCV[,"traitBWT_M:traitBWT_M.BYEAR"] + model4.1$VCV[,"traitBWT_M:traitBWT_M.units"]))
posterior.mode(DDH2_TAR_F <- model4.1$VCV[,"traitTAR_F:traitTAR_F.animal"] / (model4.1$VCV[,"traitTAR_F:traitTAR_F.animal"] + model4.1$VCV[,"traitTAR_F:traitTAR_F.MOTHER"] + model4.1$VCV[,"traitTAR_F:traitTAR_F.BYEAR"] + model4.1$VCV[,"traitTAR_F:traitTAR_F.units"]))
posterior.mode(SSH2_TAR_M <- model4.1$VCV[,"traitTAR_M:traitTAR_M.animal"] / (model4.1$VCV[,"traitTAR_M:traitTAR_M.animal"] + model4.1$VCV[,"traitTAR_M:traitTAR_M.MOTHER"] + model4.1$VCV[,"traitTAR_M:traitTAR_M.BYEAR"] + model4.1$VCV[,"traitTAR_M:traitTAR_M.units"]))
# Cross-Sex Heritability (D or S = Dam or Sire, D or S = Daughter or Son)
posterior.mode(DSH2_BWT <- DDH2_BWT_F * COR_BWT_F.BWT_M)
posterior.mode(SDH2_BWT <- SSH2_BWT_M * COR_BWT_F.BWT_M)
posterior.mode(DSH2_TAR <- DDH2_TAR_F * COR_TAR_F.TAR_M)
posterior.mode(SDH2_TAR <- SSH2_TAR_M * COR_TAR_F.TAR_M)
# G matrix
G = matrix(c( posterior.mode(VA_BWT_F), posterior.mode(COR_BWT_F.TAR_F), posterior.mode(COR_BWT_F.BWT_M), posterior.mode(COR_BWT_M.TAR_F),
posterior.mode(COV_BWT_F.TAR_F), posterior.mode(VA_TAR_F), posterior.mode(COR_BWT_F.TAR_M), posterior.mode(COR_TAR_F.TAR_M),
posterior.mode(COV_BWT_F.BWT_M), posterior.mode(COV_BWT_M.TAR_F), posterior.mode(VA_BWT_M), posterior.mode(COR_BWT_M.TAR_M),
posterior.mode(COV_BWT_F.TAR_M), posterior.mode(COV_TAR_F.TAR_M), posterior.mode(COV_BWT_M.TAR_M), posterior.mode(VA_TAR_M))
, ncol = 4, byrow=T)
colnames(G) = rownames(G) = c("BWT_F", "TAR_F", "BWT_M", "TAR_M")

### Like this:

Like Loading...