本發明涉及數據處理
技術領域：
，特別是指一種基因變異檢測方法及裝置。
背景技術：
：基因組變異檢測，這里指的是從二代測序數據的比對結果中，找出與參考基因組不同的堿基或序列片段，即單堿基變異(SNV)和插入缺失變異(INDEL)。目前被廣泛應用的10基因型模型只考慮了單堿基變異類型，插入缺失變異一般要單獨檢測，這使得現有模型的基因變異檢測不夠簡便。技術實現要素：有鑒于此，本發明的目的在于提出一種能夠同時檢測單堿基變異和插入缺失變異的基因變異檢測方法及裝置。基于上述目的本發明提供的基因變異檢測方法，包括：從基因比對結果中統計每個位點的比對信息；考慮堿基變異和插入缺失變異，創建16基因型模型；使用所述16基因型模型搜索候選變異位點；使用隨機森林對候選變異位點進行分類與篩選，并輸出篩選后的候選變異結果。在一些可選實施方式中，所述從基因比對結果中統計每個位點的比對信息，具體包括以下比對信息：堿基類型和對應的每個堿基類型的比對質量值、等位基因型及其Reads支持數量、正負鏈數量、插入缺失數量及插入序列信息，和/或，軟剪切位點數量。在一些可選實施方式中，所述考慮堿基變異和插入缺失變異，創建16基因型模型，具體包括：假設樣品是一個二倍體生物樣品，堿基類型有ATCG四種，則二倍體基因型的統計類型有{AA,AC,AG,AT,CC,CG,CT,GG,GT,TT,AX,CX,GX,TX,XX,XY},其中X和Y分別代表有最多比對reads支持和第二多reads支持的插入或缺失。在一些可選實施方式中，所述使用所述16基因型模型搜索候選變異位點，具體包括：通過貝葉斯模型計算出每個位點最大可能的基因型；將所述最大可能的基因型與參考基因組的對應位點的參考信息進行比較，得到所述候選變異位點。在一些可選實施方式中，所述使用隨機森林對候選變異位點進行分類與篩選，并輸出篩選后的候選變異結果，具體包括：定義真實變異位點和偽變異位點；建立隨機森林模型；經過隨機森林模型從所述候選變異位點中篩選得到更加可信的候選變異位點；將所述更加可信的候選變異位點以VCF格式輸出，并且直接應用于下游的分析工具。本發明的另一方面，提供了一種基因變異檢測裝置，包括：統計模塊，用于從基因比對結果中統計每個位點的比對信息；模型創建模塊，用于考慮堿基變異和插入缺失變異，創建16基因型模型；搜索模塊，用于使用所述16基因型模型搜索候選變異位點；分類與篩選模塊，用于使用隨機森林對候選變異位點進行分類與篩選，并輸出篩選后的候選變異結果。在一些可選實施方式中，所述從基因比對結果中統計每個位點的比對信息，具體包括以下比對信息：堿基類型和對應的每個堿基類型的比對質量值、等位基因型及其Reads支持數量、正負鏈數量、插入缺失數量及插入序列信息，和/或，軟剪切位點數量。在一些可選實施方式中，所述模型創建模塊，具體用于：假設樣品是一個二倍體生物樣品，堿基類型有ATCG四種，則二倍體基因型的統計類型有{AA,AC,AG,AT,CC,CG,CT,GG,GT,TT,AX,CX,GX,TX,XX,XY},其中X和Y分別代表有最多比對reads支持和第二多reads支持的插入或缺失。在一些可選實施方式中，所述搜索模塊，具體用于：通過貝葉斯模型計算出每個位點最大可能的基因型；將所述最大可能的基因型與參考基因組的對應位點的參考信息進行比較，得到所述候選變異位點。在一些可選實施方式中，所述分類與篩選模塊，具體用于：定義真實變異位點和偽變異位點；建立隨機森林模型；經過隨機森林模型從所述候選變異位點中篩選得到更加可信的候選變異位點；將所述更加可信的候選變異位點以VCF格式輸出，并且直接應用于下游的分析工具。從上面所述可以看出，本發明提供的基因變異檢測方法及裝置，通過考慮堿基變異和插入缺失變異，創建了16基因型模型，使得整體計算更加方便而且大幅提高了準確性和靈敏度；同時，利用隨機森林對檢測結果進行修正，使得檢測結果更加精確。附圖說明圖1為本發明提供的基因變異檢測方法的一個實施例的流程示意圖；圖2為本發明提供的基因變異檢測裝置的一個實施例的模塊結構示意圖；圖3為采用本發明提供的基因變異檢測方法及裝置實施例后得到的預測正確率與真實比例的對比示意圖。具體實施方式為使本發明的目的、技術方案和優點更加清楚明白，以下結合具體實施例，并參照附圖，對本發明進一步詳細說明。需要說明的是，本發明實施例中所有使用“第一”和“第二”的表述均是為了區分兩個相同名稱非相同的實體或者非相同的參量，可見“第一”“第二”僅為了表述的方便，不應理解為對本發明實施例的限定，后續實施例對此不再一一說明。基于上述目的，本發明實施例的第一個方面，提供了一種能夠同時檢測單堿基變異和插入缺失變異的基因變異檢測方法的實施例。如圖1所示，為本發明提供的基因變異檢測方法的一個實施例的流程示意圖。所述基因變異檢測方法，包括：步驟101：從基因比對結果中統計每個位點的比對信息。這里，基因比對結果可以通過任意的基因比對軟件的比對處理而得到，具體比對過程不再贅述。在一些可選實施方式中，所述步驟101——從基因比對結果中統計每個位點的比對信息，具體包括以下比對信息：堿基類型和對應的每個堿基類型的比對質量值、等位基因型及其Reads支持數量、正負鏈數量、插入缺失數量及插入序列信息，和/或，軟剪切位點數量。具體地，測序數據經過質量值校正(scorerecalibration)、序列比對(alignment)、去重復(de-duplication)和重比對(realignment)等一系列處理后，需要收集每個位點的一套詳細的統計信息以用于變異檢測分析。每個位點的統計信息如下表1：表1位點統計信息對于等位基因型及其Reads支持數量(加權和)的統計信息：在一個成功比對read中(Reads，讀長，是高通量測序中獲得的測序序列，每一個read是一段堿基序列)，每個堿基都會包含一個重校準的質量值，且質量值范圍為0到40之間。為了儲存堿基的質量值，我們為不同的質量值范圍分配相應的權重，如下表2所示：表2堿基質量值參數權重0–10[0–Weight0]011–13(Weight0–Weight1]114–17(Weight1–Weight2]218–20(Weight2–Weight3]321–40(Weight3–40]4表1中，為了將堿基質量值轉化為權重值所設定的范圍值參數，所述參數一列跟前面的堿基質量值一列的范圍是相互對應的，這里的weight0\1\2\3,分別對為10、13、17、20。每個成功配對的堿基使相應的等位基因型增加一個權重計數，如一個質量值為25的堿基A，其相應的等位基因型計數加4，若其質量值為5則計數加0。對于正負鏈數量的統計信息：依據比對結果，每個成功比對的堿基使對應的等位基因的正鏈或者負鏈計數加一。與權重計數不同，無論堿基的重校準質量值是多少，這里都增加一個計數。例如一個堿基的被多條堿基質量值小于10的reads覆蓋，它的權重計數為零而正負鏈計數則確切地反映了成功比對的reads的條數。對于插入缺失數量及插入序列信息的統計信息：若比對結果中存在插入缺失，其信息將會被記錄下來，格式為‘mI’或‘nD’，其中m和n分別表示插入和缺失的片段長度。除了不同類型的插入缺失的數量，插入的片段信息也會存儲到一個動態分配的數據結構中且高質量與低質量片段信息分別記錄在兩個計數器中。對于軟剪切位點數量的統計信息：如果比對結果中出現軟剪切位點，其數量將會被同時記錄下來。軟剪切的方向也會被記錄下來以區分頭端剪切和末端剪切。步驟102：考慮堿基變異和插入缺失變異，創建16基因型模型。對于每個位點而言，我們需要依據S1中收集的比對信息來推測該位點的真實基因型并與參考基因組作比較，從而找出那些發生變異的位點，即候選變異位點。為了實現對一個位點真實基因型的推測，首先我們需要構建相應的基因型模型。因此，在一些可選實施方式中，所述步驟102——考慮堿基變異和插入缺失變異，創建16基因型模型，具體包括以下步驟：假設樣品是一個二倍體生物樣品，堿基類型有ATCG四種，則二倍體基因型的統計類型有{AA,AC,AG,AT,CC,CG,CT,GG,GT,TT,AX,CX,GX,TX,XX,XY},其中X和Y分別代表有最多比對reads支持和第二多reads支持的插入或缺失(reads支持越多，可信度越高)。與被廣泛應用的10基因模型不同，這里提出的16基因型模型在二倍體的背景中同時考慮了堿基變異和插入缺失變異，該16基因型模型統一了A,C,G,T和INDEL(插入缺失),這個統一的模型不僅使得計算方便而且大幅提高了準確性和靈敏度。步驟103：使用所述16基因型模型搜索候選變異位點。在一些可選實施方式中，所述步驟103——使用所述16基因型模型搜索候選變異位點，具體包括以下步驟：通過貝葉斯模型計算出每個位點最大可能的基因型；將所述最大可能的基因型與參考基因組的對應位點的參考信息進行比較，得到所述候選變異位點。具體地，16基因型的后驗概率的計算，使用了貝葉斯模型：P(G|F)∝P(F|G)P(G)其中，F表示觀察到的{A,C,T,G,X,Y}各自的加權計數(weightedcount)，P(G)表示某種基因型G的先驗概率，P(F|G)表示的是基因型為G時觀察到的F的概率，P(G|F)表示的是觀察到F的基因型G的概率。一般有如下幾個原因導致我們觀察到某個位置的堿基跟參考基因組上的不一樣：測序錯誤(badbasecallorprimaryanalysis)，比對錯誤(badalignment)，基因變異(variantallele)。一般質量值校正，可以一定程度修正第1類錯誤(即測序錯誤)。這里，我們設置兩種錯誤概率：PS表示單堿基等位基因概率，PID表示插入缺失等位基因概率。一般經驗，PS設置會大于PID。如果一個錯誤(測序錯誤或比對錯誤)發生，假設：1){A,C,G,T}每種堿基被觀察到的概率相同，為PS；2){X，Y}每個被觀察到的概率相同，為PID。定義錯誤率為：Perr＝mPs+nPID其中，m為基因型G中的單個堿基{A,T,C,G}的數量，n為基因型G中{X,Y}的數量。默認的設置：PS＝0.01PID＝0.005當我們觀察到純合的基因型時，我們會期望觀察到接近100％的純合位點。當觀察到雜合的位點時，我們期望觀察到50％的兩個等位基因。為了檢測觀察到的reads覆蓋深度分布與預期匹配的好壞，我們使用雙尾費舍爾精確檢驗(Two-tailedFisher’sExactTest(FET))來檢測，計算公式如下：計算的p-value會當作某種基因型G的概率。[p-value越小表示可能性越大]。具體計算P(F|G)的過程如下：當觀察到加權計數F＝{FA,FC,FG,FT,FX,FY}，一個純合基因型G＝AA的概率的計算，表示如下：P(F|AA,Perr)＝Phom(FA)·Pe(FC,FG,FT,FX,FY)一個雜合基因型G＝CG的概率計算，表示如下：P(F|CG,Perr)＝Phet(FC,FG)·Pe(FA,FT,FX,FY)其中，Phom為觀察到純合基因型的概率：Phet為觀察到雜合基因型的概率：Pe為觀察到基因型G以外的等位基因：定義：θ表示兩個不相關單倍體單個堿基不同的頻率，ω表示兩個不相關單倍體單個插入卻是不同的頻率，ε表示轉換顛換比(Ti/Tv)。先驗概率可表示如下表3：表3默認值：θ＝0.001ω＝0.0001ε＝2.1最終輸出的基因型Gmax，為有最大后驗概率的基因型：Gmax＝argmax{P(G|F,Perr)}。至此，我們通過貝葉斯模型計算出每個位點最大可能的基因型Gmax，將這個基因型與參考基因組該位點的參考信息作比較，就能初步地得到我們想要的候選變異位點。而這些搜索到的候選變異位點還需要進一步的篩選，去除一些假陽性的變異位點，我們將在下一步使用隨機森林的模型來實現。步驟104：使用隨機森林對候選變異位點進行分類與篩選，并輸出篩選后的候選變異結果。在一些可選實施方式中，所述步驟104——使用隨機森林對候選變異位點進行分類與篩選，并輸出篩選后的候選變異結果，具體包括以下步驟：定義真實變異位點和偽變異位點；建立隨機森林模型；經過隨機森林模型從所述候選變異位點中篩選得到更加可信的候選變異位點；將所述更加可信的候選變異位點以VCF格式輸出，并且直接應用于下游的分析工具。具體地，變異分類的目的是為了給每一個檢測出來的候選變異一個更加精確的預測正確率(Probabilityofa“truesite”)，并基于這一正確率的估計值篩選出一個高準確率的變異位點的集合；這里的預測正確率，可參考表6中的預測正確率，是模型經過計算后給出預測正確的概率，模型的使用者依據這個概率來判斷一個候選變異位點是否是真實的。隨機森林是一種常用的機器學習的分類方法，我們的變異位點分類即是利用隨機森林模型來對候選變異是真實的遺傳變異(geneticvariant)而非測序及分析導致的人為誤差(artifact)的概率和變異指標之間的關系做一個連續的共變的估計，模型基于的分類依據如下：1)真實變異位點(truesites)，一般來說這些位點在SNP(SingleNucleotidePolymorphisms，單核苷酸的多態性)數據庫(如dbSNPv129,HapMap3,Omni2.5MSNPchiparrayandMills,1000Ggoldstandardindels)中呈現多態性。2)偽變異位點(falsesites)，每個候選變異位點，若5個用于偽變異位點篩選的參數指標(Strandbias；Readpositionbias；Totaldepth；Leftaveragebasequality；Rightaveragebasequality)中有3個以上落在最差的5％內，則這個位點被歸為偽變異位點。5％指的是這個候選變異在所有的由上一步貝葉斯模型檢測出的候選變異位點中落在最差的5％。參照表5，對于鏈偏差，鏈偏差值取值范圍(0,1],最差的5％則指的是鏈偏差值最小的5％的候選變異位點；對于Read位置偏差，位置偏差取值范圍[-1,1],最差的5％為絕對值最大的5％；對于各等位基因深度總和，總測序深度越深越好，越少的read支持數，可信度越差，最差5％指的是深度最少的5％；對于位點左、右側堿基平均質量值，堿基質量值取值范圍為[0,40]，值越大越好，越小越差，也就越不可信，最差5％指的是質量值最小的5％。之后，這一自適應性誤差模型就可用于變異檢測出的候選變異位點真實性的概率估算。模型訓練使用的特性如表4所示。表4模型訓練所用到的特性用于挑選偽變異位點的特性如表5所示。表5用于挑選偽變異位點的特性模型訓練詳情及結果：應用上文步驟102中介紹的16基因型模型從NA82178樣本50×150bp的雙端測序數據中搜索單堿基變異和插入缺失變異位點，再利用SNP數據庫(dbSNPv137,IndelDB,1000GandMills)從這些候選變異位點中挑選真實變異位點。這樣，我們總共得到1,813,021個“truesites”和31,588個“falsesites”。我們使用31,588個“falsesites”和26,501個隨機選取的“truesites”組成58,089個位點的訓練集合。用這個訓練集合建立了一個有96棵決策樹的隨機森林模型。模型的可靠性分析如下表6所示：表6其中，Probabilityofa“truesite”為隨機森林模型給出的變異候選位點的預測正確率，預測正確率也就是模型經過計算后給出的預測得到的正確概率，即候選變異位點為真實變異位點(truesite)的概率。模型的使用者依據這個概率來判斷一個候選變異位點是否是真實可靠的。“比例”為訓練集中“truesites”所占的真實比例，預測正確率與真實比例的對比如圖3所示。從表6和圖3可以看出我們的隨機森林模型預測的正確率與“truesites”所占的真實的比例非常接近，可以說明我們的模型可以有效的區分候選變異位點是否為真實變異位點。經過第三步的候選變異分類，我們進一步篩選出了更加可信的候選變異位點。最終的候選變異位點將以VCF(VariantCallingFile)的格式輸出，并且可以直接應用于下游的分析工具(如snpEff,VEP,GATK)和在線數據庫(如Ingenuity,GenomeTrax)。其中，所述輸出結構中還可以包括每一個變異的質量值，每一個變異的質量值計算公式如下：其中Popt(G|F)是最大的后驗概率，PsubOpt(G|F)是第二大后驗概率。一般來說，質量值q越大，這一位點的最大概率基因型的不確定性越小，Gmax也就越可信。從上述實施例可以看出，本發明實施例提供的基因變異檢測方法，通過考慮堿基變異和插入缺失變異，創建了16基因型模型，使得整體計算更加方便而且大幅提高了準確性和靈敏度；同時，利用隨機森林對檢測結果進行修正，使得檢測結果更加精確。本發明實施例的第二個方面，提供了一種基因變異檢測裝置的實施例。如圖2所示，為本發明提供的基因變異檢測裝置的一個實施例的模塊結構示意圖。所述基因變異檢測裝置，包括：統計模塊201，用于從基因比對結果中統計每個位點的比對信息；模型創建模塊202，用于考慮堿基變異和插入缺失變異，創建16基因型模型；搜索模塊203，用于使用所述16基因型模型搜索候選變異位點；分類與篩選模塊204，用于使用隨機森林對候選變異位點進行分類與篩選，并輸出篩選后的候選變異結果。從上述實施例可以看出，本發明實施例提供的基因變異檢測裝置，通過考慮堿基變異和插入缺失變異，創建了16基因型模型，使得整體計算更加方便而且大幅提高了準確性和靈敏度；同時，利用隨機森林對檢測結果進行修正，使得檢測結果更加精確。在一些可選實施方式中，所述從基因比對結果中統計每個位點的比對信息，具體包括以下比對信息：堿基類型和對應的每個堿基類型的比對質量值、等位基因型及其Reads支持數量、正負鏈數量、插入缺失數量及插入序列信息，和/或，軟剪切位點數量。在一些可選實施方式中，所述模型創建模塊202，具體用于：假設樣品是一個二倍體生物樣品，堿基類型有ATCG四種，則二倍體基因型的統計類型有{AA,AC,AG,AT,CC,CG,CT,GG,GT,TT,AX,CX,GX,TX,XX,XY},其中X和Y分別代表有最多比對reads支持和第二多reads支持的插入或缺失。在一些可選實施方式中，所述搜索模塊203，具體用于：通過貝葉斯模型計算出每個位點最大可能的基因型；將所述最大可能的基因型與參考基因組的對應位點的參考信息進行比較，得到所述候選變異位點。在一些可選實施方式中，所述分類與篩選模塊204，具體用于：定義真實變異位點和偽變異位點；建立隨機森林模型；經過隨機森林模型從所述候選變異位點中篩選得到更加可信的候選變異位點；將所述更加可信的候選變異位點以VCF格式輸出，并且直接應用于下游的分析工具。需要特別指出的是，上述裝置的實施例僅采用了所述方法的實施例來具體說明各模塊的工作過程，本領域技術人員能夠很容易想到，將這些模塊應用到所述方法的其他實施例中。當然，由于所述方法實施例中的各個步驟均可以適當地進行相互交叉、替換、增加、刪減，因此，這些合理的排列組合變換之于所述裝置也應當屬于本發明的保護范圍，并且不應將本發明的保護范圍局限在所述實施例之上。所屬領域的普通技術人員應當理解：以上任何實施例的討論僅為示例性的，并非旨在暗示本公開的范圍(包括權利要求)被限于這些例子；在本發明的思路下，以上實施例或者不同實施例中的技術特征之間也可以進行組合，步驟可以以任意順序實現，并存在如上所述的本發明的不同方面的許多其它變化，為了簡明它們沒有在細節中提供。另外，為簡化說明和討論，并且為了不會使本發明難以理解，在所提供的附圖中可以示出或可以不示出與集成電路(IC)芯片和其它部件的公知的電源/接地連接。此外，可以以框圖的形式示出裝置，以便避免使本發明難以理解，并且這也考慮了以下事實，即關于這些框圖裝置的實施方式的細節是高度取決于將要實施本發明的平臺的(即，這些細節應當完全處于本領域技術人員的理解范圍內)。在闡述了具體細節(例如，電路)以描述本發明的示例性實施例的情況下，對本領域技術人員來說顯而易見的是，可以在沒有這些具體細節的情況下或者這些具體細節有變化的情況下實施本發明。因此，這些描述應被認為是說明性的而不是限制性的。盡管已經結合了本發明的具體實施例對本發明進行了描述，但是根據前面的描述，這些實施例的很多替換、修改和變型對本領域普通技術人員來說將是顯而易見的。例如，其它存儲器架構(例如，動態RAM(DRAM))可以使用所討論的實施例。本發明的實施例旨在涵蓋落入所附權利要求的寬泛范圍之內的所有這樣的替換、修改和變型。因此，凡在本發明的精神和原則之內，所做的任何省略、修改、等同替換、改進等，均應包含在本發明的保護范圍之內。當前第1頁1 2 3