鄧繼忠， 林偉森， 甘四明， 黃華盛,3， 李　梅， 金　濟， 何明昊

(1 華南農(nóng)業(yè)大學 工程學院，廣東 廣州 510642; 2 中國林業(yè)科學研究院 熱帶林業(yè)研究所，廣東 廣州 510520; 3 廣東科技學院，廣東 東莞 523083)

?

一種二倍體片段測序中SNP檢測系統(tǒng)的構建

鄧繼忠1， 林偉森1， 甘四明2， 黃華盛1,3， 李梅2， 金濟1， 何明昊1

(1 華南農(nóng)業(yè)大學 工程學院，廣東 廣州 510642; 2 中國林業(yè)科學研究院 熱帶林業(yè)研究所，廣東 廣州 510520; 3 廣東科技學院，廣東 東莞 523083)

摘要：【目的】開發(fā)基于模式識別方法的二倍體片段測序中單核苷酸多態(tài)性(Single nucleotide polymorphism，SNP)自動檢測系統(tǒng)，提高檢測的準確性。【方法】采用LabWindows/CVI 9.0開發(fā)平臺，結合Matlab函數(shù)庫編程，以二倍體PCR片段測序的.ab1或.scf格式文件作為源數(shù)據(jù)，首先分離出堿基G、A、T和C，進行一維離散小波濾波，再對各堿基處的波形進行典型特征提取，最后運用基于反向傳播神經(jīng)網(wǎng)絡的分類器完成SNP識別和判斷?！窘Y果】系統(tǒng)界面友好、運行穩(wěn)定。SNP等級分為6級，允許用戶對可疑的SNP進行人工修正，對尾葉桉Eucalyptus urophylla的26個測序序列143個SNP的測試中檢測準確率、假陽性率和假陰性率均明顯優(yōu)于之前的類似軟件?！窘Y論】本文所構建的SNP自動檢測系統(tǒng)準確性高，不需參考序列，可用于二倍體PCR片段測序中SNP的高效檢測。

關鍵詞：二倍體; 測序; SNP檢測; 模式識別; 尾葉桉

單核苷酸多態(tài)性(Single nucleotide polymorphism，SNP)是指生物體基因組中單個堿基的顛換或轉換導致的DNA變異，是最普遍和最廣泛的基因組變異形式，如擬南芥Arabidopsisthaliana3個品系間存在超過82萬個SNP，平均約1.4 kb就存在1個SNP[1]。SNP目前已廣泛應用于生物、醫(yī)學和農(nóng)學等研究領域[2-3]。

