python PyVCF文件處理VCF文件格式實例詳解

引言

PyVCF庫的安裝

PyVCF庫的導入

PyVCF庫詳細介紹

使用實例：

_Record對象------位點信息的儲存形式

Reader對象------處理vcf文件，構建結構化信息

綜合使用：

vcf文件的全稱是variant call file，即突變識別文件，它是基因組工作流程中產生的一種文件，保存的是基因組上的突變信息。通過對vcf文件進行分析，可以得到個體的變異信息。嗯，總之，這是很重要的文件，所以怎麼處理它也顯得十分重要。它的文件信息如下：

文件的開頭是一堆以“##”開始的注釋行，包含了文件的基本信息。然後是以“#”開頭的一行，共9+n個部分，前九部分標注的是後面行每部分代表的信息，相當於表頭。後面部分是樣本名稱，可以有多個。注釋行結束後是具體的突變信息，每一行分為9+n個部分，每部分之間用制表符（‘\t’）分隔。

通常處理vcf文件時，在讀取，處理階段總是會寫很多重復代碼，核心的任務代碼很少。當然，如果僅僅是找位點的CHROM，POS，ID，REF，ALT，QUAL這幾個參數時，這樣做也可以。因為vcf格式規范，這幾個參數的結構相對簡單。但是如果處理頭文件信息，或者處理INFO，FORMAT參數時，要寫比較復雜的正則表達式，這樣做不僅繁瑣，而且容易出錯。

Python的PyVCF庫解決了這個問題，它通過正則表達式把vcf文件信息轉換成結構化的信息，簡化了vcf文件的處理過程，方便後續提取相關參數及處理。

cmd界面

或者從https://github.com/jamescasbon/PyVCF網站上下載安裝包，自行安裝。

PyVCF庫的名字為vcf，導入之後可以使用其方法對vcf文件做處理。

調用vcf.Reader類處理vcf文件，vcf文件信息就被保存到vcf_reader中了。它是一個可迭代對象，它的迭代元素都是一個_Record對象的實例，保存著非注釋行的一行信息，即變異位點的具體信息。通過它，我們可以很輕易地得到位點的詳細信息。

_Record是vcf.model中的一個對象，除了它還有_Call，_AltRecord等對象。它的基本屬性為CHROM，POS，ID，REF，ALT，QUAL，FILTER，INFO，FORMAT，也就是vcf中的一行位點信息。接下來對這些屬性一一說明：

CHROM：染色體名稱，類型為str。

POS：位點在染色體上的位置，類型為int。

ID：一般是突變的rs號，類型為str。如果是‘.’，則為None。

REF：參考基因組在該位點上的鹼基，類型為str。

ALT：在該位點的測序結果。是_AltRecord類的子類實例的列表。類型為list。_AltRecord類有4個子類，代表了突變的幾種類型：如snp，indel，structual variants等。所有的實例都可以進行比較（僅限於相等的比較，沒有大於小於之說），部分子類沒有實現str方法，也就是說不能轉成字符串。

QUAL：該位點的測序質量，類型為int或float。

FILTER：過濾信息。將FILTER列按分號分隔形成的字符串列表，類型為list。如果未給出參數則為None。

INFO：該位點的一些測試指標。將‘=’前的參數作為鍵，後面的參數作為值，構建成的字典。類型為dict。

FORMAT：基因型信息。保存vcf的FORMAT列的原始形式，類型為str。

除了這些基本屬性之外，_Record對象還有一些其他屬性：

samples：把FORMAT信息作為鍵，後面對應的信息做為值，構建成的字典（CallData對象），以及sample名稱，這兩個值組成一個Call對象，共同構成samples的一個元素。這樣就把sample和基因型信息給關聯起來，按下標訪問每一個Call對象。samples類型為list。

start：突變開始的位置

end：突變結束的位置

alleles：該位點所有的可能情況，由REF和ALT參數組成的列表（REF類型是str，ALT參數是_AltRecord對象的子類實例），類型是list。

