ovUnion {atena}R Documentation

Pre-defined overlapping mode functions

Description

The following functions control the way in which overlaps between aligned reads and annotated features are resolved when an aligned read overlaps more than one feature on the same locus:

Usage

ovUnion(reads, features, ignoreStrand, minOverlFract)

ovIntersectionStrict(reads, features, ignoreStrand, minOverlFract)

Arguments

reads

A GAlignments, GAlignmentList or a GAlignmentPairs object.

features

A GRanges object with annotated features.

ignoreStrand

(Default FALSE) A logical which defines if the strand should be taken into consideration when computing the overlap between reads and annotated features. When ignoreStrand = FALSE, an aligned read will be considered to be overlapping an annotated feature as long as they have a non-empty intersecting genomic ranges on the same strand, while when ignoreStrand = TRUE the strand will not be considered.

minOverlFract

(Default 0.2) A numeric scalar. minOverlFract is multiplied by the median read length and the resulting value is used to specify the minoverlap argument from findOverlaps from the IRanges package. When no minimum overlap is required, set minOverlFract = 0.

Details

They take the following parameters:

These functions are given to the mode parameter of the qtex() function and are similar to the functions Union() and IntersectionStrict() from the GenomicAlignments package, with the difference that instead of returning counts of reads overlapping annotated features, they return the actual overlaps, because the counting is deferred to other algorithms that follow some specific strategy when a read maps to more than one feature. For this same reason, these functions lack the inter.feature argument found in the corresponding functions from the GenomicAlignments package.

Value

A Hits object; see the Hits-class manual page.

Examples

bamfiles <- list.files(system.file("extdata", package="atena"),
                       pattern="*.bam", full.names=TRUE)
TE_annot <- readRDS(file = system.file("extdata", "Top28TEs.rds", 
                    package="atena"))
tspar <- TelescopeParam(bfl=bamfiles, teFeatures=TE_annot, 
                        singleEnd = TRUE, ignoreStrand=TRUE,
                        minOverlFract=0L)
tsSE <- qtex(tspar, mode=ovIntersectionStrict)


[Package atena version 1.0.5 Index]