Hua-Zhou/MultiResponseVarianceComponentModels.jl

Multivariate Response Variance Component Models

MultiResponseVarianceComponentModels

MultiResponseVarianceComponentModels.jl is a Julia package that allows fitting and testing multivariate response variance components linear mixed models of form

Installation

julia> ]
pkg> add https://github.com/Hua-Zhou/MultiResponseVarianceComponentModels.jl.git

Documentation

Examples

using MultiResponseVarianceComponentModels, LinearAlgebra, Random
# simulation
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)
# maximum likelihood estimation
model = MRVCModel(Y, X, V)
@timev fit!(model, verbose = true) # ~ 30 seconds
# residual maximum likelihood estimation
@timev fit!(model, reml = true, verbose = true)

model.Σ
reduce(hcat, [hcat(vec(Σ[i]), vec(model.Σ[i])) for i in 1:m])
model.Σcov # sampling variance by inverse of Fisher information matrix
diag(model.Σcov) # (binomial(d, 2) + d) * m variance/covariance parameters
model.logl

References

• 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)

See also

• 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)