Entering edit mode
Hello,
I am not able to change seqlevels (as I think I was used to) for GRanges objects.
seqlevels (gr) <- seqlevels (noexons1) Fehler in getDanglingSeqlevels(x, new2old = new2old, force = force, seqlevels(value)): won't drop seqlevels currently in use (ch8), unless you use 'force=TRUE' to also drop elements in 'x' where those seqlevels are used (e.g. with 'seqlevels(x, force=TRUE) <- new_seqlevels'). Alternatively, you can also subset 'x' first.
Trying ...
seqlevels (gr, force=TRUE) <- seqlevels (noexons1)
... leaves my with a granges object without ranges :
gr GRanges with 0 ranges and 1 metadata column: seqnames ranges strand | NIC <Rle> <IRanges> <Rle> | <numeric> --- seqlengths: chr1 chr10 chr11 chr12 chr13 chr14 ... chr6 chr7 chr8 chr9 chrX chrY NA NA NA NA NA NA ... NA NA NA NA NA NA
However, after the following procedure it works ...
seqlevels (gr) <- "chr8" seqlevels (gr) <- seqlevels (noexons1)
Can anybody explain this behaviour? What is the "normal" syntax to change seqlevels?
Thanks
Hermann
dput (gr) new("GRanges" , seqnames = new("Rle" , values = structure(1L, .Label = "ch8", class = "factor") , lengths = 4L , elementMetadata = NULL , metadata = list() ) , ranges = new("IRanges" , start = c(25951L, 91936L, 139793L, 155757L) , width = c(1L, 1L, 1L, 1L) , NAMES = NULL , elementType = "integer" , elementMetadata = NULL , metadata = list() ) , strand = new("Rle" , values = structure(3L, .Label = c("+", "-", "*"), class = "factor") , lengths = 4L , elementMetadata = NULL , metadata = list() ) , elementMetadata = new("DataFrame" , rownames = NULL , nrows = 4L , listData = structure(list(NIC = c(0.304211024770134, 0.700993422278323, 0.528789378129407, 0.528789378129407)), .Names = "NIC") , elementType = "ANY" , elementMetadata = NULL , metadata = list() ) , seqinfo = new("Seqinfo" , seqnames = "ch8" , seqlengths = NA_integer_ , is_circular = NA , genome = NA_character_ ) , metadata = list() ) seqlevels (noexons1) [1] "chr1" "chr10" "chr11" "chr12" "chr13" "chr14" "chr15" "chr16" "chr17" [10] "chr18" "chr19" "chr2" "chr20" "chr21" "chr22" "chr3" "chr4" "chr5" [19] "chr6" "chr7" "chr8" "chr9" "chrX" "chrY"
This seems to simply append a "chr" prefix to all seqlevels. This can be a problem for mitochondrial DNA, which Ensembl refers to as MT whereas UCSC refers to as chrM.
Edit: apologies, I was wrong about this. MT DNA name is handled appropriately.