[BioC] comparing seqnames

Murli Nair [guest] guest at bioconductor.org
Sat Jun 15 19:17:27 CEST 2013


Hi,
I am trying to extracts mate pairs that are mapped to a particular chromosome (say chr1). The construction of my which() function is not returning what I need. May be I am missing something in constructing it. 

bf.data= readGappedAlignments(bam_files, param=ScanBamParam(what=scanBamWhat()))
mate.pairs=table(mcols(bf.data)$qname)
onlyPairs=names(mate.pairs)[mate.pairs==2]
mappedPairs=bf.data[mcols(bf.data)$qname %in% onlyPairs]
mate1=mappedPairs[c(T,F)]
mate2=mappedPairs[c(F,T)]
isSameCzome= (seqnames(mate1)==seqnames(mate2))

chrnames=seqnames(mate1[which(seqnames(mate1[isSameCzome])=="chr1")])

With the last statement I am intending to extract mate1 that are mapped to chr1. So I expect chrnames to contain only chr1, but it returns all other chromosomes and not just for chr1.  Is there something I am missing in the way I have written the which statement?

Thanks for your help.

Cheers../Murli



 -- output of sessionInfo(): 

> sessionInfo()
R version 3.0.1 (2013-05-16)
Platform: x86_64-redhat-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=C                 LC_NAME=C                 
 [9] LC_ADDRESS=C               LC_TELEPHONE=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] VariantAnnotation_1.6.6                
 [2] org.Hs.eg.db_2.9.0                     
 [3] RSQLite_0.11.4                         
 [4] DBI_0.2-7                              
 [5] BSgenome.Hsapiens.UCSC.hg19_1.3.19     
 [6] BSgenome_1.28.0                        
 [7] TxDb.Hsapiens.UCSC.hg19.knownGene_2.9.2
 [8] GenomicFeatures_1.12.2                 
 [9] AnnotationDbi_1.22.6                   
[10] Biobase_2.20.0                         
[11] Rsamtools_1.12.3                       
[12] Biostrings_2.28.0                      
[13] GenomicRanges_1.12.4                   
[14] IRanges_1.18.1                         
[15] BiocGenerics_0.6.0                     

loaded via a namespace (and not attached):
[1] biomaRt_2.16.0     bitops_1.0-5       RCurl_1.95-4.1     rtracklayer_1.20.2
[5] stats4_3.0.1       tools_3.0.1        XML_3.96-1.1       zlibbioc_1.6.0

--
Sent via the guest posting facility at bioconductor.org.



More information about the Bioconductor mailing list