Creating A Granges Object From Rleviewslist
3
4
Entering edit mode
7.5 years ago

I am trying to create a GRanges object form package GenomicRanges from an RleViewsList but I am getting the following error >Error in normalizeSingleBracketSubscript(i, x) : invalid subscript type. Do you know why this is happening? More details about the objects and the R session below:

Thank you very much!

The RleViewsList has been created using the "coverage" function then the "slice" function to filter areas with a continuous coverage >3. A Granges object was first created with the "subsetByOverlaps" function:

>inp
GRanges with 4162112 ranges and 0 metadata columns:
seqnames                 ranges strand
<Rle>              <IRanges>  <Rle>
[1]       10       [111068, 111109]      *
[2]       10       [118972, 119013]      *
[3]       10       [119395, 119436]      *
[4]       10       [119879, 119920]      *
[5]       10       [119881, 119922]      *
...      ...                    ...    ...
[4162108]        X [154509640, 154509681]      *
[4162109]        X [154509645, 154509686]      *
[4162110]        X [154509645, 154509686]      *
[4162111]        X [154595860, 154595901]      *
[4162112]        X [154595873, 154595914]      *
---
seqlengths:
10 11 12 13 14 15 16 17 18 19  1 20  2  3  4  5  6  7  8  9  X
NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA

>cov <- coverage(inp)
>rle.i <- slice(cov, lower = 3)
>reduce(as(rle.i, "GRanges"))
>Error in normalizeSingleBracketSubscript(i, x) : invalid subscript type.

>rle.i
RleViewsList of length 21
names(21): 10 11 12 13 14 15 16 17 18 19 1 20 2 3 4 5 6 7 8 9 X
>rle.i[[1]]
Views on a 112160197-length Rle subject

views:
start       end width
[1]   1658660   1658715    56 [3 3 3 3 3 3 5 4 4 4 4 4 3 3 5 5 5 6 6 6 ...]
[2]   1658738   1658783    46 [4 4 4 4 5 5 5 5 5 6 6 6 6 6 7 7 7 7 7 7 ...]
[3]   1658819   1658868    50 [3 3 3 3 3 3 4 4 5 5 5 5 5 5 5 5 5 5 5 6 ...]
[4]   1658984   1659042    59 [3 3 4 4 4 4 4 4 4 4 5 5 5 5 5 5 5 6 6 6 ...]
[5]   1659044   1659099    56 [3 3 3 3 3 3 3 4 4 4 4 4 4 4 6 6 6 6 6 6 ...]
[6]   1659116   1659308   193 [3 3 3 3 3 3 5 5 5 5 5 5 5 5 5 5 4 5 5 5 ...]
[7]   1659330   1659377    48 [3 3 3 3 3 3 4 4 4 4 4 5 5 5 5 5 5 5 5 5 ...]
[8]   2240603   2240845   243 [ 4  4  4  4  4  4  7  7  7  7  7  7  7 ...]
[9]   2240847   2240906    60 [ 4  4  5  5  5  5  9  9 11 11 11 12 12 ...]
...       ...       ...   ... ...
[2147] 110693912 110693978    67 [3 3 3 3 3 3 4 4 4 4 4 4 4 4 4 4 5 5 5 5 ...]
[2148] 110694171 110694215    45 [4 4 4 5 4 4 4 4 4 4 4 5 5 5 5 5 5 5 5 5 ...]
[2149] 110694440 110694488    49 [4 4 5 5 5 5 5 6 6 5 5 5 5 5 5 5 5 5 5 5 ...]
[2150] 110694497 110694566    70 [3 3 3 3 3 3 3 4 4 4 5 5 5 5 5 5 5 5 5 4 ...]
[2151] 110694615 110694670    56 [3 3 3 4 4 4 4 4 4 4 4 4 4 5 6 7 7 7 7 7 ...]
[2152] 110694675 110694793   119 [3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 5 5 5 5 5 ...]
[2153] 112071553 112071610    58 [ 3  3  3  3  3  3  3  3  3  4  5  8  9 ...]
[2154] 112080357 112080432    76 [ 4  4  4  6  6  6  6  6  6  6  6  6  6 ...]
[2155] 112146168 112146211    44 [3 3 4 4 4 4 4 4 4 5 6 6 6 6 6 6 6 6 6 6 ...]

> sessionInfo()
R version 3.0.2 (2013-09-25)
Platform: x86_64-unknown-linux-gnu (64-bit)

