genomicElementUpSetR: Genomic Element data for upset plot

Description Usage Arguments Details Value Examples

View source: R/genomicElementUpSetR.R

Description

Prepare data for upset plot for genomic element distribution

Usage

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
genomicElementUpSetR(
  peaks,
  TxDb,
  seqlev,
  ignore.strand = TRUE,
  breaks = list(distal_upstream = c(-1e+05, -10000, -1, 1), proximal_upstream =
    c(-10000, -5000, -1, 1), distal_promoter = c(-5000, -2000, -1, 1), proximal_promoter
    = c(-2000, 200, -1, 0), `5'UTR` = fiveUTRsByTranscript, `3'UTR` =
    threeUTRsByTranscript, CDS = cds, exon = exons, intron = intronsByTranscript,
    gene_body = genes, immediate_downstream = c(0, 2000, 1, 1), proximal_downstream =
    c(2000, 5000, 1, 1), distal_downstream = c(5000, 1e+05, 1, 1))
)

Arguments

peaks

peak list, GRanges object or a GRangesList.

TxDb

an object of TxDb

seqlev

sequence level should be involved. Default is all the sequence levels in intersect of peaks and TxDb.

ignore.strand

logical. Whether the strand of the input ranges should be ignored or not. Default=TRUE

breaks

list. A list for labels and sets for the genomic elements. The element could be an S4 method for signature 'TxDb' or a numeric vector with length of 4. The three numbers are c(upstream point, downstream point, promoter (-1) or downstream (1), remove gene body or not (1: remove, 0: keep)).

Details

The data will be calculated by for each breaks. No precedence will be considered.

Value

list of data for plot.

Examples

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
if (interactive() || Sys.getenv("USER")=="jianhongou"){
  data(myPeakList)
  if(require(TxDb.Hsapiens.UCSC.hg19.knownGene)){
  seqinfo(myPeakList) <- 
  seqinfo(TxDb.Hsapiens.UCSC.hg19.knownGene)[seqlevels(myPeakList)]
  myPeakList <- GenomicRanges::trim(myPeakList)
  myPeakList <- myPeakList[width(myPeakList)>0]
  x <- genomicElementUpSetR(myPeakList, 
    TxDb.Hsapiens.UCSC.hg19.knownGene)
  library(UpSetR)
  upset(x$plotData, nsets=13, nintersects=NA)
  }
}

ChIPpeakAnno documentation built on April 1, 2021, 6:01 p.m.