lmebreed {lme4breeding} | R Documentation |
Fit linear or generalized linear mixed models incorporating the effects of relationships.
lmebreed(formula, data, family = NULL, REML = TRUE,
relmat = list(), addmat = list(), control = list(),
start = NULL, verbose = FALSE, subset, weights,
na.action, offset, contrasts = NULL, model = TRUE,
x = TRUE, dateWarning=TRUE, returnParams=FALSE,
rotation=FALSE, ...)
relmat |
an optional named list of relationship matrices (not the inverse).
Internally the Cholesky decomposition of those matrices will be computed.
The names of the elements must correspond to the names of grouping factors for
random-effects terms in the |
addmat |
an optional named list of customized incidence matrices.
The names of the elements must correspond to the names of grouping factors for
random-effects terms in the |
dateWarning |
an logical value indicating if you want to be warned when a new version of lme4breeding is available on CRAN. |
returnParams |
an logical value indicating if you want to only get the incidence matrices of the model. |
rotation |
an logical value indicating if you want to compute the eigen decomposition of the relationship matrix to rotate y and X and accelerate the computation. See details. |
formula |
as in |
data |
as in |
family |
as in |
REML |
as in |
control |
as in |
start |
as in |
verbose |
as in |
subset |
as in |
weights |
as in |
na.action |
as in |
offset |
as in |
contrasts |
as in |
model |
as in |
x |
as in |
... |
as in |
All arguments to this function are the same as those to the function
lmer
except relmat
and addmat
which must be
named lists. Each name must correspond to the name of a grouping factor in a
random-effects term in the formula
. The observed levels
of that factor must be contained in the rownames and columnames of the relmat.
Each relmat is the relationship matrix restricted
to the observed levels and applied to the model matrix for that term. The incidence
matrices in the addmat argument must match the dimensions of the final fit (pay
attention to missing data in responses).
It is important to remember that when you use the relmat
argument you are providing
the square root of the relationship matrix and to recover the correct BLUPs for those effects
you need to use the ranef
function which internally multiple those BLUPs the
square root of the relationship matrix one more time to recover the correct BLUPs.
The argument rotation
applies the eigen decomposition proposed by Lee and Van der Werf in 2016
and makes the genetic evaluation totally sparse leading to incredible gains in speed compared
to the classical approach. Internally, the eigen decomposition UDU' is carried in the relationship
matrix. The U matrix is then taken to the n x n level (n being the number of records), and post-multiplied
by a matrix of records presence (n x n) using the element-wise multiplication of two matrices (Schur product).
Example Datasets
The package has been equiped with several datasets to learn how to use the lme4breeding package:
* DT_halfdiallel
, DT_fulldiallel
and DT_mohring
datasets have examples to fit half and full diallel designs.
* DT_h2
to calculate heritability
* DT_cornhybrids
and DT_technow
datasets to perform genomic prediction in hybrid single crosses
* DT_wheat
dataset to do genomic prediction in single crosses in species displaying only additive effects.
* DT_cpdata
dataset to fit genomic prediction models within a biparental population coming from 2 highly heterozygous parents including additive, dominance and epistatic effects.
* DT_polyploid
to fit genomic prediction and GWAS analysis in polyploids.
* DT_gryphon
data contains an example of an animal model including pedigree information.
* DT_btdata
dataset contains an animal (birds) model.
* DT_legendre
simulated dataset for random regression model.
* DT_sleepstudy
dataset to know how to translate lme4 models to sommer models.
* DT_ige
dataset to show how to fit indirect genetic effect models.
a lmebreed
object.
Giovanny Covarrubias-Pazaran
Giovanny Covarrubias-Pazaran (2024). lme4breeding: enabling genetic evaluation in the age of genomic data. To be submitted to Bioinformatics.
Douglas Bates, Martin Maechler, Ben Bolker, Steve Walker (2015). Fitting Linear Mixed-Effects Models Using lme4. Journal of Statistical Software, 67(1), 1-48.
Lee & Van der Werf (2016). MTG2: an efficient algorithm for multivariate linear mixed model analysis based on genomic information. Bioinformatics, 32(9), 1420-1422.
data(DT_example)
DT <- DT_example
A <- A_example
ansMain <- lmebreed(Yield ~ Env + (1|Name),
relmat = list(Name = A ),
data=DT)
vc <- VarCorr(ansMain); print(vc,comp=c("Variance"))