[Bioc-devel] odd param behavior with scanBam from Rsamtools 1.14.1

Martin Morgan mtmorgan at fhcrc.org
Thu Oct 31 21:02:48 CET 2013


On 10/31/2013 12:54 PM, Jason Byars wrote:
> Using any permutation of ScanBamParam() options seems to work fine with
countBam, filterBam, etc. But, when I use a ScanBamParam object without anything
specified for "what" with scanBam, I get a length 0 list back. If all I add is
what=ScanBamWhat(), the call appears to work fine. I have tried this on both
Windows and Linux, and with a clean install of R and Bioconductor 2.13. Is this
behavior correct? Example script follows:

yes, its the way it was designed, however poorly -- you don't ask for any 
fields, so you don't get any fields. countBam and filterBam are about visiting 
records, and records exist even when none of their fields are consulted. It's like

 > data.frame(x=1:10)[,FALSE, drop=FALSE]
data frame with 0 columns and 10 rows

where there are 10 rows (records), but no columns (fields).

Martin

>
> library(Rsamtools)
>
> # param seems to work find on filterBam, but I'm getting no results on scanBam on Linux
>
> sampleFile <- "F:/nessstuff/rnaseq/mybunaligned/P61_3_MBC24-6-1-M3_mybunaligned.bam"
>
> # just picked MYB because we work on it. It's not important to show the problem.
> mybStart <- 135502453
> mybEnd <- 135540311
> mybSize <- 37859
> padding <- 10000
> flag <- scanBamFlag(isUnmappedQuery=F)
> which <- GRanges(seqname=Rle(c("chr6"), c(1)), ranges=IRanges(mybStart-padding, mybEnd+padding))
> # what <- c("rname", "strand", "pos", "qwidth", "seq")
> param <- ScanBamParam(flag=flag,which=which)
>
> countBam(sampleFile) #813k records
> countBam(sampleFile,param=param) #2397 records
>
> # 813k records loaded
> bam <- scanBam(sampleFile)[[1]]
> summary(bam)
>
> # 0 length list returned
> bam <- scanBam(sampleFile,param=param)[[1]]
> summary(bam)
>
> # specify anything for what and it works as expected #2397 records imported
> param <- ScanBamParam(flag=flag,which=which, what=scanBamWhat())
> bam <- scanBam(sampleFile,param=param)[[1]]
> summary(bam)
>
>> sessionInfo()
> R version 3.0.2 (2013-09-25)
> Platform: x86_64-w64-mingw32/x64 (64-bit)
>
> locale:
> [1] LC_COLLATE=English_United States.1252  LC_CTYPE=English_United States.1252
> [3] LC_MONETARY=English_United States.1252 LC_NUMERIC=C
> [5] LC_TIME=English_United States.1252
>
> attached base packages:
> [1] parallel  stats     graphics  grDevices utils     datasets  methods   base
>
> other attached packages:
> [1] Rsamtools_1.14.1     Biostrings_2.30.0    GenomicRanges_1.14.3 XVector_0.2.0        IRanges_1.20.4
> [6] BiocGenerics_0.8.0   BiocInstaller_1.12.0
>
> loaded via a namespace (and not attached):
> [1] bitops_1.0-6   stats4_3.0.2   tools_3.0.2    zlibbioc_1.8.0
>
>
>
> 	[[alternative HTML version deleted]]
>
> _______________________________________________
> Bioc-devel at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/bioc-devel
>


-- 
Computational Biology / Fred Hutchinson Cancer Research Center
1100 Fairview Ave. N.
PO Box 19024 Seattle, WA 98109

Location: Arnold Building M1 B861
Phone: (206) 667-2793



More information about the Bioc-devel mailing list