EMM2.cpp {RAINBOWR}

Equation of mixed model for one kernel, EMMA-based method (inplemented by Rcpp)

Description

This function solves single-kernel linear mixed model by EMMA (efficient mixed model association; Kang et al., 2008) approach.

Usage

EMM2.cpp(
  y,
  X = NULL,
  ZETA,
  eigen.G = NULL,
  eigen.SGS = NULL,
  tol = NULL,
  optimizer = "nlminb",
  traceInside = 0,
  REML = TRUE,
  bounds = c(1e-09, 1e+09),
  SE = FALSE,
  return.Hinv = FALSE
)

Arguments

y

A n \times 1 vector. A vector of phenotypic values should be used. NA is allowed.

X

A n \times p matrix. You should assign mean vector (rep(1, n)) and covariates. NA is not allowed.

ZETA

A list of variance (relationship) matrix (K; m \times m) and its design matrix (Z; n \times m) of random effects. You can use only one kernel matrix. For example, ZETA = list(A = list(Z = Z, K = K)) Please set names of list "Z" and "K"!

eigen.G

A list with

$values

Eigen values

$vectors

Eigen vectors

The result of the eigen decompsition of G = ZKZ'. You can use "spectralG.cpp" function in RAINBOWR. If this argument is NULL, the eigen decomposition will be performed in this function. We recommend you assign the result of the eigen decomposition beforehand for time saving.

eigen.SGS

A list with

$values

Eigen values

$vectors

Eigen vectors

The result of the eigen decompsition of SGS, where S = I - X(X'X)^{-1}X', G = ZKZ'. You can use "spectralG.cpp" function in RAINBOWR. If this argument is NULL, the eigen decomposition will be performed in this function. We recommend you assign the result of the eigen decomposition beforehand for time saving.

tol

The tolerance for detecting linear dependencies in the columns of G = ZKZ'. Eigen vectors whose eigen values are less than "tol" argument will be omitted from results. If tol is NULL, top 'n' eigen values will be effective.

optimizer

The function used in the optimization process. We offer "optim", "optimx", and "nlminb" functions.

traceInside

Perform trace for the optimzation if traceInside >= 1, and this argument shows the frequency of reports.

REML

You can choose which method you will use, "REML" or "ML". If REML = TRUE, you will perform "REML", and if REML = FALSE, you will perform "ML".

bounds

Lower and Upper bounds of the parameter lambda. If the updated parameter goes out of this range, the parameter is reset to the value in this range.

SE

If TRUE, standard errors are calculated.

return.Hinv

If TRUE, the function returns the inverse of H = ZKZ' + \lambda I where \lambda = \sigma^2_e / \sigma^2_u. This is useful for GWAS.

Value

$Vu

Estimator for \sigma^2_u

$Ve

Estimator for \sigma^2_e

$beta

BLUE(\beta)

$u

BLUP(u)

$LL

Maximized log-likelihood (full or restricted, depending on method)

$beta.SE

Standard error for \beta (If SE = TRUE)

$u.SE

Standard error for u^*-u (If SE = TRUE)

$Hinv

The inverse of H = ZKZ' + \lambda I (If return.Hinv = TRUE)

References

Kang, H.M. et al. (2008) Efficient Control of Population Structure in Model Organism Association Mapping. Genetics. 178(3): 1709-1723.