MMEst | R Documentation |
This is the main function of the MM4LMM package. It performs inference in a variance component mixed model using a Min-Max algorithm. Inference in multiple models (e.g. for GWAS analysis) can also be performed.
MMEst(Y, Cofactor = NULL, X = NULL, formula=NULL, VarList, ZList = NULL, Method = "Reml",
Henderson=NULL, Init = NULL, CritVar = 0.001, CritLogLik = 0.001,
MaxIter = 100, NbCores = 1, Verbose = TRUE)
Y |
A vector of response values. |
Cofactor |
An incidence matrix corresponding to fixed effects common to all models to be adjusted. If |
X |
An incidence matrix or a list of incidence matrices corresponding to fixed effects specific to each model. If |
formula |
A formula object specifying the fixed effect part of all models separated by + operators. To specify an interaction between |
VarList |
A list of covariance matrices associated with random and residual effects. |
ZList |
A list of incidence matrices associated with random and residual effects (default is |
Method |
The method used for inference. Available methods are "Reml" (Restricted Maximum Likelihood) and "ML" (Maximum Likelihood). |
Henderson |
If |
Init |
A vector of initial values for variance parameters (default is |
CritVar |
Value of the criterion for the variance components to stop iteration. (see Details) |
CritLogLik |
Value of the criterion for the log-likelihood to stop iteration. (see Details) |
MaxIter |
Maximum number of iterations per model. |
NbCores |
Number of cores to be used. |
Verbose |
A boolean describing if messages have to be printed (TRUE) or not (FALSE). Default is TRUE. |
If X
is NULL
, the following model is fitted:
Y = X_C \beta_C + \sum_{k=1}^K Z_k u_k
with X_C
the matrix provided in Cofactor
, \beta_C
the unknown fixed effects, Z_k
the incidence matrix provided for the kth component of ZList
and u_k
the kth vector of random effects. If ZList
is unspecified, all incidence matrices are assumed to be the Identity matrix. Random effects are assumed to follow a Gaussian distribution with mean 0 and covariance matrix R_k \sigma_k^2
, where R_k
is the kth correlation matrix provided in VarList
.
If X
is not NULL
, the following model is fitted for each i:
Y = X_C \beta_C + X_{[i]} \beta_{[i]} + \sum_{k=1}^K Z_k u_k
where X_{[i]}
is the incidence matrix corresponding to the ith component (i.e. column if X
is a matrix, element otherwise) of X
, and \beta_{[i]}
is the (unknow) fixed effect associated to X_{[i]}
.
All models are fitted using the MM algorithm. If Henderson
=TRUE
, at each step the quantities required for updating the variance components are computed using the Mixed Model Equation (MME) trick. See Johnson et al. (1995) for details.
The result is a list where each element corresponds to a fitted model. Each element displays the following:
Beta |
Estimated values of |
Sigma2 |
Estimated values of |
VarBeta |
Variance matrix of |
LogLik (Method) |
The value of the (restricted, if |
NbIt |
The number of iterations required to reach the optimum |
Method |
The method used for the inference |
attr |
An integer vector with an entry for each element of |
Factors |
Names of each term in the formula |
F. Laporte and T. Mary-Huard
Laporte, F., Charcosset, A., & Mary-Huard, T. (2022). Efficient ReML inference in variance component mixed models using a Min-Max algorithm. PLOS Computational Biology, 18(1), e1009659.
Johnson, D. L., & Thompson, R. (1995). Restricted maximum likelihood estimation of variance components for univariate animal models using sparse matrix techniques and average information. Journal of dairy science, 78(2), 449-456.
Hunter, D. R., & Lange, K. (2004). A tutorial on MM algorithms. The American Statistician, 58(1), 30-37.
Zhou, H., Hu, L., Zhou, J., & Lange, K. (2015). MM algorithms for variance components models. arXiv preprint arXiv:1509.07426.
require('MM4LMM')
#### Example 1: variance component analysis, 1 model
data(VarianceComponentExample)
DataHybrid <- VarianceComponentExample$Data
KinF <- VarianceComponentExample$KinshipF
KinD <- VarianceComponentExample$KinshipD
##Build incidence matrix for each random effect
Zf <- t(sapply(as.character(DataHybrid$CodeFlint), function(x)
as.numeric(rownames(KinF)==x)))
Zd <- t(sapply(as.character(DataHybrid$CodeDent), function(x)
as.numeric(rownames(KinD)==x)))
##Build the VarList and ZList objects
VL = list(Flint=KinF , Dent=KinD , Error = diag(1,nrow(DataHybrid)))
ZL <- list(Flint=Zf , Dent=Zd , Error = diag(1,nrow(DataHybrid)))
##Perform inference
#A first way to call MMEst
ResultVA <- MMEst(Y=DataHybrid$Trait , Cofactor = matrix(DataHybrid$Trial)
, ZList = ZL , VarList = VL)
length(ResultVA)
print(ResultVA)
#A second way to call MMEst (same result)
Formula <- as.formula('~ Trial')
ResultVA2 <- MMEst(Y=DataHybrid$Trait , Cofactor = DataHybrid,
formula = Formula
, ZList = ZL , VarList = VL)
length(ResultVA2)
print(ResultVA2)
#### Example 2: Marker Selection with interaction between Cofactor and X matrix
Formula <- as.formula('~ Trial+Xeffect+Xeffect:Trial')
ResultVA3 <- MMEst(Y=DataHybrid$Trait , Cofactor = DataHybrid,
X = VarianceComponentExample$Markers,
formula = Formula
, ZList = ZL , VarList = VL)
length(ResultVA3)
print(ResultVA3[[1]])
#### Example 3: QTL detection with two variance components
data(QTLDetectionExample)
Pheno <- QTLDetectionExample$Phenotype
Geno <- QTLDetectionExample$Genotype
Kinship <- QTLDetectionExample$Kinship
##Build the VarList object
VLgd <- list(Additive=Kinship , Error=diag(1,length(Pheno)))
##Perform inference
ResultGD <- MMEst(Y=Pheno , X=Geno
, VarList=VLgd , CritVar = 10e-5)
length(ResultGD)
print(ResultGD[[1]])
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.