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 β_C + ∑_{k=1}^K Z_k u_k
with X_C the matrix provided in Cofactor
, β_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 σ_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 β_C + X_{[i]} β_{[i]} + ∑_{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 β_{[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 β_C and β_{i} |
Sigma2 |
Estimated values of σ_1^2,...,σ_K^2 |
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.