本發(fā)明涉及生物信息學(xué)領(lǐng)域，更特別地，涉及一種估算四倍體物種基因組的二倍化程度的方法。

背景技術(shù)：

四倍體基因組是正常二倍體基因組在自然情況或人工操作過(guò)程中，經(jīng)自我加倍產(chǎn)生；或者兩個(gè)親緣關(guān)系相對(duì)較近物種通過(guò)自然或人工雜交并染色體加倍產(chǎn)生。前者往往稱之為同源四倍體，后者稱之為異源四倍體。四倍體物種的細(xì)胞核內(nèi)包含組染色體，可以組成兩套二倍體基因組。在進(jìn)化過(guò)程中，四倍體基因組中的兩套二倍體基因組會(huì)逐漸趨于形成更多的差異，向兩個(gè)方向發(fā)展，使得同源四倍體異源化最終成為二倍體基因組。這個(gè)過(guò)程可以稱之為四倍體基因組的二倍化。

當(dāng)前基因組學(xué)研究中，并沒(méi)有專門針對(duì)四倍體基因組二倍化程度進(jìn)行直接分析的方法。在對(duì)四倍體基因組進(jìn)行組裝分析時(shí)，將四倍體基因組中兩套基因組分別或同時(shí)組裝出來(lái)，通過(guò)組裝出來(lái)的contigs或染色體水平的基因組內(nèi)互相比較估算兩套基因組間的差異。如果兩套基因組間差異相對(duì)較大，處于高度二倍化狀態(tài)時(shí)，這種方法可以得到相對(duì)較好的結(jié)果。但由于四倍體基因組的特性，如果兩套基因組間差異相對(duì)較小，那么往往不能得到理想的組裝效果，不能分開兩套基因組，因此這套方法對(duì)于同源四倍體效果較差。如果采用將兩套基因組分別組裝的方案，首先需要找到兩套基因組對(duì)應(yīng)的親本物種，隨后分別組裝。由于異源四倍體是由兩親本物種雜交后又歷經(jīng)進(jìn)化過(guò)程產(chǎn)生，親本物種的基因組與異源四倍體中的基因組存在差異，這套方法會(huì)有一定程度的誤差。

因此，需要一種直接估算四倍體基因組的二倍化程度并進(jìn)行量化的方法。

技術(shù)實(shí)現(xiàn)要素：

為解決以上問(wèn)題，本發(fā)明提供了一種估算四倍體物種基因組的二倍化程度的方法，其特征在于，包括以下步驟：

s1：對(duì)所述四倍體物種的基因組進(jìn)行二代測(cè)序，得到四倍體基因組測(cè)序數(shù)據(jù)；

s2：將所述四倍體基因組測(cè)序數(shù)據(jù)與二倍體基因組測(cè)序數(shù)據(jù)進(jìn)行比較，估算所述四倍體物種基因組的二倍化程度，所述二倍體基因組測(cè)序數(shù)據(jù)為所述四倍體物種的近緣二倍體物種基因組的二代測(cè)序數(shù)據(jù)。

本發(fā)明以分析四倍體基因組二倍化程度為目標(biāo)的直接分析方法，首次提出了二倍化率的概念用于量化四倍體基因組的二倍化程度，不依賴于目標(biāo)的基因組序列，具有成本低，速度快，成功率高等優(yōu)勢(shì)。

在一個(gè)實(shí)施方案中，s1和s2中所述二代測(cè)序?yàn)閕llumina測(cè)序。

在另一個(gè)實(shí)施方案中，所述二倍體基因組測(cè)序數(shù)據(jù)通過(guò)測(cè)序得到，或?yàn)橐延械臏y(cè)序數(shù)據(jù)。

在另一個(gè)實(shí)施方案中，所述四倍體基因組測(cè)序數(shù)據(jù)的測(cè)序深度不小于100x。

在另一個(gè)實(shí)施方案中，所述二倍體基因組測(cè)序數(shù)據(jù)的測(cè)序深度不小于30x。

在另一個(gè)實(shí)施方案中，s2包括：

s21：獲得所述二倍體基因組測(cè)序數(shù)據(jù)；

s22：對(duì)所述四倍體基因組測(cè)序數(shù)據(jù)和所述二倍體基因組測(cè)序數(shù)據(jù)進(jìn)行分析處理，分別得到四倍體基因組k-mer集合和二倍體基因組k-mer集合；

