collectResiduals: Residual hooks and collection methods

Description Usage Arguments Value Functions Total residual types Partial residuals See Also Examples

View source: R/zlmHooks.R

Description

After each gene is fit, a hook function can optionally be run and the output saved. This allows extended computations to be done using the fitted model, without keeping it in memory. Here this is used to calculate various residuals, though in some cases they can be done using only the information contained in the ZlmFit-class.

Usage

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13

Arguments

x

ZlmFit-class

sca

SingleCellAssay object to which the residuals should be added

newLayerName

character name of the assay layer

effectRegex

a regular expression naming columns of the design corresponding to Z_0. Generally these should be the treatment effects of interest.

Value

copy of sca with new layer

Functions

Total residual types

Each component of the model contributes several flavors of residual, which can be combined in various fashions. The discrete residual can be on the response scale (thus subtracting the predicted probability of expression from the 0/1 expression value). Or it can be a deviance residual, revealing something about the log-likelihood.

Partial residuals

It's also possible to consider partial residuals, in which the contribution of a particular covariate is added back into the model.

See Also

zlm

Examples

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
data(vbetaFA)
svbeta <- subset(vbetaFA, ncells==1)
svbeta <- svbeta[freq(svbeta)>.4,]
window <- function(x1) lapply(assays(x1), function(x2) x2[1:3, 1:6])
#total residuals of the response
z1 <- zlm(~ Stim.Condition, svbeta, hook=discrete_residuals_hook)
window(collectResiduals(z1, svbeta))
z2 <- zlm(~ Stim.Condition, svbeta, hook=continuous_residuals_hook)
window(collectResiduals(z2, svbeta))
z3 <- zlm(~ Stim.Condition, svbeta, hook=combined_residuals_hook)
window(collectResiduals(z3, svbeta))
#partial residuals
colData(svbeta)$ngeneson <- colMeans(assay(svbeta)>0)
z5 <- zlm(~ Stim.Condition + ngeneson, svbeta)
partialScore(z5, 'Stim.Condition')

MAST documentation built on Nov. 8, 2020, 8:19 p.m.