MultiResponseVarianceComponentModels.jl is a Julia package that allows fitting and testing multivariate response variance components linear mixed models of form
julia> ]
pkg> add https://github.com/Hua-Zhou/MultiResponseVarianceComponentModels.jl.git
using MultiResponseVarianceComponentModels, LinearAlgebra, Random
# simulation
begin
Random.seed!(1234)
n = 1_000 # n of observations
d = 4 # n of responses
p = 10 # n of covariates
m = 3 # n of variance components
X = rand(n, p)
B = rand(p, d)
V = [zeros(n, n) for _ in 1:m] # kernel matrices
Σ = [zeros(d, d) for _ in 1:m] # variance components
for i in 1:m
Vi = randn(n, n)
copy!(V[i], Vi' * Vi)
Σi = randn(d, d)
copy!(Σ[i], Σi' * Σi)
end
Ω = zeros(n * d, n * d) # overall nd-by-nd covariance matrix Ω
for i = 1:m
Ω += kron(Σ[i], V[i])
end
Ωchol = cholesky(Ω)
Y = X * B + reshape(Ωchol.L * randn(n * d), n, d)
end
# maximum likelihood estimation
model = MRVCModel(Y, X, V)
@timev fit!(model) # ~ 30 seconds
# residual maximum likelihood estimation
model = MRVCModel(Y, X, V; reml = true)
@timev fit!(model) # ~ 30 seconds
# variance components estimates
model.Σ
# comparison of true values and estimates
reduce(hcat, [hcat(vech(Σ[i]), vech(model.Σ[i])) for i in 1:m])
# sampling variance by inverse of Fisher information matrix
model.Σcov
diag(model.Σcov) # m * (binomial(d, 2) + d) parameters
# log-likelihood
model.logl
- H. Zhou, L. Hu, J. Zhou, and K. Lange: MM algorithms for variance components models (2019) (link)
- M. Kim: Gene regulation in the human brain and the biological mechanisms underlying psychiatric disorders (2022) (link)
- J. Kim, J. Shen, A. Wang, D.V. Mehrotra, S. Ko, J.J. Zhou, and H. Zhou: VCSEL: Prioritizing SNP-set by penalized variance component selection (2021) (link)
- L. Hu, W. Lu, J. Zhou, and H. Zhou: MM algorithms for variance component estimation and selection in logistic linear mixed models (2019) (link)
- J.J. Zhou, T. Hu, D. Qiao, M.H. Cho, and H. Zhou: Boosting gene mapping power and efficiency with efficient exact variance component tests of single nucleotide polymorphism sets (2016) (link)