傳統(tǒng)的PCR片段測序中，測序廠家所提供的軟件只能識別各波峰位置的單個堿基，而對二倍體個體內(nèi)存在SNP的雙峰處的低峰堿基不能有效檢測。已有一些第三方軟件可用于二倍體物種PCR片段測序的SNP自動檢測，如Mutation Surveyor (http://www.softgenetics.com/MutationSurveyor.html)、novoSNP[4]和PolyPhred[5]等。但是，這些軟件均需參考序列，具有局限性，不能有效用于一些序列的分析，如表達序列標簽(Expressed sequence tag，EST)和候選基因的測序中的內(nèi)含子區(qū)域，而且操作上較為繁瑣。這些軟件的檢測準確性均偏低[6]。此外，PolyPhred[5]還需要至少8個測序文件以保證足夠的準確性，不利于少量樣品的檢測。針對二倍體PCR片段測序中如何進行SNP的自動、有效檢測，本文基于模式識別的方法構建了無需參考序列進行SNP自動檢測的計算機系統(tǒng)，并驗證了系統(tǒng)的準確性。

1材料與方法

1.1系統(tǒng)開發(fā)平臺

本研究基于LabWindows/CVI (以下簡稱CVI)平臺，以CVI 9.0作為系統(tǒng)前臺進行主界面的設計和基本功能菜單的實現(xiàn)。CVI是基于C語言的虛擬儀器軟件開發(fā)平臺，其功能強大，具有靈活的交互式編程方法和豐富的用戶接口資源[7-8]，且運行速度快、界面控件豐富。

采用MATLAB 2012函數(shù)庫編程實現(xiàn)模式識別的分類器構建。MATLAB具有編程簡單、仿真能力強和易于擴展移植等優(yōu)點[9]，提供內(nèi)部神經(jīng)網(wǎng)絡工具箱，且可外接用戶自行開發(fā)的模式識別算法或軟件包，有利于用戶進行各種模式識別算法的測試與建模。本文主要通過MATLAB實現(xiàn)模式識別的分類器結構的建立，利用反向傳播神經(jīng)網(wǎng)絡(Back-propagation neural network, BPNN)進行權值和閾值的訓練，再將訓練好的BPNN結構移植至CVI環(huán)境用作SNP分類器。

系統(tǒng)的運行環(huán)境為Windows XP Professional、Windows 7或Window 8。

1.2系統(tǒng)功能模塊和界面的設計

系統(tǒng)主要包括測序數(shù)據(jù)導入、堿基分離、噪聲處理、SNP識別、數(shù)據(jù)顯示與存儲、人工校正等功能模塊。測序數(shù)據(jù)導入可打開擴展名為.ab1或.scf的測序文件。堿基分離是根據(jù)測序數(shù)據(jù)中G、A、T和C 4種堿基的不同熒光顏色而將它們分別提取出來。噪聲處理通過一維離散小波變換濾除噪聲，便于后續(xù)的數(shù)據(jù)特征提取。SNP識別主要根據(jù)測序峰圖的波形特征，結合二倍體內(nèi)SNP表現(xiàn)為雙峰的現(xiàn)象，提取各堿基位置的典型波形特征作為BPNN的輸入向量，輸出SNP屬性分數(shù)，再參照判別標準給出SNP等級。數(shù)據(jù)顯示可在運行窗口直觀顯示測序峰圖和序列以及SNP的位置、雙峰堿基、屬性分數(shù)和等級，識別結果可存為.txt文件。人工校正允許用戶對誤判和漏判的SNP進行刪除、添加和更改堿基等手動操作。系統(tǒng)的工作流程見圖1。

圖1　SNP檢測流程圖

系統(tǒng)界面包括主操作區(qū)、原始峰圖顯示區(qū)、SNP識別結果區(qū)和人工校正區(qū)等4個功能區(qū)域(圖2)。

主操作區(qū)位于系統(tǒng)界面的上部，如圖2所示，包含面板選項卡、菜單欄、快捷鍵和信息提示欄。面板選項卡可提供多個檢測面板，各檢測面板的序列分析數(shù)據(jù)相互獨立，以便進行多個測序文件的檢測。菜單欄包括文件、編輯、濾波選擇、識別方法和幫助等一級菜單，各一級菜單可下拉顯示二級菜單。文件的二級菜單包括原始測序文件導入、SNP檢測結果保存等。編輯的二級菜單包括X軸和Y軸縮放以及隱藏/顯示堿基G、A、T和C。濾波選擇的二級菜單包括小波變換等方法。識別方法的二級菜單包括BPNN等多種方法?？旖萱I工具欄列出常用功能，包括文件打開、SNP檢測(按默認的濾波和模式識別方法)、增加面板、關閉當前面板、峰圖和結果顯示窗口X軸和Y軸縮放等。信息提示欄顯示文件路徑等。

原始峰圖顯示區(qū)位于界面右中區(qū)域，見圖2，不同顏色代表堿基G、A、T和C顯示測序的峰圖，可通過主操作區(qū)的快捷按鈕選擇性地顯示一種或幾種堿基類型，也可進行X軸和Y軸的縮放。

SNP識別結果區(qū)位于界面下部區(qū)域，見圖2，界面右下區(qū)域顯示SNP識別結果，包括堿基位置、堿基類型、雙峰SNP、SNP分數(shù)和SNP等級等。界面左下區(qū)域列表總結SNP檢測的結果，點擊某一行可與右下窗口相應位置的SNP進行快速定位。

人工校正區(qū)位于界面左中區(qū)域(圖2)，允許用戶對任一位置堿基進行更改、刪除和增加等操作。

圖2　SNP檢測系統(tǒng)的界面

1.3系統(tǒng)主要功能的實現(xiàn)

1.3.1測序數(shù)據(jù)的導入和堿基分離本系統(tǒng)可導入.ab1或.scf格式的測序文件，根據(jù)文件內(nèi)容分離出各種堿基。以 .ab1文件為例說明。按照美國Applied Biosystems公司的ab1f文件格式， .ab1文件是由一個指定的文件夾存儲文件信息，包括文件名稱、堿基類型和堿基數(shù)量等。文件7 ～ 34 字節(jié)是‘DirEntry’結構體，指向文件夾所在的位置，其結構為：

Struct DirEntry{ sInt32 name;sInt32 number;sInt16 elementtype;sInt16 elementsize;sInt32 numelements; sInt32 datasize; sInt32 dataoffset;sInt32 datahandle;}，

其中，number顯示通道編號(對應4種堿基)，elementsize顯示堿基所占空間大小，numelements顯示堿基數(shù)量， dataoffset顯示堿基位置的偏移量。根據(jù)這些信息即可將各通道的堿基和熒光數(shù)據(jù)提取出來。

1.3.2濾波的實現(xiàn)分離出的堿基熒光信號通常含有噪聲，影響SNP分析的準確性，可通過小波分解進行噪聲濾除。小波分解具有良好的時頻分析特性，在時域和頻域內(nèi)都具有突顯信號局部特征的能力，并且能夠在多種分辨率下觀測數(shù)據(jù)，是噪聲濾除的有效方法[10-11]。

本系統(tǒng)調(diào)用了CVI中Signal Processing Toolkit 7.0.2工具箱所含的小波工具包的函數(shù)，一維離散小波的濾波過程如下：分離4個通道的堿基數(shù)據(jù);讀取各通道的堿基數(shù)據(jù)長度，獲得濾波參數(shù);調(diào)用SptAllocCoeffWFBD函數(shù)，分配小波濾波工具包的濾波參數(shù)結構;調(diào)用SptReadCoeffWFBD函數(shù)，讀取小波工具包參數(shù);調(diào)用SptDiscreteWaveletTransform (filterData, dataLenght, analysisFilter, SptSymmetric, NULL, NULL, scales, shift, dwt, &outputSize, lengths)函數(shù)，進行一維離散小波分解，濾除后的數(shù)據(jù)存儲于變量dwt中，用于后續(xù)分析。

1.3.3SNP識別分類器的構建利用BPNN構建SNP識別分類器。BPNN是一種反向傳播且能修正誤差的多層映射模式，由輸入層、隱含層和輸出層及各層神經(jīng)元之間連接組成[12]。SNP識別可以視作從波形特征輸入到SNP類型輸出的非線性映射問題。

基于文獻[13]，選擇了SNP處雙峰的波峰距離、高度比值和起伏度比值作為特征參數(shù)，歸一化處理后作為BPNN輸入層的特征向量，即輸入層的神經(jīng)元數(shù)為3。BPNN輸出層為SNP類別的診斷結果，可能是雙峰SNP位點或單峰非SNP位點，即輸出層神經(jīng)元為1;輸出層中設定范圍0 ～ 100為1個位點的屬性分數(shù)，并結合屬性分數(shù)和周邊雜峰將SNP可能性分為6個等級，等級越小越可能為SNP。

MATLAB神經(jīng)網(wǎng)絡工具箱中調(diào)用BPNN構建SNP分類器的過程具體如下：調(diào)用net=newff[input, output, 10, (tansig,purelin),trainlm]，初始建立適于SNP識別的BPNN分類器;調(diào)用net.trainParam.goal=0.01，設置誤差，當網(wǎng)絡精度達到誤差則停止訓練，訓練好的權值和閾值存于MATLAB中net結構體的net.IW、net.LW和net.b元素;調(diào)用net=train (net, input, output)，CVI中重新塑造BPNN即可用于SNP檢測。

1.4數(shù)據(jù)處理

因SNP檢測可能存在假陽性和假陰性，人工校正的堿基刪除、添加和更改的結果可直接顯示在當前面板。SNP檢測結果的保存格式為.txt文件。

利用尾葉桉Eucalyptusurophylla26個測序文件的143個SNP進行系統(tǒng)測試，并與Mutation Surveyor (http://www.softgenetics.com/MutationSurveyor.html)、novoSNP[4]和PolyPhred[5]的結果進行比較，分析軟件的有效性。所用26個測序文件為尾葉桉SNP開發(fā)的EST重測序文件[14]。

2結果與分析

2.1系統(tǒng)主要功能的實現(xiàn)

系統(tǒng)導入測序文件后能正確地進行堿基判讀，判讀結果與測序結果較為一致。圖3顯示了本系統(tǒng)對1個序列的判讀與測序結果一致。但對測序質(zhì)量較差或者雜峰較多的區(qū)域，測序結果亦不準確，與本系統(tǒng)的判讀結果可能有異。圖3B表明本系統(tǒng)能可靠地識別雙峰堿基位置的SNP，155 bp處的SNP判讀為T/A，分值為87，等級為1。并且，其他較低的雜峰被有效排除，未被誤判為SNP。

b：本系統(tǒng)判讀的序列，155 bp處雙峰、SNP為T/A

人工校正可對特定位置的第1和/或第2堿基進行修改、刪除和添加操作。圖4顯示了對174 bp處誤讀的‘T’修改為‘G’，表明系統(tǒng)可有效地進行人工校正。

數(shù)據(jù)存儲的.txt文件見圖5，第1行為文件名，數(shù)字代表顯示堿基位置，識別序列上面一行為堿基序列，下面一行顯示識別的SNP的第2堿基。測序前20個堿基一般是較雜亂的峰，SNP可靠性低。

b：校正后174 bp處更正為‘G’

圖5　SNP檢測結果存為.txt文件

2.2SNP檢測的有效性

對于測試的尾葉桉26個測序文件，對比另外3個軟件，本系統(tǒng)具有較高的準確性和誤檢率。表1顯示了不同軟件在高和低的標準下對143個SNP的識別準確率、漏判(假陰性)率和誤判(假陽性)率，不管是在高還是低的標準下，本文所開發(fā)的系統(tǒng)均有最高的SNP判讀準確率以及最低的假陰性率和假陽性率，PolyPhred除外，因其未檢測到任何SNP。

表1不同軟件識別SNP的結果對比1)

Tab.1Comparison of software performance in SNP detection

判別標準SNP準確性DiSNPindelnovoSNPMutationSurveyorPolyPhredSNP數(shù)量判別率/%SNP數(shù)量判別率/%SNP數(shù)量判別率/%SNP數(shù)量判別率/%高準確13896.5128.43021.000漏判53.513191.611379.0143100誤判4524.61758.63150.800低準確14198.63423.83121.700漏判21.410976.211478.3143100誤判7542.98972.43150.000

1) DiSNPindel、novoSNP、Mutation Surveyor和PolyPhred高的標準分別為等級2、分值13、中敏感性和等級2[6]，低的標準分別為等級4、分值6、高敏感性和等級4準確SNP的判別率=準確SNP數(shù)量/(準確SNP數(shù)量+漏判SNP數(shù)量)/×100%,漏判SNP的百分比=漏判SNP數(shù)量/(準確SNP數(shù)量+漏判SNP數(shù)量)/×100%,誤判SNP的百分比=誤判SNP數(shù)量/(誤判SNP數(shù)量+準確SNP數(shù)量)×100%。

