utils {GenomicRanges}R Documentation

seqlevels utility functions

Description

Rename or subset the seqlevels in a GenomicRanges, GRangesList or GappedAlignments object.

Usage

  ## S4 method for signature 'GenomicRanges,character'
keepSeqlevels(x, value, ...) 
  ## S4 method for signature 'GenomicRanges,character'
renameSeqlevels(x, value, ...) 

Arguments

x

The GenomicRanges, GRangesList or GappedAlignments object in which the seqlevels will be removed or renamed.

value

For keepSeqlevels, a GRanges, GRangesList or GappedAlignments or character vector. x is subset on the seqlevels in value.

For renameSeqlevels, a named character vector where the names are the ‘old’ and the values are the ‘new’ seqlevels.

...

Arguments passed to other functions.

Details

Many operations on GRanges objects require the seqlevels to match before a comparison can be made (e.g., findOverlaps(type="within")). keepSeqlevels and renameSeqlevels are convenience functions for subsetting and renaming the seqlevels of these objects.

keepSeqlevels subsets x based on the seqlevels provided in value. If value does not match any seqlevels in x the original x is returned. When x is a GRangesList, there may be multiple chromosomes in a single list element. If not all chromosomes are specified in value, a reduced list element is returned. All empty list elements are dropped. See examples.

renameSeqlevels renames the seqlevels in x to those provided in value. value is a named character vector where the names are the ‘old’ seqlevel names and the values are the ‘new’. If no names in value match the seqlevels in x the original x is returned.

Value

The x object with seqlevels removed or renamed. If x has no seqlevels (empty object) or no replacement values match the current seqlevels in x the unchanged x is returned.

Author(s)

Valerie Obenchain vobencha@fhcrc.org

See Also

VCF

Examples

gr1 <- GRanges(seqnames = c("chr1", "chr2"),
               ranges = IRanges(c(7,13), width = 3),
               strand = c("+", "-"), score = 3:4, GC = c(0.3, 0.5))
gr2 <- GRanges(seqnames = c("chr1", "chr1", "chr2", "chr3", "chr3"),
               ranges = IRanges(c(1, 4, 8, 9, 16), width=5),
               strand = "-", score = c(3L, 2L, 5L, 6L, 2L), 
               GC = c(0.4, 0.1, 0.55, 0.20, 0.10))
gr3 <- GRanges(seqnames = c("CHROM4", "CHROM4"), 
               ranges = IRanges(c(20, 45), width=6), 
               strand = "+", score = c(2L, 5L), GC = c(0.30, 0.45))
 
## GRanges :
gr3_rename <- renameSeqlevels(gr3, c(CHROM4="chr4"))
gr3_rename
 
gr2_subset_chr <- keepSeqlevels(gr2, c("chr1", "chr2"))
gr2_subset_gr <- keepSeqlevels(gr2, gr1)
identical(gr2_subset_chr, gr2_subset_gr)
 
## GRangesList :
grl1 <- GRangesList("gr1" = gr1, "gr2" = gr2, "gr3" = gr3)
grl2 <- GRangesList("gr1" = gr1, "gr2" = gr2, "gr3" = gr3_rename)
grl1_rename <- renameSeqlevels(grl1, c(CHROM4="chr4"))
identical(grl1_rename, grl2)
 
grl1_subset <- keepSeqlevels(grl1, "chr3")
 
## GappedAlignments :
library(Rsamtools)
galn_file <- system.file("extdata", "ex1.bam", package="Rsamtools")
galn <- readGappedAlignments(galn_file)
 
galn_rename <- renameSeqlevels(galn, c(seq2="chr2"))
galn_subset <- keepSeqlevels(galn_rename, gr1)
galn_subset

## See ?VCF for examples using renameSeqlevels and keepSeqlevels with 
## VCF class objects

[Package GenomicRanges version 1.12.1 Index]