This note instructs how to reproduce the genetics (Table 3) simulation example in the manuscript. This simulation compares the performance among MM, EM, FS (Fisher scoring) and lme4 under model: $$ \mathbf{y} \sim \text{Normal}(\mathbf{X} \beta, \sigma_a^2 \Phi + \sigma_e^2 \mathbf{I}). $$ Here the covariates $\mathbf{X}$ and kinship matrix $\Phi$ are from the example in Option 29 of the Mendel software. In simulation $\beta = [33, 6.0]$ and $\sigma_e^2$ is fixed at 2.0. $\sigma_a = r \cdot \sigma_e^2$, where $r$ varies in {0, 0.05, 0.1, 1, 10, 20}.
Comparison metrics include:
versioninfo()
To run this document, we need Julia packages
Pkg.add("JLD")
Pkg.add("RCall")
and R packages MixedModels
and lme4
install.packages("lme4")
install.packages("devtools")
devtools::install_github("jrklasen/MixedModels")
The script mendel29_sim.jl
will run the simulations with $r = \sigma_a^2 / \sigma_e^2$ varying in {0, 0.05, 0.1, 1, 10, 20}, and nsim=50
replicates. It saves the results into a file mendel29_sim_result.jld
:
tic()
include("mendel29_sim.jl")
toc()
To obtain Table 3 in manuscript:
include("mendel29_sim_print.jl")