qtl2 (version 0.32)

scan1blup: Calculate BLUPs of QTL effects in scan along one chromosome

Description

Calculate BLUPs of QTL effects in scan along one chromosome, with a single-QTL model treating the QTL effects as random, with possible allowance for covariates and for a residual polygenic effect.

Usage

scan1blup(
  genoprobs,
  pheno,
  kinship = NULL,
  addcovar = NULL,
  nullcovar = NULL,
  contrasts = NULL,
  se = FALSE,
  reml = TRUE,
  tol = 0.000000000001,
  cores = 1,
  quiet = TRUE
)

Value

An object of class "scan1coef": a matrix of estimated regression coefficients, of dimension positions x number of effects. The number of effects is n_genotypes + n_addcovar + (n_genotypes-1)*n_intcovar. May also contain the following attributes:

• SE - Present if se=TRUE: a matrix of estimated standard errors, of same dimension as coef.

• sample_size - Vector of sample sizes used for each phenotype

Arguments

genoprobs

Genotype probabilities as calculated by calc_genoprob().

pheno

A numeric vector of phenotype values (just one phenotype, not a matrix of them)

kinship

Optional kinship matrix, or a list of kinship matrices (one per chromosome), in order to use the LOCO (leave one chromosome out) method.

addcovar

An optional numeric matrix of additive covariates.

nullcovar

An optional numeric matrix of additional additive covariates that are used under the null hypothesis (of no QTL) but not under the alternative (with a QTL). This is needed for the X chromosome, where we might need sex as a additive covariate under the null hypothesis, but we wouldn't want to include it under the alternative as it would be collinear with the QTL effects. Only used if kinship is provided but hsq is not, to get estimate of residual heritability.

contrasts

An optional numeric matrix of genotype contrasts, size genotypes x genotypes. For an intercross, you might use cbind(mu=c(1,0,0), a=c(-1, 0, 1), d=c(0, 1, 0)) to get mean, additive effect, and dominance effect. The default is the identity matrix.

se

If TRUE, also calculate the standard errors.

reml

If reml=TRUE, use REML to estimate variance components; otherwise maximum likelihood.

tol

Tolerance value for convergence of linear mixed model fit.

cores

Number of CPU cores to use, for parallel calculations. (If 0, use parallel::detectCores().) Alternatively, this can be links to a set of cluster sockets, as produced by parallel::makeCluster().

quiet

If FALSE, print message about number of cores used when multi-core.

Details

For each of the inputs, the row names are used as individual identifiers, to align individuals.

If kinship is provided, the linear mixed model accounts for a residual polygenic effect, with a the polygenic variance estimated under the null hypothesis of no (major) QTL, and then taken as fixed as known in the scan to estimate QTL effects.

If contrasts is provided, the genotype probability matrix, \(P\), is post-multiplied by the contrasts matrix, \(A\), prior to fitting the model. So we use \(P \cdot A\) as the \(X\) matrix in the model. One might view the rows of A^-1 as the set of contrasts, as the estimated effects are the estimated genotype effects pre-multiplied by A^-1.

References

Haley CS, Knott SA (1992) A simple regression method for mapping quantitative trait loci in line crosses using flanking markers. Heredity 69:315--324.

Kang HM, Zaitlen NA, Wade CM, Kirby A, Heckerman D, Daly MJ, Eskin E (2008) Efficient control of population structure in model organism association mapping. Genetics 178:1709--1723.

Robinson GK (1991) That BLUP is a good thing: The estimation of random effects. Statist Sci 6:15--32.

Examples

# read data
iron <- read_cross2(system.file("extdata", "iron.zip", package="qtl2"))
iron <- iron[,c(7,19)] # reduce to chr 7 and 19

# insert pseudomarkers into map
map <- insert_pseudomarkers(iron$gmap, step=1)

# calculate genotype probabilities
probs <- calc_genoprob(iron, map, error_prob=0.002)
probs[["7"]] <- probs[["7"]][,,1:5] # reduce to very small number

# convert to allele probabilities
aprobs <- genoprob_to_alleleprob(probs)

# grab phenotypes and covariates; ensure that covariates have names attribute
pheno <- iron$pheno[,1]
covar <- match(iron$covar$sex, c("f", "m")) # make numeric
names(covar) <- rownames(iron$covar)

# calculate BLUPs of coefficients for chromosome 7
blup <- scan1blup(aprobs[,"7"], pheno, addcovar=covar)

# leave-one-chromosome-out kinship matrix for chr 7
kinship7 <- calc_kinship(probs, "loco")[["7"]]

# calculate BLUPs of coefficients for chromosome 7, adjusting for residual polygenic effect
blup_pg <- scan1blup(aprobs[,"7"], pheno, kinship7, addcovar=covar)