Pysam

Python
#1

Добрый день!
Подскажите, может кто знает как справиться с этой проблемой.
Я работаю с Pysam, с функцией pileup, аналогичная функция есть в инструменте samtools, для того, чтобы сделать вывод в таком же виде, как в samtools (что мне и надо), нужно поставить у pileupcolumn.get_query_sequences значения mark_matches=True, add_indels=True.
И тогда по идее там должны выводиться запятые/точки, совпало ли с референсом (правильным значением в позиции) или нет.

mark_matches (bool) – If True, output bases matching the reference as “,” or “.” for forward and reverse strand, respectively. This mark requires the reference sequence. If no reference is present, this option is ignored.
То есть, чтобы это заработало, необходима reference sequence, та самая референсная последовательность, однако как и где её можно указать, нигде не написано, может кто-нибудь знает, где путь к ней надо прописать?

Вот мой код, если что:

samfile = pysam.AlignmentFile(sampathsort, "rb")

for pileupcolumn in samfile.pileup("chr1",1000,12000):
    print ("\ncoverage at base %s = %s" %
           (pileupcolumn.pos, pileupcolumn.n))
    for query_seq in pileupcolumn.get_query_sequences(mark_matches=True, add_indels=True):
        print(query_seq)
#2

Тогда видимо надо смотреть исходники.

Но поиском по reference_sequence нашлась только эта функция просто принимающая уже прочитанные данные

и несколько функций рядом про создание этой последовательности из каких-то данных.

Так что похоже оно должно передаваться где-то на этапе чтения PileupColumn из файла.

Может этот параметр AlignmentFile?

reference_filename (string) – Path to a FASTA-formatted reference file. Valid only for CRAM files. When reading a CRAM file, this overrides both $REF_PATH and the URL specified in the header (UR tag), which are normally used to find the reference.

#3

Спасибо, нет, этот параметр работает только для CRAM файлов, а у меня BAM файл.

#4

Тогда может для них просто не реализовано еще это.
Можно попробовать спросить разработчиков. Issues · pysam-developers/pysam · GitHub