calculateStats: Calculate F-statistics at base pair resolution from a loaded...

In lcolladotor/derfinder: Annotation-agnostic differential expression analysis of RNA-seq data at base-pair resolution via the DER Finder approach

Calculate F-statistics at base pair resolution from a loaded BAM files

Description

After defining the models of interest (see makeModels) and pre-processing the data (see preprocessCoverage), use calculateStats to calculate the F-statistics at base-pair resolution.

Usage

calculateStats(coveragePrep, models, lowMemDir = NULL, ...)

Arguments

coveragePrep	A list with ⁠$coverageProcessed⁠, ⁠$mclapplyIndex⁠, and ⁠$position⁠ normally generated using preprocessCoverage.
models	A list with ⁠$mod⁠ and ⁠$mod0⁠ normally generated using makeModels.
lowMemDir	The directory where the processed chunks are saved when using preprocessCoverage with a specified lowMemDir.
...	Arguments passed to other methods and/or advanced arguments. Advanced arguments: verbose If TRUE basic status updates will be printed along the way. scalefac This argument is passed to fstats.apply and should be the same as the one used in preprocessCoverage. Default: 32. method Has to be either 'Matrix' (default), 'Rle' or 'regular'. See details in fstats.apply. adjustF A single value to adjust that is added in the denominator of the F-stat calculation. Useful when the Residual Sum of Squares of the alternative model is very small. Default: 0. Passed to define_cluster.

Value

A numeric Rle with the F-statistics per base pair that passed the cutoff.

Author(s)

Leonardo Collado-Torres

See Also

makeModels, preprocessCoverage

Examples

## Collapse the coverage information
collapsedFull <- collapseFullCoverage(list(genomeData$coverage),
    verbose = TRUE
)

## Calculate library size adjustments
sampleDepths <- sampleDepth(collapsedFull, probs = c(0.5), verbose = TRUE)

## Build the models
group <- genomeInfo$pop
adjustvars <- data.frame(genomeInfo$gender)
models <- makeModels(sampleDepths, testvars = group, adjustvars = adjustvars)

## Preprocess the data
prep <- preprocessCoverage(genomeData,
    cutoff = 0, scalefac = 32,
    chunksize = 1e3, colsubset = NULL
)

## Run the function
fstats <- calculateStats(prep, models, verbose = TRUE, method = "regular")
fstats
## Not run: 
## Compare vs pre-packaged F-statistics
library("testthat")
expect_that(fstats, is_equivalent_to(genomeFstats))

## End(Not run)

lcolladotor/derfinder documentation built on Dec. 17, 2024, 4:53 p.m.