transcriptsByOverlaps {GenomicFeatures}R Documentation

Extract genomic features from an object based on their by genomic location

Description

Generic functions to extract genomic features for specified genomic locations. This page documents the methods for TranscriptDb objects only.

Usage

transcriptsByOverlaps(x, ranges,
                      maxgap = 0L, minoverlap = 1L,
                      type = c("any", "start", "end"), ...)
## S4 method for signature 'TranscriptDb'
transcriptsByOverlaps(x, ranges,
                      maxgap = 0L, minoverlap = 1L,
                      type = c("any", "start", "end"),
                      columns = c("tx_id", "tx_name"))

exonsByOverlaps(x, ranges,
                maxgap = 0L, minoverlap = 1L,
                type = c("any", "start", "end"), ...)
## S4 method for signature 'TranscriptDb'
exonsByOverlaps(x, ranges,
                maxgap = 0L, minoverlap = 1L,
                type = c("any", "start", "end"),
                columns = "exon_id")

cdsByOverlaps(x, ranges,
              maxgap = 0L, minoverlap = 1L,
              type = c("any", "start", "end"), ...)
## S4 method for signature 'TranscriptDb'
cdsByOverlaps(x, ranges,
              maxgap = 0L, minoverlap = 1L,
              type = c("any", "start", "end"),
              columns = "cds_id")

Arguments

x

A TranscriptDb object.

...

Arguments to be passed to or from methods.

ranges

A GRanges object to restrict the output.

type

How to perform the interval overlap operations of the ranges. See the findOverlaps manual page in the GRanges package for more information.

maxgap

A non-negative integer representing the maximum distance between a query interval and a subject interval.

minoverlap

Ignored.

columns

Columns to include in the output. See ?transcripts for the possible values.

Details

These functions subset the results of transcripts, exons, and cds function calls with using the results of findOverlaps calls based on the specified ranges.

Value

a GRanges object

Author(s)

P. Aboyoun

See Also

Examples

  txdb <- loadDb(system.file("extdata", "UCSC_knownGene_sample.sqlite",
                                   package="GenomicFeatures"))
  gr <- GRanges(seqnames = rep("chr1",2),
                ranges = IRanges(start=c(500,10500), end=c(10000,30000)),
                strand = strand(rep("-",2)))
  transcriptsByOverlaps(txdb, gr)

[Package GenomicFeatures version 1.14.2 Index]