從bam文件中提取讀取位置

Question

我有一個包含多個SNP的vcf文件，現在我想看看，這些SNP是否均勻分佈在從中獲得SNP的bam文件的讀取中。具體而言，我想繪製讀取位置上的SNP數量。 我想知道是否有這樣做的工具，或者我是否必須自己寫一個腳本。如果是這樣，R中是否有一個包可以做到這一點（我已經習慣了R但是沒有太多的perl經驗）？

Answer 1

還不確定是什麼在讀取位置的SNP「是指，但你可以閱讀與R/Bioconductor包和功能VariantAnnotation :: readVcf的VCF，並使用基因組座標查詢與Rsamtools :: countBam巴姆文件，使用ScanBamParam。未經測試，沿

## first-time installation 
source("http://bioconductor.org/biocLite.R") 
biocLite(c("VariantAnnotation", "Rsamtools")) 

線安裝相關的包，然後

library(VariantAnnotation) # also loads Rsamtools 
snps = readVcf("/some/file.vcf") 
param = ScanBamParam(which=rowData(vcf)) 
reads = countBam("/some/file.bam", param=param) 

最好的方式來實現，這可能取決於你有多少個SNP有興趣了很多東西。我」 d建議您使用預發佈的R-2.15 alpha版，因爲您將獲得更新的Bioconductor套件。這些包有廣泛的護身符（vignette(package="VariantAnnotation")和知識的人在Bioconductor的mailing list，以及平時的幫助頁面?readVcf