3結論與討論

本文基于CVI 9.0開發(fā)平臺，結合MATLAB函數(shù)庫，運用特征提取和模式識別等技術構建了二倍體片段測序的SNP自動檢測系統(tǒng)。本系統(tǒng)開發(fā)中結合了CVI運行速度快和用戶接口資源豐富以及MATLAB函數(shù)庫多樣和用戶可自行擴充函數(shù)等優(yōu)點，降低了系統(tǒng)開發(fā)的難度，并且檢測系統(tǒng)界面友好、易于擴充。系統(tǒng)各模塊功能的完整實現(xiàn)和測試的有效性也證明相關技術的合理運用。

高通量的新一代測序技術已經(jīng)興起，但傳統(tǒng)的PCR片段測序仍非常重要[15]。測序廠家提供的軟件只能識別各序列位置的最高峰所對應的堿基，因此雙峰位置的低峰堿基的識別需要第三方軟件。另外，測序也常出現(xiàn)雜峰，且雜峰可能高于堿基峰而導致測序廠家的軟件不能正確讀序。本文針對這一問題研發(fā)了SNP自動檢測系統(tǒng)，可以有效實現(xiàn)雙峰判讀和雜峰濾除等功能，為二倍體PCR片段測序的SNP判讀提供了有力工具。

