callPatterns: Extract posterior probabilities (or combinatorial patterns)...

View source: R/callPaterns.R

callPatternsR Documentation

Extract posterior probabilities (or combinatorial patterns) associated with differential regions

Description

Given results from epigraHMM's differential peak caller, this function will output either posterior probabilities or combinatorial patterns associated with the mixture components of the embedded mixture model.

Usage

callPatterns(
  object,
  peaks,
  hdf5 = metadata(object)$output,
  type = "all",
  fdr = NULL,
  pattern = NULL,
  ranges = NULL
)

Arguments

object

an epigraHMMDataSet

peaks

a GRanges object with differential peaks from 'callPeaks'

hdf5

a character with the location of the epigraHMM HDF5 output file

type

a character string that defines which output will be givem (see details; default is 'all')

fdr

the desired fdr thresholding level to define combinatorial patterns

pattern

a string that explicitly specifies the combinatorial pattern to be output

ranges

a GRanges object with the genomic ranges to subset the output

Details

The output of 'callPatterns' is always restricted to genomic windows intersecting peaks.

If ‘type = ’all'‘, all windows’ posterior probabilities associated with all differential combinatorial patterns are returned. If ‘type = ’fdr'', users must also specify the input argument 'pattern' and this function will output windows wich are associated with the given 'pattern' that pass a particular fdr threshold level. If ‘type = ’max'', this function will output the combinatorial pattern which has the maximal posterior probability for each window. If ‘type = ’ranges'', the windows that are output are restricted to those that intersect the 'ranges' input argument.

Value

A GRanges object with metadata

Author(s)

Pedro L. Baldoni, pedrobaldoni@gmail.com

References

https://github.com/plbaldoni/epigraHMM

Examples

# Creating dummy object
countData <- cbind(rbind(matrix(rnbinom(1e2, mu = 1, size = 10), ncol = 1),
                         matrix(rnbinom(1e2, mu = 10, size = 5), ncol = 1),
                         matrix(rnbinom(1e2, mu = 1, size = 10), ncol = 1),
                         matrix(rnbinom(1e2, mu = 10, size = 5), ncol = 1),
                         matrix(rnbinom(1e2, mu = 1, size = 10), ncol = 1),
                         matrix(rnbinom(1e2, mu = 1, size = 10), ncol = 1),
                         matrix(rnbinom(1e2, mu = 1, size = 10), ncol = 1)),
                   rbind(matrix(rnbinom(1e2, mu = 1, size = 10), ncol = 1),
                         matrix(rnbinom(1e2, mu = 1, size = 10), ncol = 1),
                         matrix(rnbinom(1e2, mu = 1, size = 10), ncol = 1),
                         matrix(rnbinom(1e2, mu = 10, size = 5), ncol = 1),
                         matrix(rnbinom(1e2, mu = 1, size = 10), ncol = 1),
                         matrix(rnbinom(1e2, mu = 10, size = 5), ncol = 1),
                         matrix(rnbinom(1e2, mu = 1, size = 10), ncol = 1)))

colData <- data.frame(condition = c('A','B'), replicate = c(1,1))
rowRanges <- GenomicRanges::GRanges('chrA',
                     IRanges::IRanges(start = seq(1,by = 500,
                     length.out = nrow(countData)),width = 500))

object <- epigraHMMDataSetFromMatrix(countData,colData,rowRanges = rowRanges)

# Initializing
object <- initializer(object,controlEM())

# Running epigraHMM
object <- epigraHMM(object,controlEM(),type = 'differential',dist = 'nb')

# Calling peaks
peaks <- callPeaks(object = object,
                  hdf5 = S4Vectors::metadata(object)$output,
                  method = 'viterbi')
                  
# Extracting posterior probabilities
patterns <- callPatterns(object = object,peaks = peaks,type = 'max')


plbaldoni/epigraHMM documentation built on Oct. 15, 2023, 7:53 p.m.