compHists {MMDiff2} | R Documentation |
This function computes histograms at pre-defined regions (peaks)
from mapped fragments, i.e. fragment counts at genomic position. Note,
in contrast to genomic coverage or density maps, this function uses a single
position per fragment (usually its center) rather than the whole extend of
the fragment. This results in a significant increase in resolution.
The parameter whichPos
determines whether
fragment centers, start or end positions should be considered
('Center','Left','Right').
Results are stored as a list in the Hists
slot of the DBAmmd Object,
with one entry per peak. For each peak i, a (n x L_i) matrix is generated,
where n is the number of samples and L_i is the number of bins used to cover
the extend of the peak. Note, L_i varies between peaks of different lengths.
compHists(MD, bin.length = 20, whichPos = "Center")
MD |
DBAmmd Object. This Object can be created using |
bin.length |
size of binning window (in bp) (DEFAULT: 20) |
whichPos |
specifies which relative positions of mapped fragments should to be considered. Can be one of: 'Left.p', 'Right.p', 'Right.p' and 'Left.n': Start and end positions of fragments mapping to positive or negative strand, respectively ('Right.p' and 'Left.n' are not available for single-end reads). Additionally inferred positions: 'Center.n','Center.p','Center','Left','Right'. (DEFAULT: 'Center') |
DBAmmd object with updated slot Hists
DBAmmd
, getPeakReads
,
estimateFragmentCenters
, plotPeak
,
## Example using a small data set provided with this package: data("MMD") bin.length <- 20 MMD.1 <- compHists(MMD,bin.length) # use \code{plotPeak()} to plot indivdual peaks: Peak.id <- '6' plotPeak(MMD.1, Peak.id=Peak.id) # or explicitly using the histograms: H <- Hists(MMD.1, whichPos='Center') Sample <- 'WT.AB2' Peak.idx <- match(Peak.id, names(Regions(MMD.1))) plot(H[[Peak.idx]][Sample,],t='l') # add peak cooridnates: Peak <- Regions(MMD.1)[Peak.idx] meta <- metaData(MMD.1) PeakBoundary <- meta$AnaData$PeakBoundary x.coords <- as.integer(colnames(H[[Peak.idx]])) + start(Peak) - PeakBoundary plot(x.coords,H[[Peak.idx]]['WT.AB2',],t='l', xlab=names(H)[Peak.idx], ylab='counts', main=Sample)