####=========================================####
#### For CRAN time limitations most lines in the
#### examples are silenced with one '#' mark,
#### remove them and run the examples using
#### command + shift + C |OR| control + shift + C
####=========================================####
# data(DT_gryphon)
# DT <- DT_gryphon
# A <- A_gryphon
# P <- P_gryphon
# #### look at the data
# head(DT)
# #### fit the model with no fixed effects (intercept only)
# 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
# Ai <- solve(A + diag(1e-4,ncol(A),ncol(A)))
# Ai <- as(as(as( Ai, "dMatrix"), "generalMatrix"), "CsparseMatrix")
# ####====================####
# #### multivariate model ####
# #### 2 traits ####
# ####====================####
# head(DT)
# DT2 <- stackTrait(DT, traits = c("BWT","TARSUS"))
# head(DT2$long)
#
# # #### fit the multivariate model with no fixed effects (intercept only)
# mix2 <- mmes(valueS~trait, henderson=FALSE, # using direct inversion
# random=~vsm(usm(trait),ism(ANIMAL),Gu=A),
# rcov=~vsm(dsm(trait),ism(units)),
# data=DT2$long)
# summary(mix2)$varcomp
# cov2cor(mix2$theta[[1]])
Run the code above in your browser using DataLab