coverage.density: Coverage density plot

Description Usage Arguments Details Value Author(s) See Also Examples

View source: R/coverage.density.R

Description

Visualization of target coverage density for one or more samples.

Usage

1
coverage.density(coveragelist, normalized = TRUE, legend, main, xlab, col, lwd, lty, xlim, ylim, ...)

Arguments

coveragelist

Output of function coverage.target, where option perBase had to be set to TRUE, i.e. a list with elements coverageTarget and avgTargetCoverage. Or, when density of several samples shall be visualized, a list with respective outputs of coverage.target.

normalized

if TRUE, densities of normalized coverages will be shown; original coverages otherwise

legend

legend text. If missing, names of coveragelist will be taken. If NULL, no legend will be drawn.

main

main title

xlab

x-axis label

col

line color(s)

lwd

line width(s)

lty

line style(s)

xlim, ylim

x- and y-axis coordinate ranges

...

further graphical parameters passed to plot

Details

If normalized = TRUE, the function calculates normalized coverages: per-base coverages divided by average coverage over all targeted bases. Normalized coverages are not dependent on the absolute quantity of reads and are hence better comparable between different samples or even different experiments.

Value

Line plot(s) showing densities.

Author(s)

Manuela Hummel m.hummel@dkfz.de

See Also

coverage.target, covered.k, coverage.hist, coverage.uniformity, coverage.correlation, coverage.plot

Examples

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
## get reads and targets
exptPath <- system.file("extdata", package="TEQC")
readsfile <- file.path(exptPath, "ExampleSet_Reads.bed")
reads <- get.reads(readsfile, idcol=4, skip=0)
targetsfile <- file.path(exptPath, "ExampleSet_Targets.bed")
targets <- get.targets(targetsfile, skip=0)

## calculate per-base coverages
Coverage <- coverage.target(reads, targets, perBase=TRUE)

## coverage density
coverage.density(Coverage)

Example output

Loading required package: BiocGenerics
Loading required package: parallel

Attaching package: 'BiocGenerics'

The following objects are masked from 'package:parallel':

    clusterApply, clusterApplyLB, clusterCall, clusterEvalQ,
    clusterExport, clusterMap, parApply, parCapply, parLapply,
    parLapplyLB, parRapply, parSapply, parSapplyLB

The following objects are masked from 'package:stats':

    IQR, mad, sd, var, xtabs

The following objects are masked from 'package:base':

    Filter, Find, Map, Position, Reduce, anyDuplicated, append,
    as.data.frame, cbind, colMeans, colSums, colnames, do.call,
    duplicated, eval, evalq, get, grep, grepl, intersect, is.unsorted,
    lapply, lengths, mapply, match, mget, order, paste, pmax, pmax.int,
    pmin, pmin.int, rank, rbind, rowMeans, rowSums, rownames, sapply,
    setdiff, sort, table, tapply, union, unique, unsplit, which,
    which.max, which.min

Loading required package: IRanges
Loading required package: S4Vectors
Loading required package: stats4

Attaching package: 'S4Vectors'

The following object is masked from 'package:base':

    expand.grid

Loading required package: Rsamtools
Loading required package: GenomeInfoDb
Loading required package: GenomicRanges
Loading required package: Biostrings
Loading required package: XVector

Attaching package: 'Biostrings'

The following object is masked from 'package:base':

    strsplit

Loading required package: hwriter
[1] "read 19546 sequenced reads"
[1] "read 50 (non-overlapping) target regions"
Warning message:
system call failed: Cannot allocate memory 

TEQC documentation built on Nov. 8, 2020, 8:07 p.m.