[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