aggregateTagClusters: Aggregate TCs across all samples

aggregateTagClustersR Documentation

Aggregate TCs across all samples

Description

Aggregates tag clusters (TCs) across all CAGE datasets within the CAGEr object to create a referent set of consensus clusters.

Usage

aggregateTagClusters(
  object,
  tpmThreshold = 5,
  excludeSignalBelowThreshold = TRUE,
  qLow = NULL,
  qUp = NULL,
  maxDist = 100,
  useMulticore = FALSE,
  nrCores = NULL
)

## S4 method for signature 'CAGEr'
aggregateTagClusters(
  object,
  tpmThreshold = 5,
  excludeSignalBelowThreshold = TRUE,
  qLow = NULL,
  qUp = NULL,
  maxDist = 100,
  useMulticore = FALSE,
  nrCores = NULL
)

Arguments

object

A CAGEr object

tpmThreshold

Ignore tag clusters with normalized signal ⁠< tpmThreshold⁠ when constructing the consensus clusters.

excludeSignalBelowThreshold

When TRUE the tag clusters with normalized signal ⁠< tpmThreshold⁠ will not contribute to the total CAGE signal of a consensus cluster. When set to FALSE all TCs that overlap consensus clusters will contribute to the total signal, regardless whether they pass the threshold for constructing the clusters or not.

qLow, qUp

Set which "lower" (or "upper") quantile should be used as 5' (or 3') boundary of the tag cluster. If NULL the start (for qLow) or end (for qUp) position of the TC is used.

maxDist

Maximal length of the gap (in base-pairs) between two tag clusters for them to be part of the same consensus clusters.

useMulticore

Logical, should multicore be used (supported only on Unix-like platforms).

nrCores

Number of cores to use when useMulticore = TRUE. Default (NULL) uses all detected cores.

Details

Since the tag clusters (TCs) returned by the CTSS clustering functions function are constructed separately for every CAGE sample within the CAGEr object, they can differ between samples in both their number, genomic coordinates, position of dominant TSS and overall signal. To be able to compare all samples at the level of clusters of TSSs, TCs from all CAGE datasets are aggregated into a single set of consensus clusters. First, TCs with signal ⁠>= tpmThreshold⁠ from all CAGE datasets are selected, and their 5' and 3' boundaries are determined based on provided qLow and qUp parameter (or the start and end coordinates, if they are set to NULL). Finally, the defined set of TCs from all CAGE datasets is reduced to a non-overlapping set of consensus clusters by merging overlapping TCs and TCs ⁠<= maxDist⁠ base-pairs apart. Consensus clusters represent a referent set of promoters that can be further used for expression profiling or detecting "shifting" (differentially used) promoters between different CAGE samples.

Value

Returns the object in which the experiment consensusClusters will be occupied by a RangedSummarizedExperiment containing the cluster coordinates as row ranges, and their expression levels in the counts and normalized assays. These genomic ranges are returned by the consensusClustersGR function and the whole object can be accessed with the consensusClustersSE function. The CTSS ranges of the tagCountMatrix experiment will gain a cluster column indicating which cluster they belong to. Lastly, the number of CTSS outside clusters will be documented in the outOfClusters column data.

Author(s)

Vanja Haberle

Charles Plessy

See Also

Other CAGEr object modifiers: CTSStoGenes(), CustomConsensusClusters(), annotateCTSS(), cumulativeCTSSdistribution(), distclu(), getCTSS(), normalizeTagCount(), paraclu(), quantilePositions(), quickEnhancers(), resetCAGEexp(), summariseChrExpr()

Other CAGEr clusters functions: CTSScumulativesTagClusters(), CustomConsensusClusters(), consensusClustersDESeq2(), consensusClustersGR(), cumulativeCTSSdistribution(), distclu(), paraclu(), plotInterquantileWidth(), quantilePositions(), tagClustersGR()

Examples


consensusClustersGR(exampleCAGEexp)
ce <- aggregateTagClusters( exampleCAGEexp, tpmThreshold = 50
                          , excludeSignalBelowThreshold = FALSE, maxDist = 100)
consensusClustersGR(ce)

ce <- aggregateTagClusters( exampleCAGEexp, tpmThreshold = 50
                          , excludeSignalBelowThreshold = TRUE, maxDist = 100)
consensusClustersGR(ce)

ce <- aggregateTagClusters( exampleCAGEexp, tpmThreshold = 50
                          , excludeSignalBelowThreshold = TRUE, maxDist = 100
                          , qLow = 0.1, qUp = 0.9)
consensusClustersGR(ce)


charles-plessy/CAGEr documentation built on Aug. 2, 2024, 4:35 p.m.