Файл BAM: получение всех прочтений на определенной позиции с помощью pysam

Question

У меня есть BAM-файл с определенной позицией чтения 520817 (как видно из 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()

Answer 1

Проверьте опцию IGV View ->Preference->Alignment, некоторые опции "filter xxxx" (дубликаты, вторичные выравнивания, низкое качество) могут изменить вывод.

Обычно pysam не накапливает чтения с флагами BAM_FUNMAP, BAM_FSECONDARY, BAM_FQCFAIL, BAM_FDUP, поэтому убедитесь, что параметры просмотра IGV совпадают с параметрами в pysam.AlignmentFile.pileup. В противном случае они могут генерировать разные результаты.