fit1 {qtl2} | R Documentation |
Fit single-QTL model at a single position
Description
Fit a single-QTL model at a single putative QTL position and get detailed results about estimated coefficients and individuals contributions to the LOD score.
Usage
fit1(
genoprobs,
pheno,
kinship = NULL,
addcovar = NULL,
nullcovar = NULL,
intcovar = NULL,
weights = NULL,
contrasts = NULL,
model = c("normal", "binary"),
zerosum = TRUE,
se = TRUE,
hsq = NULL,
reml = TRUE,
blup = FALSE,
...
)
Arguments
genoprobs |
A matrix of genotype probabilities, individuals x genotypes.
If NULL, we create a single intercept column, matching the individual IDs in |
pheno |
A numeric vector of phenotype values (just one phenotype, not a matrix of them) |
kinship |
Optional kinship matrix. |
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. |
intcovar |
An optional numeric matrix of interactive covariates. |
weights |
An optional numeric vector of positive weights for the
individuals. As with the other inputs, it must have |
contrasts |
An optional numeric matrix of genotype contrasts, size
genotypes x genotypes. For an intercross, you might use
|
model |
Indicates whether to use a normal model (least
squares) or binary model (logistic regression) for the phenotype.
If |
zerosum |
If TRUE, force the genotype or allele coefficients
sum to 0 by subtracting their mean and add another column with
the mean. Ignored if |
se |
If TRUE, calculate the standard errors. |
hsq |
(Optional) residual heritability; used only if
|
reml |
If |
blup |
If TRUE, fit a model with QTL effects being random, as in |
... |
Additional control parameters; see Details; |
Details
For each of the inputs, the row names are used as individual identifiers, to align individuals.
If kinship
is absent, Haley-Knott regression is performed.
If kinship
is provided, a linear mixed model is used, with a
polygenic effect estimated under the null hypothesis of no (major)
QTL, and then taken as fixed as known in the genome scan.
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.
The ...
argument can contain several additional control
parameters; suspended for simplicity (or confusion, depending on
your point of view). tol
is used as a tolerance value for linear
regression by QR decomposition (in determining whether columns are
linearly dependent on others and should be omitted); default
1e-12
. maxit
is the maximum number of iterations for
converence of the iterative algorithm used when model=binary
.
bintol
is used as a tolerance for converence for the iterative
algorithm used when model=binary
. eta_max
is the maximum value
for the "linear predictor" in the case model="binary"
(a bit of a
technicality to avoid fitted values exactly at 0 or 1).
Value
A list containing
-
coef
- Vector of estimated coefficients. -
SE
- Vector of estimated standard errors (included ifse=TRUE
). -
lod
- The overall lod score. -
ind_lod
- Vector of individual contributions to the LOD score (not provided ifkinship
is used). -
fitted
- Fitted values. -
resid
- Residuals. Ifblup==TRUE
, onlycoef
andSE
are included at present.
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.
See Also
pull_genoprobpos()
, find_marker()
Examples
# read data
iron <- read_cross2(system.file("extdata", "iron.zip", package="qtl2"))
# insert pseudomarkers into map
map <- insert_pseudomarkers(iron$gmap, step=5)
# calculate genotype probabilities
probs <- calc_genoprob(iron, map, error_prob=0.002)
# 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)
# scan chromosome 7 to find peak
out <- scan1(probs[,"7"], pheno, addcovar=covar)
# find peak position
max_pos <- max(out, map)
# genoprobs at max position
pr_max <- pull_genoprobpos(probs, map, max_pos$chr, max_pos$pos)
# fit QTL model just at that position
out_fit1 <- fit1(pr_max, pheno, addcovar=covar)