Hi Xiaobei,

Thanks for reporting this bug! The problem was with how the "*" strand 
was assigned to a GRangesList when ignore.strand = TRUE. This is now 
fixed in GenomicRanges in devel 1.9.45 and release 1.8.12 in release. 
Both versions will be available through biocLite() tomorrow morning PST 
or through svn immediately.

Valerie


On 08/10/2012 09:07 AM, xiaobei.zhou@uzh.ch wrote:
> Hello everyone,
>
> Maybe "I've found an error running summarizeOverlaps()" ...
>
> This is what the example in summarizeOverlaps() does:
> -------------------------------------------------------------------------------------------------------------------------
>> Ā 
>> ## paired-end readsĀ 
>> ##
>> Ā 
>> library("TxDb.Dmelanogaster.UCSC.dm3.ensGene")
>> exbygene<- exonsBy(TxDb.Dmelanogaster.UCSC.dm3.ensGene, "gene")
>> fl<- system.file("extdata", "untreated3_chr4.bam",Ā 
>> Ā  Ā  package="pasillaBamSubset")Ā 
>> Ā 
>> ## Paired-end reads are stored in a GappedAlignmentPairs object
>> reads<- readGappedAlignmentPairs(fl)
>> Ā 
>> res<- summarizeOverlaps(exbygene, reads, "Union")
>> stopifnot(length(assays(res)$counts) == length(exbygene))
>> Ā 
>> -------------------------------------------------------------------------------
>
> I run exactly the same example but with ignore.strand = TRUE.Ā 
> Then it gives an error:
>
>
>>> library("TxDb.Dmelanogaster.UCSC.dm3.ensGene")
>>> exbygene<- exonsBy(TxDb.Dmelanogaster.UCSC.dm3.ensGene, "gene")
>>> fl<- system.file("extdata", "untreated3_chr4.bam",Ā
>> + Ā  Ā  package="pasillaBamSubset")Ā
>>> Ā 
>>> ## Paired-end reads are stored in a GappedAlignmentPairs object
>>> reads<- readGappedAlignmentPairs(fl)
>> Ā
>>>> res<- summarizeOverlaps(exbygene, reads, "Union", ignore.strand = TRUE)
>>> Error in .local(x, ..., value) :Ā
>> Ā  replacement 'value' is not an AtomicList with the same elementLengths as 'x'
>>>> traceback()
>>> 7: stop("replacement 'value' is not an AtomicList with the same ",Ā
>> Ā  Ā  Ā  Ā "elementLengths as 'x'")
>>> 6: .local(x, ..., value)
>>> 5: `strand<-`(`*tmp*`, value = "*")
>>> 4: `strand<-`(`*tmp*`, value = "*")
>>> 3: .dispatchOverlaps(grglist(reads), features, mode, ignore.strand)
>>> 2: summarizeOverlaps(exbygene, reads, "Union", ignore.strand = TRUE)
>>> 1: summarizeOverlaps(exbygene, reads, "Union", ignore.strand = TRUE)
>> Ā 
>> Ā
>>>> Ā sessionInfo()
>>> R version 2.15.1 (2012-06-22)
>>> Platform: x86_64-unknown-linux-gnu (64-bit)
>> Ā
>>> locale:
>>> [1] LC_CTYPE=en_CA.UTF-8 Ā  Ā  Ā  LC_NUMERIC=C Ā  Ā  Ā  Ā  Ā  Ā  Ā 
>>> [3] LC_TIME=en_CA.UTF-8 Ā  Ā  Ā  Ā LC_COLLATE=en_CA.UTF-8 Ā  Ā 
>>> [5] LC_MONETARY=en_CA.UTF-8 Ā  Ā LC_MESSAGES=en_CA.UTF-8 Ā Ā 
>>> [7] LC_PAPER=C Ā  Ā  Ā  Ā  Ā  Ā  Ā  Ā  LC_NAME=C Ā  Ā  Ā  Ā  Ā  Ā  Ā  Ā Ā 
>>> [9] LC_ADDRESS=C Ā  Ā  Ā  Ā  Ā  Ā  Ā  LC_TELEPHONE=C Ā  Ā  Ā  Ā  Ā  Ā 
>>> [11] LC_MEASUREMENT=en_CA.UTF-8 LC_IDENTIFICATION=C Ā  Ā  Ā Ā
>> Ā
>>> attached base packages:
>>> [1] stats Ā  Ā  graphics Ā grDevices utils Ā  Ā  datasets Ā methods Ā  base Ā  Ā Ā
>> Ā
>>> other attached packages:
>>> [1] Rsamtools_1.9.25 Ā  Ā  Ā  Ā  Ā  Ā  Ā  Ā  Ā  Ā  Ā  Ā Ā 
>>> [2] Biostrings_2.25.8 Ā  Ā  Ā  Ā  Ā  Ā  Ā  Ā  Ā  Ā  Ā  Ā 
>>> [3] rtracklayer_1.17.15 Ā  Ā  Ā  Ā  Ā  Ā  Ā  Ā  Ā  Ā  Ā 
>>> [4] TxDb.Dmelanogaster.UCSC.dm3.ensGene_2.7.1
>>> [5] GenomicFeatures_1.9.28 Ā  Ā  Ā  Ā  Ā  Ā  Ā  Ā  Ā Ā 
>>> [6] AnnotationDbi_1.19.28 Ā  Ā  Ā  Ā  Ā  Ā  Ā  Ā  Ā  Ā 
>>> [7] Biobase_2.17.6 Ā  Ā  Ā  Ā  Ā  Ā  Ā  Ā  Ā  Ā  Ā  Ā  Ā Ā 
>>> [8] GenomicRanges_1.9.42 Ā  Ā  Ā  Ā  Ā  Ā  Ā  Ā  Ā  Ā Ā 
>>> [9] IRanges_1.15.29 Ā  Ā  Ā  Ā  Ā  Ā  Ā  Ā  Ā  Ā  Ā  Ā  Ā 
>>> [10] BiocGenerics_0.3.0 Ā  Ā  Ā  Ā  Ā  Ā  Ā  Ā  Ā  Ā  Ā Ā
>> Ā
>>> loaded via a namespace (and not attached):
>>> [1] biomaRt_2.13.2 Ā bitops_1.0-4.1 Ā BSgenome_1.25.6 DBI_0.2-5 Ā  Ā  Ā 
>>> [5] RCurl_1.91-1 Ā  Ā RSQLite_0.11.1 Ā stats4_2.15.1 Ā  tools_2.15.1 Ā Ā 
>>> [9] XML_3.9-4 Ā  Ā  Ā  zlibbioc_1.3.0Ā
> -----------------------------------------------------------------------------------------
> Ā 
>
> Help much appreciated.
> Ā 
> Best,
> Xiaobei
>
>
> =
> 	[[alternative HTML version deleted]]
>
>
>
> _______________________________________________
> Bioc-devel@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/bioc-devel


	[[alternative HTML version deleted]]