locale:
[1] LC_CTYPE=en_US.UTF-8       LC_NUMERIC=C
[3] LC_TIME=en_US.UTF-8        LC_COLLATE=en_US.UTF-8
[5] LC_MONETARY=en_US.UTF-8    LC_MESSAGES=en_US.UTF-8
[7] LC_PAPER=en_US.UTF-8       LC_NAME=C
[11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C

attached base packages:
[1] parallel  stats     graphics  grDevices utils     datasets  methods
[8] base

other attached packages:
[1] rtracklayer_1.22.0     GenomicFeatures_1.14.0 AnnotationDbi_1.24.0
[4] Biobase_2.22.0         Rsamtools_1.14.1       Biostrings_2.30.0
[7] GenomicRanges_1.14.1   XVector_0.2.0          IRanges_1.20.0
[10] BiocGenerics_0.8.0

loaded via a namespace (and not attached):
[1] biomaRt_2.18.0  bitops_1.0-6    BSgenome_1.30.0 DBI_0.2-7
[5] RCurl_1.95-4.1  RSQLite_0.11.4  stats4_3.0.2    tools_3.0.2
[9] XML_3.98-1.1    zlibbioc_1.8.0

r bioconductor coverage • 3.6k views
1
Entering edit mode
7.5 years ago

The answer is: "Works for me!"

using a small dataset from example(GRanges). Would it be possible to save and give us the complete GRanges object inp? And could you test the following code?

 example(GRanges)
mycov = coverage(longGR)
rle.i = slice(mycov, lower=3)
gr = as (rle.i, "GRanges")
reduce(gr)


Output:

 GRanges with 6 ranges and 0 metadata columns:
seqnames     ranges strand
<Rle>  <IRanges>  <Rle>
[1]     chr1 [  2,  10]      *
[2]     chr1 [103, 110]      *
[3]     chr2 [  4,  10]      *
[4]     chr2 [106, 114]      *
[5]     chr3 [  8,  12]      *
[6]     chr3 [110, 117]      *
---
seqlengths:
chr1 chr2 chr3
NA   NA   NA


My versions are a bit older, possibly a bug was introduced recently?

> sessionInfo()
R version 3.0.1 (2013-05-16)
Platform: x86_64-apple-darwin10.8.0 (64-bit)

locale:
[1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8

attached base packages:
[1] parallel  grid      stats     graphics  grDevices utils     datasets  methods   base

other attached packages:
[1] GenomicRanges_1.12.4 IRanges_1.18.1       BiocGenerics_0.6.0   cluster_1.14.4       gplots_2.11.3        MASS_7.3-26
[7] KernSmooth_2.23-10   caTools_1.14         gdata_2.12.0.2       gtools_2.7.1         biomaRt_2.16.0

loaded via a namespace (and not attached):
[1] bitops_1.0-5   RCurl_1.95-4.1 splines_3.0.1  stats4_3.0.1   tools_3.0.1    XML_3.95-0.2

1
Entering edit mode
7.5 years ago

Thank you very much! With the following versions it works also for me:

R version 3.0.2 (2013-09-25)
Platform: x86_64-pc-linux-gnu (64-bit)

locale:
[1] LC_CTYPE=en_GB.UTF-8       LC_NUMERIC=C
[3] LC_TIME=en_GB.UTF-8        LC_COLLATE=en_GB.UTF-8
[5] LC_MONETARY=en_GB.UTF-8    LC_MESSAGES=en_GB.UTF-8
[7] LC_PAPER=en_GB.UTF-8       LC_NAME=C
[11] LC_MEASUREMENT=en_GB.UTF-8 LC_IDENTIFICATION=C

attached base packages:
[1] parallel  stats     graphics  grDevices utils     datasets  methods
[8] base

other attached packages:
[1] rtracklayer_1.20.4     GenomicFeatures_1.12.4 AnnotationDbi_1.22.6
[4] Biobase_2.20.1         Rsamtools_1.12.4       Biostrings_2.28.0
[7] GenomicRanges_1.12.5   IRanges_1.18.4         BiocGenerics_0.6.0

loaded via a namespace (and not attached):
[1] biomaRt_2.16.0  bitops_1.0-6    BSgenome_1.28.0 DBI_0.2-7
[5] RCurl_1.95-4.1  RSQLite_0.11.4  stats4_3.0.2    tools_3.0.2
[9] XML_3.98-1.1    zlibbioc_1.6.0


So it should be a recently introduced bug. Interesting...

0
Entering edit mode

Could you please try the small example code I posted too? This can serve as a reproducible example when submitting a bug report.

1
Entering edit mode
7.5 years ago

Confirmed, I am getting the error only with the newer versions. The complete session info of the versions which give the error is:

 > reduce(as(rle.i, "GRanges"))
Error in reduce(as(rle.i, "GRanges")) :
error in evaluating the argument 'x' in selecting a method for function 'reduce': Error in normalizeSingleBracketSubscript(i, x) : invalid subscript type
> sessionInfo()
R version 3.0.2 (2013-09-25)
Platform: x86_64-unknown-linux-gnu (64-bit)

locale:
[1] LC_CTYPE=en_US.UTF-8       LC_NUMERIC=C
[3] LC_TIME=en_US.UTF-8        LC_COLLATE=en_US.UTF-8
[5] LC_MONETARY=en_US.UTF-8    LC_MESSAGES=en_US.UTF-8
[7] LC_PAPER=en_US.UTF-8       LC_NAME=C
[11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C

attached base packages:
[1] parallel  stats     graphics  grDevices utils     datasets  methods
[8] base

other attached packages:
[1] rtracklayer_1.22.0     GenomicFeatures_1.14.0 AnnotationDbi_1.24.0
[4] Biobase_2.22.0         Rsamtools_1.14.1       Biostrings_2.30.0
[7] GenomicRanges_1.14.1   XVector_0.2.0          IRanges_1.20.0
[10] BiocGenerics_0.8.0

loaded via a namespace (and not attached):
[1] biomaRt_2.18.0  bitops_1.0-6    BSgenome_1.30.0 DBI_0.2-7
[5] RCurl_1.95-4.1  RSQLite_0.11.4  stats4_3.0.2    tools_3.0.2
[9] XML_3.98-1.1    zlibbioc_1.8.0


Is there any way to report this bug? I think it has to do with the GenomiRangesn package version.

1
Entering edit mode

Yes, please send a report to the mailling list bioconductor@stat.math.ethz.ch giving the reproducible code example that can be run from example data only (without needing to share your data). You can also point to this thread for documentation.