_Record對象方法：

對象之間比較大小方法：根據染色體名稱和位置信息比較。Python3只實現了‘=’和‘<’的比較。

迭代方法：對samples裡的元素進行迭代。

字符串方法：只返回CHROM，POS，REF，ALT四列信息。

genotype(name)方法，和samples按下標訪問不同，這個方法提供按sample名稱進行訪問的功能。

add_format(fmt)， add_filter(flt)， add_info(info, value=True)：給相應的屬性添加元素。

get_hom_refs()：拿到samples中該位點未突變的所有sample，返回列表。

get_hom_alts()：拿到samples中該位點100%突變的所有sample，返回列表。

get_hets()：拿到samples中該位點基因型為雜合的所有sample，返回列表。

get_unknown()：拿到samples中該位點基因型未知的所有sample，返回列表。

_Record對象還有很多有用的方法屬性：

num_called：該位點已識別的sample數目。

call_rate：已識別的sample數目占sample總數的比例。

num_hom_ref，num_hom_alt，num_het，num_unknown：四種基因型的數量

aaf：所有sample等位基因的頻率（即除開REF），返回列表。

heterozygosity：該位點的雜合度，0.5為雜合突變，0為純合突變。

var_type：突變類型，包括‘snp’，‘indel’，‘sv’（structural variant），‘unknown’。

var_subtype：更加細化的突變類型，如‘indel’包括‘del’，‘ins’，‘unknown’。

is_snp，is_indel，is_sv，is_transition，is_deletion：判斷突變是不是snp，indel，sv，transition，deletion等等。

在讀vcf文件時，總共有六個參數可供選擇，如上圖所示。

fsock：目標文件的文件對象，可以用open(文件名)得到這個文件對象。

filename：文件名，當fsock和filename同時存在時，優先考慮fsock。

compressed：是否要解壓，不提供參數時由程序自行判斷（以文件名是否以.gz結尾判斷是否需要解壓）。
prepend_chr：在保存染色體名稱時，是否加前綴‘chr’，默認不加，如果vcf文件的染色體名稱本來沒有前綴‘chr’，可設置為True，自動加上。

strict_whitespace：是否嚴格以制表符‘\t’分隔數據。True則表示嚴格按制表符分，False表示可以夾雜空格。

encoding：文件編碼。

頭文件信息主要保存在Reader對象的屬性中，包括alts，contigs，filters，formats，infos，metadata。

alts使用實例：

其他的屬性用法類似。

Reader對象實現了兩個方法：

next()：獲得下一行的數據，也就是返回下一個_Record對象。可以顯式調用next()得到下一行數據，也可以直接迭代Reader對象，它會自動調用next()函數以獲得下一行數據。

fetch(chrom，start=None，end=None)：返回chrom染色體從start+1到end坐標的所有突變位點。不給end，就返回chrom染色體從start+1到末尾的所有突變位點；

start和end都不給，就返回chrom染色體所有的突變位點。這個方法需要用另一個第三方Python模塊pysam來建立文件索引，如果沒有安裝這個模塊，會導致錯誤。

另外，使用這個方法之後，它會將對象的可迭代范圍改成fetch()得到的突變位點，所以用這個方法，原來的迭代進度就失效了。

這個庫還有一個Writer對象，在此就不詳細介紹了，因為大部分對vcf文件的處理都可以用上面兩個對象的知識搞定。

這個庫實現的所有功能，都可以自己寫代碼實現，而且實現方法比較簡單。之所以要用這個庫來處理vcf文件，是因為這個庫考慮的東西可能比我們自己了解的更多，其實現也可能比我們自己的代碼更加完備合理。

還有，重復造車總歸是不好的。

以上就是python PyVCF文件處理VCF文件格式實例詳解的詳細內容，更多關於python PyVCF處理VCF格式的資料請關注軟件開發網其它相關文章！