與功能類似的軟件Mutation Surveyor (http://www.softgenetics.com/MutationSurveyor.html)、novoSNP[4]和PolyPhred[5]相比，本系統(tǒng)不需參考序列，既為EST和候選基因的測序中內(nèi)含子無參考序列的SNP檢測提供了便利，也減少了輸入?yún)⒖夹蛄械姆爆嵅僮?并且，本文所開發(fā)的檢測系統(tǒng)具有最高的SNP判讀準確率以及最低的假陰性率和假陽性率，可用于二倍體PCR片段測序中SNP的高效檢測。

實際應用中，Mutation Surveyor和novoSNP可能假陽性率極高，如對15個基因171個SNP的檢測中的假陽性SNP分別為3 728和505個，顯著多于正確判讀的SNP[16]。另一方面，任何檢測系統(tǒng)的準確率都嚴重依賴于測序質(zhì)量，實際應用時需要優(yōu)化測序體系[4]。

參考文獻：

[1]OSSOWSKI S, SCHNEEBERGER K, CLARK R M, et al. Sequencing of natural strains ofArabidopsisthalianawit short reads[J]. Genome Res， 2008，18(12)：2024-2033．

[2]唐立瓊，肖層林，王偉平. SNP分子標記的研究及其應用進展[J]. 中國農(nóng)學通報， 2012，28(12)：154-158．

[3]許家磊，王宇，后猛，等. SNP檢測方法的研究進展[J]. 分子植物育種， 2015，13(2)：475-482．

[4]WECKX S，DEL-FAVERO J，RADEMAKERS R，et al. novoSNP, a novel computational tool for sequence variation discovery[J]. Genome Res， 2005，15(3)：436-442．

[5]MATTHEW S，JAMES S，ROBERTSON P D，et al. Automating sequence-based detection and genotyping of single nucleotide polymorphisms (SNPs) from diploid samples[J]. Nat Genet， 2006，38(3)：375-381．

[6]DENG J Z，HUANG H S，YU X L，et al. DiSNPindel: Improved intra-individual SNP and InDel detection in direct amplicon sequencing of a diploid[J]. BMC Bioinformatics， 2015，16：343．

[7]仇志平，李樹軍. LabWindows/CVI虛擬儀器軟件在測試領域中的應用[J]. 計算機工程與設計，2007，28(22)：5544-5548．

[8]劉君華. 虛擬程序編程語言LabWindows/CVI編程[M]. 北京：電子工業(yè)出版社，2001．

[9]肖偉，劉忠，曾新勇，等. MATLAB 程序設計與應用[M]. 北京：清華大學出版社，2005．

[10]BUI T D，CHEN G. Translation-invariant denoising using multiwavelets[J]. IEEE Trans Sig Proc，1998，46(12)：3414-3420．

[11]PAN Q，ZHANG P，DAI G，et al. Two denoising methods by wavelet transform[J]. IEEE Trans Sig Proc，1999，47(12)：3401-3406．

[12]MCKEOWN J J，STELLA F，HALL G. Some numerical aspects of the training problem for feed-forward neural nets[J].Neural Netw，1997，10(9)：1455-1463.

[13]黃華盛. 基于模式識別的二倍體個體內(nèi)SNP和InDel自動檢測[D]. 廣州：華南農(nóng)業(yè)大學，2014．

[14]YU X，GUO Y，ZHANG X，et al. Integration of EST-CAPS markers into genetic maps of andEucalgptusurophyllaandE.tereticornisand their alignment withE.grandisgenome sequence[J]. Silvae Genet， 2012，61(6)：247-255．

[15]STUDER A，ZHAO Q，ROSS-IBARRA J，et al. Identification of a functional transposon insertion in the maize domestication genetb1[J]. Nat Genet，2011，43(11)：1160-1163．

[16]NGAMPHIW C，KULAWONGANUNCHAI S，ASSAWAMAKIN A，et al. VarDetect: A nucleotide sequence variation exploratory tool[J]. BMC Bioinformatics，2008，9(S12)：9．

【責任編輯霍歡】

Development of an automatic system for SNP detection in diploid fragment sequencing

DENG Jizhong1, LIN Weisen1, GAN Siming2, HUANG Huasheng1,3, LI Mei2, JIN Ji1, HE Minghao1

(1 College of Engineering, South China Agricultural University, Guangzhou 510642, China;2 Research Institute of Tropical Forestry, Chinese Academy of Forestry, Guangzhou 510520, China;3 Guangdong University of Science and Technology, Dongguan 523083, China)

Abstract：【Objective】 This study aims to develop a pattern-recognition based system for automatic single nucleotide polymorphism (SNP)detection in diploid fragment sequencing and improve the detection accuracy. 【Method】The LabWindows/CVI 9.0 platform and Matlab environment were combined for analyzing .ab1 or .scf files generated in diploid PCR fragment sequencing. Firstly， four bases G, A, T and C were separated for eliminating noise through one-dimensional discrete wavelet filtering, following with extraction of typical features of each base position (peak) from a fluorescence curve. A classifier based on back-propagation neural network was then used for SNP recognition and diagnosis. 【Result】 This established system was characterized by friendly interface, stable operation and manual modification accessibility. It classified the SNP reliability into six grades. Performance test with 143 SNPs of 26 sequencing fragments from Eucalyptus urophylla demonstrated that our system outperformed three previously reported software packages in detecting accuracy, false positive and false negative rates. 【Conclusion】 Our system has a high rate of accuracy without the need for a reference sequence. It could be used for efficient SNP detection in diploid PCR fragment sequencing.

Key words：diploid; sequencing; SNP detection;pattern recognition; Eucalyptus urophylla

中圖分類號：TP391.4; Q523.8

文獻標志碼：A

文章編號：1001- 411X(2016)03- 0115- 06

基金項目：863計劃(2013AA102705);國家自然科學基金(31270702， 31070592)

作者簡介：鄧繼忠(1963—)，男，副教授，博士，E-mail: jz-deng@scau.edu.cn；通信作者：甘四明(1970—)，男，研究員，博士，E-mail：Siming_Gan@126.com

收稿日期：2015- 08- 26優(yōu)先出版時間：2016-04-15

優(yōu)先出版網(wǎng)址：http://www.cnki.net/kcms/detail/44.1110.s.20160415.1554.002.html

鄧繼忠， 林偉森， 甘四明，等.一種二倍體片段測序中SNP檢測系統(tǒng)的構建[J].華南農(nóng)業(yè)大學學報，2016,37(3):115- 120.