enhancer (version 1.1.0)

DT_gryphon: Gryphon data from the Journal of Animal Ecology

Description

This is a dataset that was included in the Journal of animal ecology by Wilson et al. (2010; see references) to help users understand how to use mixed models with animal datasets with pedigree data.

The dataset contains 3 elements:

gryphon; variables indicating the animal, the mother of the animal, sex of the animal, and two quantitative traits named 'BWT' and 'TARSUS'.

pedi; dataset with 2 columns indicating the sire and the dam of the animals contained in the gryphon dataset.

A; additive relationship matrix formed using the 'getA()' function used over the pedi dataframe.

Usage

data("DT_gryphon")

Format

The format is: chr "DT_gryphon"

References

Wilson AJ, et al. (2010) An ecologist's guide to the animal model. Journal of Animal Ecology 79(1): 13-26.

Giovanny Covarrubias-Pazaran (2024). lme4breeding: enabling genetic evaluation in the age of genomic data. To be submitted to Bioinformatics.

Douglas Bates, Martin Maechler, Ben Bolker, Steve Walker (2015). Fitting Linear Mixed-Effects Models Using lme4. Journal of Statistical Software, 67(1), 1-48.

Examples

data(DT_gryphon)
DT <- DT_gryphon
A <- A_gryphon
P <- P_gryphon
#### look at the data
head(DT)

# \donttest{

############## sommer #################

if(requireNamespace("sommer")){
library(sommer)
mix1 <- mmes(BWT~1,
             random=~vsm(ism(ANIMAL),Gu=A),
             rcov=~units,
             data=DT)
summary(mix1)$varcomp

## mmes algorithm uses the inverse of the relationship matrix
## if you select henderson=TRUE
Ai <- solve(A + diag(1e-4,ncol(A),ncol(A)))
Ai <- as(as(as( Ai,  "dMatrix"), "generalMatrix"), "CsparseMatrix")
attr(Ai, 'inverse')=TRUE
####====================####
#### multivariate model ####
####     2 traits       ####
####====================####
head(DT)

traits <- c("BWT","TARSUS")
DT[,traits] <- apply(DT[,traits],2,scale)
DTL <- reshape(DT[,c("ANIMAL","MOTHER","BYEAR","SEX", traits)],
               idvar = c("ANIMAL","MOTHER","BYEAR","SEX"),
               varying = traits,
               v.names = "value", direction = "long",
               timevar = "trait", times = traits )
DTL <- DTL[with(DTL, order(trait,ANIMAL)), ]
head(DTL)

# #### fit the multivariate model with no fixed effects (intercept only)
mix2 <- mmes(value~trait, # henderson=TRUE,
             random=~vsm(usm(trait),ism(ANIMAL),Gu=A),
             rcov=~vsm(dsm(trait),ism(units)),
             data=DTL)
summary(mix2)$varcomp
cov2cor(mix2$theta[[1]])

}

############## lme4breeding #################
if(requireNamespace("lme4breeding")){
library(lme4breeding)
## fit the model with no fixed effects (intercept only)
mix1 <- lmeb(BWT~ (1|ANIMAL),
                 relmat = list(ANIMAL=A),
                 data=DT)
vc <- VarCorr(mix1); print(vc,comp=c("Variance"))
sigma(mix1)^2 # error variance
BLUP <- ranef(mix1, condVar=TRUE)
condVAR <- lapply(BLUP, function(x){attr(x, which="postVar")}) # take sqrt() for SEs


### multi-trait model
traits <- c("BWT","TARSUS")
for(iTrait in traits){DT[,iTrait] <- scale(imputev(DT[,iTrait]))}
DTL <- reshape(DT[,c("ANIMAL", traits)], idvar = "ANIMAL", varying = traits,
               v.names = "value", direction = "long",
               timevar = "trait", times = traits )
DTL <- DTL[with(DTL, order(trait)), ]
head(DTL)

system.time(
  mix <- lmeb(value ~ (0+trait|ANIMAL),
                  relmat = list(ANIMAL=A),
                  # rotation = TRUE, 
                  data=DTL)
)
vc <- VarCorr(mix); print(vc,comp=c("Variance"))
cov2cor(vc$ANIMAL)
sigma(mix)^2 # error variance

}

# }