pysam使用心得

pysam以及Samtools使用心得

0. pysam簡介

pysam 是一個基於 htslib 的 C++ API 進行封裝的 python 模塊,實現了對 SAM / BAM / CRAM 文件的便捷操作,可以簡化 bam 文件處理的代碼複雜度，同時也可以處理 VCF / BCF 等其他文件

1. pysam 的安裝

pysam 已經在 pypi 中包含，可以直接使用 命令安裝

2. pysam 的使用

2.1 統計 bam 文件中總 reads 數目以及比對到參考基因組上的 reads 數目

利用 對象，我們可以使用 方法快速統計其中的比對片段數目，然而由於存在一條 read 比對到基因組中多個位置的情況，因此直接統計的比對片段數目與 reads 數是不同的

可以使用 參數，使用callback函數對比對片段進行篩選，具體細節請參考pysam API

–

select a call-back to ignore reads when counting. It can be either a string with the following values:

-

skip reads in which any of the following flags are set: BAM_FUNMAP, BAM_FSECONDARY, BAM_FQCFAIL, BAM_FDUP

-

uses every single read

使用結合自定義callback函數，我們可以方便地統計總 reads 數目以及比對到基因組上的reads數

Note: 默認情況下， 只統計比對上的片段，因此統計 total reads 數目時需要指定 參數，才可以包含未比對上的片段信息

下面是一個使用 pysam 統計 total reads 與 mapped reads 數目的實常式序：

2.2 pysam 多進程處理 bam 文件

pysam中多進程處理需要注意一些事情

- 如何多進程處理bam文件

通過方法，可以提取出比對到指定染色體或者指定區域的比對片段，因此可以多進程同時處理多個染色體上的比對結果

- 注意事項

利用同一個句柄同時進行操作可能發生衝突，因此每個進程需要重新打開自己的新句柄，並使用參數，使得句柄之間相互獨立，可以避免文件讀取時發生問題

對象不能pickle化，因此無法在多進程中作為其他函數的參數，需要自己構建tuple/list用於參數傳遞

如果每個子進程都需要相同的輔助數據，可以使用initializer進行全局對象的共享，效率比為每個進程提供相同參數要高很多

3. pysam與samtools結合使用

在打開bam文件時，一般需要bam文件已經排序並且已建立好index

可以使用以及多進程

需要注意，多進程運行之後生成的index文件的時間戳會出現異常，在程序里需要等待幾秒後再使用pysam繼續處理，否則會提示

喜歡這篇文章嗎？立刻分享出去讓更多人知道吧！