lmm.restricted.likelihood {gaston}

Likelihood of a linear mixed model

Description

Compute the Restricted or the Full Likelihood of a linear mixed model.

Usage

 
 lmm.restricted.likelihood(Y, X = matrix(1, nrow = length(Y)), K, tau, s2)
 lmm.profile.restricted.likelihood(Y, X = matrix(1, nrow = length(Y)), K, h2) 

Arguments

Y

Phenotype vector

X

Covariable matrix

K

A positive definite matrix or a list of such matrices

tau

Value(s) of parameter(s) \tau

s2

Value of parameter \sigma^2

h2

Value(s) of heritability

Details

Theses function respectively compute the Restricted and the Profile Likelihood under the linear mixed model

Y = X\beta + \omega_1 + \ldots + \omega_k + \varepsilon

with \omega_i \sim N(0,\tau_i K_i) for i \in 1, \dots,k and \varepsilon \sim N(0,\sigma^2 I_n) .

The variance matrices K_1, ..., K_k, are specified through the parameter K. The parameter tau should be a vector of length k.

The function lmm.restricted.likelihood computes the restricted likelihood for the given values of \tau and \sigma^2. Whenever k = 1, it is similar to lmm.diago.likelihood(tau, s2, Y = Y, X = X, eigenK = eigen(K)) which should be prefered (with a preliminary computation of eigen(K)).

The function lmm.profile.restricted.likelihood computes a profile restricted likelihood: the values of \tau and \sigma^2 which maximizes the likelihood are computed under the constraint {\tau \over \tau + \sigma^2 } = h^2 , and the profiled likelihood value for these parameters is computed. Whenever k = 1, it is similar to lmm.diago.likelihood(h2 = h2, Y = Y, X = X, eigenK = eigen(K)).

Value

The restricted likelihood value.

Author(s)

Hervé Perdry and Claire Dandine-Roulland

See Also

lmm.diago.likelihood, lmm.diago, lmm.aireml

Examples

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

# Compute Genetic Relationship Matrix and its eigen decomposition
K <- GRM(x)
eiK <- eigen(K)

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

# compute restricted likelihood for tau = 0.2 and s2 = 0.8
lmm.restricted.likelihood(y, K=K, tau = 0.2, s2 = 0.8)

# compute profile restricted likelihood for h2 = 0.2
lmm.profile.restricted.likelihood(y, K=K, h2 = 0.2)

# identity with the values computed with the diagonalisation trick
lmm.diago.likelihood(tau = 0.2, s2 = 0.8, Y = y, eigenK = eiK)
lmm.diago.likelihood(h2 = 0.2, Y = y, eigenK = eiK)