s23：分別統(tǒng)計(jì)備所述四倍體基因組k-mer集合和所述二倍體基因組k-mer集合中的k-mer總數(shù)，并以k-mer的出現(xiàn)頻數(shù)為橫坐標(biāo)，k-mer的種類數(shù)縱坐標(biāo)分別制備所述四倍體基因組k-mer集合的k-mer種類數(shù)頻數(shù)分布圖和所述二倍體基因組k-mer集合的k-mer種類數(shù)頻數(shù)分布圖(例如分別統(tǒng)計(jì)在四倍體基因組和二倍體基因組中出現(xiàn)頻率在1到1000次的k-mer的種類數(shù))，并以所述k-mer種類數(shù)頻數(shù)分布圖中的第一波谷前k-mer為錯(cuò)誤k-mer；

s24：根據(jù)所述四倍體基因組k-mer集合的k-mer種類數(shù)頻數(shù)分布圖和所述二倍體基因組k-mer集合的k-mer種類數(shù)頻數(shù)分布圖分別計(jì)算所述四倍體基因組的序列重復(fù)率和雜合度，以及所述二倍體基因組的序列重復(fù)率；

s25：根據(jù)所述四倍體基因組的序列重復(fù)率和雜合度，以及所述二倍體基因組的序列重復(fù)率計(jì)算所述四倍體基因組的二倍化率，計(jì)算公式如下：

公式ii：

d：四倍體基因組二倍化率

α：四倍體基因組序列重復(fù)率

β：二倍體基因組序列重復(fù)率

k：四倍體基因組雜合度。

基于基因組k-mer的頻數(shù)以及種類數(shù)的分布估算基因組特征。該方法僅需要一定覆蓋度的基本二代高通量數(shù)據(jù)就可以完成，且該過(guò)程對(duì)測(cè)序文庫(kù)的種類及插入片段類型無(wú)要求且不需要進(jìn)行基因組組裝，所以受基因組復(fù)雜度影響極小。相比較，通過(guò)基因組序列估算方法受限于基因組序列，尤其是高度復(fù)雜的同源四倍體基因組序列的組裝目前是一個(gè)世界性難題。往往需要構(gòu)建多種測(cè)序文庫(kù)，采用各種方法測(cè)序，并制定復(fù)雜的組裝策略，而往往得不到質(zhì)量足夠好的基因組序列進(jìn)行后續(xù)分析。

在另一個(gè)實(shí)施方案中，s24通過(guò)以下方法計(jì)算所述四倍體基因組的序列重復(fù)率和所述二倍體基因組的序列重復(fù)率：在所述四倍體基因組k-mer集合的k-mer種類數(shù)頻數(shù)分布圖中，以第一雜合峰2x處為主峰位置，以主峰1.8x處為界限，出現(xiàn)頻率大于該界限的k-mer為所述四倍體基因組的重復(fù)k-mer；所述二倍體基因組k-mer集合的k-mer種類數(shù)頻數(shù)分布圖中，以主峰后1.8x處為界限，出現(xiàn)頻率大于該界限的k-mers為重復(fù)k-mer，并根據(jù)公式i分別計(jì)算所述四倍體基因組的序列重復(fù)頻率和所述二倍體基因組的序列重復(fù)頻率

公式i：

r：基因組序列重復(fù)率

nkspecies：非重復(fù)k-mer種類數(shù)

nkfrequency：非重復(fù)k-mer頻數(shù)

ekmer：錯(cuò)誤k-mer數(shù)

akmer：總k-mer數(shù)；

并以所述四倍體基因組k-mer集合的k-mer種類數(shù)頻數(shù)分布圖中第一雜合峰計(jì)算其基因組雜合度?；蚪M雜合度的計(jì)算方法是現(xiàn)有技術(shù)，在本發(fā)明中本著突出重點(diǎn)的原則不做贅述。

在另一個(gè)實(shí)施方案中，s22中，通過(guò)以下方法處理所述四倍體基因組測(cè)序數(shù)據(jù)和所述二倍體基因組測(cè)序數(shù)據(jù)：

s221：過(guò)濾掉所述四倍體基因組測(cè)序數(shù)據(jù)和所述二倍體基因組測(cè)序數(shù)據(jù)中低質(zhì)量堿基和/或短于一定長(zhǎng)度的讀序；

