scan_genomic_contigs {spiky} | R Documentation |
The default workflow for spiky is roughly as follows:
scan_genomic_contigs(bam, spike, param = NULL, ...)
bam |
the BAM or CRAM file |
spike |
the spike-in reference database (e.g. data(spike)) |
param |
a ScanBamParam object specifying which reads to count (NULL) |
... |
additional arguments to pass to scanBamFlag() |
Identify and quantify the spike-in contigs in an experiment.
Fit a model for sequence-based abundance artifacts using the spike-ins.
Quantify raw fragment abundance on genomic contigs, and adjust per step 2.
scan_genomic_contigs addresses the first half of step 3. The assumption is
that anything which isn't a spike contig, is a genomic contig. This isn't
necessarily true, so the user can also supply a ScanBamParam object for the
param
argument and restrict scanning to whatever contigs they wish, which
also allows for non-default MAPQ, pairing, and quality filters.
a CompressedGRangesList with bin- and spike-level coverage
Rsamtools::ScanBamParam
library(Rsamtools) data(spike, package="spiky") fl <- system.file("extdata", "ex1.bam", package="Rsamtools", mustWork=TRUE) scan_genomic_contigs(fl, spike=spike) # will warn user about spike contigs sb <- system.file("extdata", "example.spike.bam", package="spiky", mustWork=TRUE) scan_genomic_contigs(sb, spike=spike) # will warn user about genomic contigs