lmm.diago: Linear mixed model fitting with the diagonalization trick

View source: R/lmm_diago.r

lmm.diagoR Documentation

Linear mixed model fitting with the diagonalization trick

Description

Estimate the parameters of a linear mixed model, using the "diagonalization trick".

Usage

 lmm.diago(Y, X = matrix(1, nrow=length(Y)), eigenK, p = 0, 
                  method = c("newton", "brent"), min_h2 = 0, max_h2 = 1,
                  verbose = getOption("gaston.verbose", TRUE), 
                  tol = .Machine$double.eps^0.25) 

Arguments

Y

Phenotype vector

X

Covariable matrix

eigenK

Eigen decomposition of K (a positive symmetric matrix)

p

Number of Principal Components included in the mixed model with fixed effect

method

Optimization method to use

min_h2

Minimum admissible value

max_h2

Maximum admissible value

verbose

If TRUE, display information on the function actions

tol

Accuracy of estimation

Details

Estimate the parameters of the following linear mixed model, computing the restricted likelihood as in lmm.diago.likelihood, and using either a Newton algorithm, or Brent algorithm as in optimize:

Y = (X|PC)\beta + \omega + \varepsilon

with \omega \sim N(0,\tau K) and \varepsilon \sim N(0,\sigma^2 I_n) .

The matrix K is given through its eigen decomposition, as produced by eigenK = eigen(K, symmetric = TRUE). The matrix (X|PC) is the concatenation of the covariable matrix X and of the first p eigenvectors of K, included in the model with fixed effects.

Value

If the parameter p is a scalar, a list with following elements :

sigma2

Estimate of the model parameter \sigma^2

tau

Estimate(s) of the model parameter(s) \tau_1, \dots, \tau_k

Py

Last computed value of vector Py (see reference)

BLUP_omega

BLUPs of random effects

BLUP_beta

BLUPs of fixed effects \beta (only the components corresponding to X)

Xbeta

Estimate of (X|PC)\beta

varbeta

Variance matrix for \beta estimates (only the components corresponding to X)

varXbeta

Participation of fixed effects to variance of Y

p

Number of Principal Components included in the linear mixed model with fixed effect

If the paramer p is a vector of length > 1, a list of lists as described above, one for each value in p.

Author(s)

Hervé Perdry and Claire Dandine-Roulland

See Also

lmm.diago.likelihood, lmm.aireml, optimize

Examples

# Load data
data(AGT)
x <- as.bed.matrix(AGT.gen, AGT.fam, AGT.bim)

# Compute Genetic Relationship Matrix
K <- GRM(x)

# eigen decomposition of K
eiK <- eigen(K)

# simulate a phenotype
set.seed(1)
y <- 1 + lmm.simu(tau = 1, sigma2 = 2, eigenK = eiK)$y

# Estimations
R <- lmm.diago(Y = y, eigenK = eiK, p = c(0,10))
str(R)

gaston documentation built on May 29, 2024, 7:33 a.m.