| est_herit | R Documentation |
Estimate the heritability of a set of traits via a linear mixed model, with possible allowance for covariates.
est_herit(
pheno,
kinship,
addcovar = NULL,
weights = NULL,
reml = TRUE,
cores = 1,
...
)
pheno |
A numeric matrix of phenotypes, individuals x phenotypes. |
kinship |
A kinship matrix. |
addcovar |
An optional numeric matrix of additive covariates. |
weights |
An optional numeric vector of positive weights for the
individuals. As with the other inputs, it must have |
reml |
If true, use REML; otherwise, use maximimum likelihood. |
cores |
Number of CPU cores to use, for parallel calculations.
(If |
... |
Additional control parameters (see details). |
We fit the model y = X \beta + \epsilon where
\epsilon is multivariate normal with mean 0 and covariance
matrix \sigma^2 [h^2 (2 K) + I] where
K is the kinship matrix and I is the identity matrix.
If weights are provided, the covariance matrix becomes
\sigma^2 [h^2 (2 K) + D] where
D is a diagonal matrix with the reciprocal of the weights.
For each of the inputs, the row names are used as individual identifiers, to align individuals.
If reml=TRUE, restricted maximum likelihood (reml) is used
to estimate the heritability, separately for each phenotype.
Additional control parameters include tol for the tolerance
for convergence, quiet for controlling whether messages will
be display, max_batch for the maximum number of phenotypes
in a batch, and check_boundary for whether the 0 and 1
boundary values for the estimated heritability will be checked
explicitly.
A vector of estimated heritabilities, corresponding to the
columns in pheno. The result has attributes "sample_size",
"log10lik" and "resid_sd".
# read data
iron <- read_cross2(system.file("extdata", "iron.zip", package="qtl2"))
# insert pseudomarkers into map
map <- insert_pseudomarkers(iron$gmap, step=1)
# calculate genotype probabilities
probs <- calc_genoprob(iron, map, error_prob=0.002)
# kinship matrix
kinship <- calc_kinship(probs)
# grab phenotypes and covariates; ensure that covariates have names attribute
pheno <- iron$pheno
covar <- match(iron$covar$sex, c("f", "m")) # make numeric
names(covar) <- rownames(iron$covar)
# perform genome scan
hsq <- est_herit(pheno, kinship, covar)
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.