s222：將所述四倍體基因組測(cè)序數(shù)據(jù)和所述二倍體基因組測(cè)序數(shù)據(jù)分化成k-mer，分別得到四倍體基因組k-mer集合和二倍體基因組k-mer集合。

進(jìn)一步地，s221中，所述低質(zhì)量堿基讀序?yàn)樾蛄袃啥速|(zhì)量值小于20的讀序，所述短于一定長(zhǎng)度的讀序?yàn)樾蛄锌傞L(zhǎng)小于50的讀序。

附圖說(shuō)明

圖1為本發(fā)明方法的流程圖。

具體實(shí)施方式

以下結(jié)合實(shí)例對(duì)本發(fā)明的原理和特征進(jìn)行描述，所舉實(shí)例只用于解釋本發(fā)明，并非用于限定本發(fā)明的范圍。

本發(fā)明的方法的流程示意圖如圖1所示。

以某種魚類(因論文未發(fā)表，為保密起見，暫時(shí)不公開其種屬)為例，該魚類在自然種群中包括四倍體基因組類型和二倍體基因組類型。其四倍體基因組預(yù)估大小約2.4g，二倍體基因組預(yù)估大小約1.2g。我們基于本方法估算該四倍體基因組物種的兩套二倍體基因組間是否發(fā)生了分化以及分化程度有多少。具體實(shí)施過(guò)程如下：

1)分別對(duì)兩個(gè)物種進(jìn)行建庫(kù)測(cè)序，均建立插入片段300-350bp的illuminahiseq文庫(kù)并進(jìn)行pe150高通量測(cè)序。四倍體物種總共測(cè)得約280g數(shù)據(jù)，測(cè)序深度約117x，二倍體物種總共測(cè)得約52g數(shù)據(jù)，測(cè)序深度約43x。

2)使用htqc軟件的ht-trim和ht-filter模塊采用默認(rèn)參數(shù)分別對(duì)兩組數(shù)據(jù)進(jìn)行堿基質(zhì)量過(guò)濾和讀序過(guò)濾，約過(guò)濾掉0.02％的數(shù)據(jù)，整體測(cè)序深度不變。

3)使用jellyfishcount首先以k-mer＝17計(jì)算兩組數(shù)據(jù)的所有k-mer類型及頻數(shù)；使用jellyfishstats統(tǒng)計(jì)并獲得兩數(shù)據(jù)的k-mer總數(shù)，四倍體數(shù)據(jù)有250,217,368,293個(gè)k-mer，二倍體數(shù)據(jù)有46,513,565,383個(gè)k-mer；使用jellyfishhisto繪制以出現(xiàn)頻率為橫坐標(biāo)，k-mer的種類數(shù)為縱坐標(biāo)統(tǒng)計(jì)k-mer種類數(shù)頻數(shù)分布。

4)二倍體k-mer種類數(shù)頻數(shù)分布中主峰位于k-mer頻數(shù)為40處，因此頻數(shù)大于40×1.8＝72的k-mer為重復(fù)k-mer，重復(fù)k-mer有2,3765,979,213個(gè)。其第一波谷為頻數(shù)＝7處，因此錯(cuò)誤k-mer有976,614,636個(gè)。通過(guò)公式(1)計(jì)算得其基因組序列重復(fù)約52％。

5)四倍體k-mer種類數(shù)頻數(shù)分布中第一雜合峰位于k-mer頻數(shù)46處，因此頻數(shù)大于46×2×1.8≈166的k-mer為重復(fù)k-mer，重復(fù)k-mer有202,733,541,319個(gè)。其第一波谷為頻數(shù)＝11處，因此錯(cuò)誤k-mer有7,879,758,786個(gè)。通過(guò)公式計(jì)算的其基因組序列重復(fù)約84％。同過(guò)其第一雜合峰處的k-mer數(shù)計(jì)算可得雜合度約1.1％。

6)同過(guò)公式(2)代入以上計(jì)算得到的兩基因組重復(fù)序列含量以及四倍體基因組雜合度可得該四倍體物種基因組二倍化率約32％。

以上所述僅為本發(fā)明的較佳實(shí)施例，并不用以限制本發(fā)明，凡在本發(fā)明的精神和原則之內(nèi)，所作的任何修改、等同替換、改進(jìn)等，均應(yīng)包含在本發(fā)明的保護(hù)范圍之內(nèi)。