特定の位置に 520817 読み取りの BAM ファイルがあります (IGV で見られるように)。ただし、pysam を使用して特定の位置の読み取り名と関連するヌクレオチドを取得すると、その量ははるかに得られません (約 7000 回の読み取りしか取得できません)。その位置のヌクレオチドが参照ゲノムと異なる場合にのみ読み取りが得られると思います。回避策があるので、すべての読み取りを取得できますか? 私はバイオインフォマティクスから始めています...だから、私を助けるためにもっと必要なものを教えてください!
どうもありがとうございました!
これは私のコードです:
import pysam
import csv
import sys
#---Get a table with in the first column: read-ID; second column: SNP-location; third column: nucleotide---#
mybam = pysam.AlignmentFile("file.bam", "rb")
w = csv.writer(open("snp.csv", "wb"), delimiter=",")
w.writerow(["Read", "Loc", "Nucl"])
for pileupcolumn in mybam.pileup('chr6', 29911198,29911199):
print ("\ncoverage at base %s = %s" %
(pileupcolumn.pos, pileupcolumn.n))
for pileupread in pileupcolumn.pileups:
if not pileupread.is_del:
if pileupcolumn.pos == 29911198:
w.writerow((pileupread.alignment.query_name, 29911198, pileupread.alignment.query_sequence[pileupread.query_position]))
print ('\tbase in read %s = %s' % (pileupread.alignment.query_name, pileupread.alignment.query_sequence[pileupread.query_position]))
mybam.close()