黃體浩，李俊青，趙海勇

聊城大學(xué)計(jì)算機(jī)學(xué)院，山東 聊城 252000

拷貝數(shù)變異（copy number variation，CNV）是指DNA拷貝數(shù)偏離正?？截悢?shù)時(shí)基因組序列發(fā)生擴(kuò)增或缺失，是基因組結(jié)構(gòu)變異的一種重要形式。生物醫(yī)學(xué)研究估計(jì)，人類(lèi)基因組中包含1 000 多個(gè)CNV 區(qū)域，約占基因序列的1%。CNV 與腫瘤疾病的產(chǎn)生和發(fā)展有著極其重要的聯(lián)系，也是引發(fā)眾多復(fù)雜腫瘤疾病的主要原因之一，因此發(fā)現(xiàn)原癌基因和抑癌基因尤為關(guān)鍵的一步就是準(zhǔn)確檢測(cè)目標(biāo)基因組的拷貝數(shù)變異。

利用下一代測(cè)序（next-generation sequencing，NGS）技術(shù)檢測(cè)拷貝數(shù)變異的策略主要分為四種：基于讀深（read depth，RD）、基于雙端比對(duì)、基于分裂讀段、基于序列組裝。其中最常用的拷貝數(shù)檢測(cè)方法是基于RD的方法，Xie等[1]通過(guò)利用病例樣本與對(duì)照樣本之間的讀段數(shù)目（read count，RC）比例建立統(tǒng)計(jì)模型，并根據(jù)比例差異識(shí)別CNV；Yoon等[2]采用基于事件測(cè)試的方法在個(gè)體基因組中檢測(cè)到的缺失和擴(kuò)增跨多個(gè)基因組進(jìn)行檢查，以確定個(gè)體之間的多態(tài)性；Ivakhno 等[3]采用離散小波變換對(duì)QC 過(guò)濾的RC 進(jìn)行去噪，并使用隱馬爾可夫模型（hidden markov model，HMM）進(jìn)行分段。早期基于RD的拷貝數(shù)檢測(cè)方法過(guò)度依賴RD信號(hào)，Abyzov等[4]采用均值漂移算法對(duì)單個(gè)樣本的RD信號(hào)進(jìn)行處理，并以此統(tǒng)計(jì)模型建立檢測(cè)CNV；Xi等[5]將最小化貝葉斯信息準(zhǔn)則算法擴(kuò)展應(yīng)用到多個(gè)樣本中，以檢測(cè)重復(fù)出現(xiàn)的CNV；Miller等[6]建立服從負(fù)二項(xiàng)分布的RC統(tǒng)計(jì)模型以減少均方根誤差；Boeva 等[7]對(duì)GC 含量進(jìn)行建模，并將模型用于歸一化RD 信號(hào)；Gusnanto 等[8]通過(guò)損失函數(shù)糾正RC比例以消除GC含量偏差的影響。然而上述方法只能檢測(cè)出顯著性的CNV，且不能準(zhǔn)確定位CNV 的斷點(diǎn)位置。反向傳播（back propagation，BP）神經(jīng)網(wǎng)絡(luò)擁有強(qiáng)大的非線性映射能力，可用來(lái)解決拷貝數(shù)變異檢測(cè)中多個(gè)特征的聯(lián)合作用，然而B(niǎo)P 神經(jīng)網(wǎng)絡(luò)在實(shí)際應(yīng)用中仍存在諸多缺點(diǎn)，如局部最優(yōu)、收斂速度慢、過(guò)擬合等。

基于上述研究，本文提出CNV-GABP（copy number variation detection of BP neural network based on genetic algorithm）檢測(cè)方法。方法分為四個(gè)部分：特征值計(jì)算與歸一化、BP神經(jīng)網(wǎng)絡(luò)結(jié)構(gòu)確定、遺傳算法優(yōu)化和BP神經(jīng)網(wǎng)絡(luò)預(yù)測(cè)。該方法的優(yōu)點(diǎn)在于：（1）充分考慮了多個(gè)特征之間的聯(lián)合作用；（2）在BP 神經(jīng)網(wǎng)絡(luò)基礎(chǔ)上，使用遺傳算法（genetic algorithm，GA）優(yōu)化BP 神經(jīng)網(wǎng)絡(luò)模型，提高模型的泛化能力；（3）優(yōu)化后的BP 神經(jīng)網(wǎng)絡(luò)用以解決多個(gè)特征之間的聯(lián)合作用，進(jìn)一步提高了檢測(cè)準(zhǔn)確率和精度，可有效識(shí)別特征不明顯的CNV。

1 相關(guān)工作

近年來(lái)，研究者綜合考慮RD信號(hào)和其他因素，提出了很多CNV檢測(cè)算法，Wang等[9-11]采用基于HMM的方法，分別利用期望最大化（expectation-maximization，EM）算法和置信傳播算法估計(jì)模型的最優(yōu)參數(shù)；Onsongo等[12]采用隨機(jī)森林算法結(jié)合機(jī)器學(xué)習(xí)識(shí)別CNV；Johansson 等[13]選擇與分析樣本最具相似度的對(duì)照樣本，歸一化后分析兩個(gè)樣本間的平均覆蓋率差異；Dharanipragada等[14]使用覆蓋深度方法檢測(cè)CNV；Chen等[15]采用最大懲罰似然估計(jì)的方法識(shí)別CNV和檢測(cè)CNV邊界；Melcher等[16]和Jugas 等[17]分別采用局部多項(xiàng)式回歸擬合算法和中位數(shù)方法進(jìn)行GC偏差矯正，然后利用全變分模型和循環(huán)二元分割算法進(jìn)行分割讀??；Yuan等[18]采用基于統(tǒng)計(jì)混合模型的貝葉斯方法來(lái)推斷拷貝數(shù)的擴(kuò)增和缺失；Xiao 等[19]采用基于高斯混合模型的聚類(lèi)策略消除假陽(yáng)性，并在聚類(lèi)步驟中引入EM算法優(yōu)化模型參數(shù)；Li等[20]建立高斯混合模型，并用先大窗后小窗的分割方式進(jìn)行分割；Yuan等[21]采用孤立森林算法將RD信號(hào)生成孤立樹(shù)集合，并使用全變分模型進(jìn)行去噪；Li等[22]在RD信號(hào)的基礎(chǔ)上引入了連續(xù)位置的相關(guān)性作為統(tǒng)計(jì)量以提高對(duì)CNV 的識(shí)別能力；Yuan 等[23]為每個(gè)分段分配一個(gè)離群點(diǎn)因子，采用局部異常因子算法識(shí)別CNV；Yuan等[24]利用腫瘤純度和RD信號(hào)建立非線性模型，采用窮舉搜索策略計(jì)算最優(yōu)解；Zhao 等[25]利用BP 神經(jīng)網(wǎng)絡(luò)構(gòu)建多特征訓(xùn)練算法以預(yù)測(cè)CNV；Yuan 等[26]利用潛變量模型重建觀察到的信號(hào)，并對(duì)其進(jìn)行分層處理。上述方法不僅僅依賴RD信號(hào)，減少了由于讀段數(shù)目分布不均和定位模糊產(chǎn)生的影響，提高了靈敏度和特異性，但是仍然忽略了多個(gè)特征之間聯(lián)合作用的影響，并且由于測(cè)序誤差、GC 含量偏差等影響，一些特征較不明顯的CNV 無(wú)法被識(shí)別。

2 CNV-GABP檢測(cè)方法

本文提出的CNV-GABP方法是一種基于多個(gè)特征的檢測(cè)方法，用于從NGS 數(shù)據(jù)中準(zhǔn)確檢測(cè)CNV。目前很多方法也逐漸考慮基因組位置之間的相關(guān)性，例如文獻(xiàn)[4]和文獻(xiàn)[18]都將基因組位置之間的相關(guān)性綜合到統(tǒng)計(jì)量中來(lái)提高RD信號(hào)的檢測(cè)效率，但是這些方法忽略了特征值之間的相互作用。CNV-GABP將GC含量、RD 信號(hào)、基因組相鄰位置之間的相關(guān)性以及比對(duì)質(zhì)量這四個(gè)特征融入到基因組序列的分析中，并通過(guò)訓(xùn)練BP 神經(jīng)網(wǎng)絡(luò)使四個(gè)特征相互作用。為了克服BP 神經(jīng)網(wǎng)絡(luò)全局搜索能力差、易陷入局部最優(yōu)等缺點(diǎn)，本文使用遺傳算法優(yōu)化BP 神經(jīng)網(wǎng)絡(luò)模型，增強(qiáng)模型的泛化能力，提高方法的檢測(cè)性能。

2.1 算法流程

CNV-GABP 方法工作流程如圖1 所示。首先輸入一個(gè)參考基因組和一個(gè)測(cè)序樣本，然后對(duì)輸入數(shù)據(jù)進(jìn)行預(yù)處理，最后經(jīng)過(guò)四個(gè)步驟實(shí)現(xiàn)CNV 的預(yù)測(cè)：（1）定義與CNV相關(guān)的四個(gè)特征值，計(jì)算四個(gè)特征值并歸一化；（2）利用四個(gè)特征值構(gòu)建基于遺傳算法優(yōu)化的BP 神經(jīng)網(wǎng)絡(luò)；（3）訓(xùn)練BP神經(jīng)網(wǎng)絡(luò)，其中標(biāo)記為CNV的數(shù)據(jù)從仿真數(shù)據(jù)和真實(shí)數(shù)據(jù)中采樣獲取；（4）根據(jù)訓(xùn)練完成的BP 神經(jīng)網(wǎng)絡(luò)模型預(yù)測(cè)每個(gè)窗的CNV 狀態(tài)（正常、擴(kuò)增或缺失）。

2.2 數(shù)據(jù)預(yù)處理與特征值量化

CNV-GABP使用比對(duì)算法BWA來(lái)處理格式為fasta的參考文件和格式為fastq 的序列樣本文件，以生成格式為bam 的比對(duì)文件。使用samtools 軟件從bam 文件中獲取RC 文件，RD 信號(hào)是根據(jù)RC 文件計(jì)算得來(lái)的，為方便計(jì)算RD 信號(hào)，CNV-GABP 將提取的RC 文件分成連續(xù)非重疊且大小相等的窗，可根據(jù)需求自由設(shè)置窗的大小，本文中窗的大小為1 000 bp，取每個(gè)窗內(nèi)的平均RC值作為該窗的RD值。由于GC含量偏差對(duì)RD信號(hào)產(chǎn)生噪聲，通常如文獻(xiàn)[21]中的方法會(huì)預(yù)先對(duì)GC 含量偏差進(jìn)行校正，以期降低GC含量偏差的影響，這種方法適用于一些特定的情況，但在大多數(shù)應(yīng)用場(chǎng)景下并未取得十分顯著的效果。因此CNV-GABP并未同常規(guī)方法般預(yù)先對(duì)每個(gè)窗的GC 含量偏差進(jìn)行校正，而是將GC含量作為該窗的特征值之一?；蚪M相鄰位置之間的相關(guān)性和比對(duì)質(zhì)量在一定程度上反應(yīng)出該區(qū)域的CNV狀態(tài)，因此CNV-GABP 也將這兩個(gè)因素作為特征值。其中基因組相鄰位置之間的相關(guān)性值參考文獻(xiàn)[25]，計(jì)算方式如式（1）：

其中Ci表示基因組第i個(gè)位置與相鄰位置的相關(guān)性值，Ri表示表示第i個(gè)窗的RD 值，n表示窗的總數(shù)，且j

每個(gè)窗的量化方式如式（2）：

其中，Gi表示GC 含量，Qi表示比對(duì)質(zhì)量，Ci和Ri的定義由式（1）給出。

由于四個(gè)特征值的表示標(biāo)準(zhǔn)不同，為消除特征值之間的量綱影響，CNV-GABP 對(duì)特征值進(jìn)行歸一化處理。以比對(duì)質(zhì)量為例，歸一化方式如式（3）：

2.3 BP神經(jīng)網(wǎng)絡(luò)構(gòu)建

BP神經(jīng)網(wǎng)絡(luò)是一種按照誤差逆向傳播算法訓(xùn)練的多層前饋式神經(jīng)網(wǎng)絡(luò)，以其優(yōu)越的非線性映射能力和柔性的網(wǎng)絡(luò)結(jié)構(gòu)被廣泛應(yīng)用于分類(lèi)問(wèn)題，BP 神經(jīng)網(wǎng)絡(luò)算法的學(xué)習(xí)過(guò)程包括信號(hào)的正向傳播與誤差的反向回傳。正向傳播時(shí)，數(shù)據(jù)從輸入層傳入，經(jīng)隱含層處理，傳向輸出層，若輸出層輸出與期望輸出不一致，則將誤差作為調(diào)整信號(hào)反向回傳，調(diào)整神經(jīng)元之間的權(quán)值矩陣，減小誤差。本文用歸一化之后的特征值組成的四元組Xi=(Ri,Gi,Qi,Ci)作為輸入構(gòu)建BP 神經(jīng)網(wǎng)絡(luò)。如圖2表示該網(wǎng)絡(luò)的拓?fù)浣Y(jié)構(gòu)。

圖2 BP神經(jīng)網(wǎng)絡(luò)拓?fù)浣Y(jié)構(gòu)Fig.2 Topology of BP neural network

其中ω表示輸入層與隱含層、隱含層與輸出層之間的權(quán)重，θ表示偏差，X表示神經(jīng)元輸入向量，即歸一化之后的特征值組成的四元組。除此之外，學(xué)習(xí)速率與收斂性密切相關(guān)，是神經(jīng)網(wǎng)絡(luò)的重要參數(shù)之一，本文使用BP神經(jīng)網(wǎng)絡(luò)的默認(rèn)值0.1作為該網(wǎng)絡(luò)的學(xué)習(xí)速率。

2.4 基于遺傳算法的BP神經(jīng)網(wǎng)絡(luò)

遺傳算法是一種啟發(fā)式隨機(jī)搜索算法，能夠在不需要誤差函數(shù)梯度信息的情況下學(xué)習(xí)近似最優(yōu)解，具有良好的全局搜索能力和隱并行性，適合復(fù)雜的非線性問(wèn)題求解。遺傳算法的主要步驟包括編碼、種群初始化、選擇、交叉和變異。使用遺傳算法優(yōu)化BP 神經(jīng)網(wǎng)絡(luò)主要是對(duì)神經(jīng)網(wǎng)絡(luò)各項(xiàng)參數(shù)的調(diào)整，本文使用遺傳算法優(yōu)化BP神經(jīng)網(wǎng)絡(luò)的權(quán)值和閾值。

基于遺傳算法優(yōu)化BP神經(jīng)網(wǎng)絡(luò)算法的基本思想是將BP神經(jīng)網(wǎng)絡(luò)算法與遺傳算法相結(jié)合，當(dāng)BP算法訓(xùn)練網(wǎng)絡(luò)的收斂速度較慢時(shí)，將BP 神經(jīng)網(wǎng)絡(luò)每個(gè)隱含層節(jié)點(diǎn)的閾值和權(quán)重作為遺傳算法的輸入信息，并對(duì)它們進(jìn)行編碼以生成染色體，然后使用遺傳算法的選擇算子、交叉算子和變異算子來(lái)生成新的后代作為BP算法的初始值，最后繼續(xù)使用BP 算法訓(xùn)練網(wǎng)絡(luò)?；谶z傳算法優(yōu)化的BP神經(jīng)網(wǎng)絡(luò)算法流程如圖3所示。

圖3 基于遺傳算法的BP神經(jīng)網(wǎng)絡(luò)程序流程圖Fig.3 Flow chart of BP neural network based on genetic algorithm

算法分為兩個(gè)部分：遺傳算法優(yōu)化部分和BP 神經(jīng)網(wǎng)絡(luò)部分。

（1）遺傳算法優(yōu)化部分。首先對(duì)種群進(jìn)行初始化，初始化種群由BP 神經(jīng)網(wǎng)絡(luò)的初始權(quán)值和閾值組成，隨后執(zhí)行選擇、交叉和變異操作，將種群中執(zhí)行完操作的個(gè)體代入BP神經(jīng)網(wǎng)絡(luò)進(jìn)行迭代以尋找最優(yōu)個(gè)體。具體內(nèi)容如下：

步驟1 初始化種群。本文采用實(shí)數(shù)編碼，種群中的每個(gè)個(gè)體都是由輸入層與隱含層之間的權(quán)值、隱含層閾值、隱含層與輸出層之間的權(quán)值、輸出層閾值組成的實(shí)數(shù)串。

步驟2 計(jì)算適應(yīng)度函數(shù)。本文將期望輸出與預(yù)測(cè)輸出之間的誤差絕對(duì)值作為個(gè)體的適應(yīng)度值F，計(jì)算方式如下：

其中m表示神經(jīng)網(wǎng)絡(luò)輸出節(jié)點(diǎn)總數(shù)，di表示神經(jīng)網(wǎng)絡(luò)第i個(gè)節(jié)點(diǎn)的期望輸出，pi表示神經(jīng)網(wǎng)絡(luò)第i個(gè)節(jié)點(diǎn)的預(yù)測(cè)輸出。

步驟3 選擇操作。本文采用輪盤(pán)賭法作為選擇策略，選擇適應(yīng)度好的個(gè)體組成新的種群，選擇算子的計(jì)算方式如下：

其中Fi表示個(gè)體i的適應(yīng)度值，N為種群中個(gè)體的總數(shù)。

步驟4 交叉操作。在種群中隨機(jī)選擇兩個(gè)個(gè)體按照一定概率進(jìn)行交叉操作，兩個(gè)父代個(gè)體Pk和Pl在第j個(gè)點(diǎn)位的實(shí)數(shù)交叉方式如下：

其中α是介于（0，1）之間得隨機(jī)數(shù)。

步驟5 變異操作。在種群中隨機(jī)選擇一個(gè)個(gè)體按照一定概率進(jìn)行變異操作，個(gè)體Pi在第j個(gè)點(diǎn)位變異方式如下：

其中Pmax表示個(gè)體的上界，Pmin表示個(gè)體的下界，v 表示當(dāng)前迭代次數(shù)，Vmax表示最大進(jìn)化次數(shù)，μ是介于（0，1）之間的隨機(jī)數(shù)。

步驟6 檢查新的個(gè)體是否符合最優(yōu)個(gè)體標(biāo)準(zhǔn)，若不符合，則返回步驟2；若符合，則將最優(yōu)個(gè)體分割為BP神經(jīng)網(wǎng)絡(luò)的權(quán)值和閾值，調(diào)整網(wǎng)絡(luò)參數(shù)，重復(fù)學(xué)習(xí)訓(xùn)練，直至達(dá)到所需精度或?qū)W習(xí)次數(shù)的上限為止。

（2）BP 神經(jīng)網(wǎng)絡(luò)部分。將經(jīng)過(guò)遺傳算法優(yōu)化得到的初始權(quán)值和閾值等個(gè)體，代入BP神經(jīng)網(wǎng)絡(luò)，輸入訓(xùn)練數(shù)據(jù)，進(jìn)行BP 神經(jīng)網(wǎng)絡(luò)訓(xùn)練。利用遺傳算法檢驗(yàn)個(gè)體的適應(yīng)度，無(wú)法達(dá)到適應(yīng)度標(biāo)準(zhǔn)的個(gè)體進(jìn)行選擇、交叉和變異操作，重新進(jìn)行BP神經(jīng)網(wǎng)絡(luò)訓(xùn)練，直至個(gè)體滿足適應(yīng)度標(biāo)準(zhǔn)或滿足精度要求為止。具體內(nèi)容如下：

步驟1 初始化網(wǎng)絡(luò)。X為神經(jīng)元輸入向量，n為輸入層節(jié)點(diǎn)數(shù)；Y為神經(jīng)元輸出向量，m為輸出層節(jié)點(diǎn)數(shù)。

其中輸入向量根據(jù)式（1）進(jìn)行量化，輸出向量分別代表拷貝數(shù)變異的3種狀態(tài)。

步驟2 計(jì)算隱含層輸出。隱含層節(jié)點(diǎn)數(shù)為p，輸入層與隱含層之間的權(quán)值為ωij，閾值為λa，隱含層輸出值Ha的計(jì)算方式如下：

其中f為激勵(lì)函數(shù)，計(jì)算方式如式（4）所示。

步驟3 計(jì)算輸出層輸出。隱含層與輸出層之間的權(quán)值為ωjk，閾值為λb，輸出層輸出值Hb的計(jì)算方式如下：

步驟4 計(jì)算誤差并驗(yàn)證。誤差e表示為預(yù)測(cè)輸出Hb與期望輸出Y之差，計(jì)算方式如下：

假定e0為允許的誤差最大值，若max(em)>e0，則調(diào)整BP神經(jīng)網(wǎng)絡(luò)參數(shù)，重新進(jìn)行學(xué)習(xí)訓(xùn)練，直至全部訓(xùn)練數(shù)據(jù)滿足誤差標(biāo)準(zhǔn)。

步驟5 更新權(quán)值和閾值。根據(jù)上述步驟計(jì)算出的誤差e更新ωij和ωjk，計(jì)算方式如下：

其中η為學(xué)習(xí)速率，介于（0，1）之間。閾值修正如下：

BP 神經(jīng)網(wǎng)絡(luò)算法的結(jié)果還取決于網(wǎng)絡(luò)的原始狀態(tài)，在局部搜索中使用BP 算法更為有效，GA 則擅長(zhǎng)全局搜索。因此，本文使用GA 優(yōu)化初始權(quán)值和閾值，并在解空間中找到一些更好的搜索空間，然后使用BP 算法在這些小的解空間中搜索最優(yōu)解。通常使用混合算法的效果要優(yōu)于僅使用BP神經(jīng)網(wǎng)絡(luò)算法的效果。

3 實(shí)驗(yàn)結(jié)果與分析

為防止細(xì)胞污染等外部因素造成的誤差，本文在不同的腫瘤純度和覆蓋深度中選擇樣本作為訓(xùn)練數(shù)據(jù)。利用經(jīng)過(guò)上述步驟訓(xùn)練完成的BP 神經(jīng)網(wǎng)絡(luò)模型，對(duì)測(cè)試數(shù)據(jù)集進(jìn)行預(yù)測(cè)。對(duì)于每個(gè)窗，算法輸出3個(gè)從tansig函數(shù)映射的0 到1 之間的值，分別表示CNV 這3 種狀態(tài)的概率，并根據(jù)最大的概率確定該窗的CNV 狀態(tài)。例如，代表擴(kuò)增的概率大于代表其他兩種CNV 狀態(tài)的概率，則該窗的CNV狀態(tài)為擴(kuò)增。

本文中BP神經(jīng)網(wǎng)絡(luò)結(jié)構(gòu)為3層，利用式（4）給出的激活函數(shù)，輸入層節(jié)點(diǎn)數(shù)為4，隱含層節(jié)點(diǎn)數(shù)為5，輸出層節(jié)點(diǎn)數(shù)為3，共4×5+5×3=35個(gè)權(quán)值，5+3=8個(gè)閾值，因此需要優(yōu)化43 個(gè)參數(shù)。交叉概率pc=0.2，變異概率pm=0.1，初始種群規(guī)模S=10。

3.1 仿真數(shù)據(jù)研究

為了驗(yàn)證CNV-GABP 方法的合理性和可靠性，首先使用模擬數(shù)據(jù)進(jìn)行實(shí)驗(yàn)。本文從基因組變異模擬方法中選擇了一種綜合模擬方法IntSIM[27]來(lái)生成測(cè)序數(shù)據(jù)。腫瘤純度設(shè)置為0.2、0.3和0.4，覆蓋深度設(shè)置為4x和6x，以生成不同配置的數(shù)據(jù)集。從標(biāo)準(zhǔn)參考基因組中選擇一條染色體Chr21作為IntSIM輸入序列，并針對(duì)每種配置模擬50個(gè)樣本，共計(jì)300個(gè)樣本。將CNV-GABP與4 種現(xiàn)有方法（CNAnator[4]、FREEC[7]、MFCNV[25]、GROM_RD[28]）在靈敏度、精度和F1評(píng)分上進(jìn)行比較，指標(biāo)的計(jì)算方式參考文獻(xiàn)[16]。每個(gè)指標(biāo)均取50 個(gè)樣本的測(cè)量均值，5種方法的性能比較結(jié)果如圖4所示。

圖4 仿真數(shù)據(jù)中5種方法的性能比較Fig.4 Performance comparison of five methods in simulation data

由圖4可以看出，腫瘤純度增加的同時(shí)，5種方法的檢測(cè)性能都在提高。例如腫瘤純度0.2、覆蓋深度均為4x 時(shí)GROM_RD 的靈敏度約為0.45，而腫瘤純度0.4 時(shí)約為0.56；5種方法F1 評(píng)分均值在腫瘤純度0.2、覆蓋深度均為4x 時(shí)約為0.59，而在腫瘤純度0.4 時(shí)約為0.71。就精度而言，在6類(lèi)樣本中GROM_RD獲得2個(gè)最高值，CNVnator 獲得1 個(gè)最高值，CNV-GABP 獲得3 個(gè)最高值；就靈敏度和F1 評(píng)分而言，在6類(lèi)樣本中CNV-GABP均獲得最高值。

檢測(cè)CNV邊界的準(zhǔn)確性同樣是衡量算法可靠性的指標(biāo)之一，邊界偏差越小，CNV 檢測(cè)方法性能越高。MFCNV在邊界偏差上取得了良好的效果，本文將CNVGABP與之比較，結(jié)果如圖5所示。

圖5 兩種方法的邊界偏差比較Fig.5 Boundary bias comparison of two methods

由圖5 可以看出，隨著腫瘤純度的增加，兩種方法的邊界偏差都隨之減小。實(shí)驗(yàn)結(jié)果表明，CNV-GABP在六類(lèi)樣本中的邊界偏差均小于MFCNV。因此，CNVGABP在模擬數(shù)據(jù)中的性能表現(xiàn)優(yōu)于其他4種算法。

CNV-GABP 的性能優(yōu)勢(shì)主要由三方面原因?qū)е拢阂环矫?，CNV-GABP在不同的腫瘤純度和測(cè)序深度中選擇大量訓(xùn)練數(shù)據(jù)進(jìn)行訓(xùn)練，提高了容錯(cuò)率；另一方面，其他方法只考慮到多個(gè)特征的邊際效應(yīng)，而CNV-GABP通過(guò)BP神經(jīng)網(wǎng)絡(luò)解決多個(gè)特征之間的聯(lián)合作用；最后，CNV-GABP 使用遺傳算法優(yōu)化BP 神經(jīng)網(wǎng)絡(luò)模型，增強(qiáng)了模型的魯棒性。

3.2 真實(shí)數(shù)據(jù)研究

為驗(yàn)證CNV-GABP 的有效性，應(yīng)用該算法分析3個(gè)實(shí)際測(cè)序樣本：NA19238、NA19239 和NA19240，數(shù)據(jù)可以在1 000 基因組項(xiàng)目網(wǎng)站http：//www.internationalgome.org/上獲取。本文將5 種檢測(cè)方法分別用于這3個(gè)樣本，檢測(cè)結(jié)果如表1 所示。基因組變異數(shù)據(jù)庫(kù)（http：//DGV.tcag.ca/）提供已確認(rèn)的CNV，可量化這5種方法的性能。

表1 5種方法檢測(cè)的CNV個(gè)數(shù)比較Table 1 Number comparison of detected copy number variations by 5 methods

由表1 展示的結(jié)果可以看出，在NA19238 樣本中，CNV-GABP 檢測(cè)的CNV 個(gè)數(shù)少于CNAnator，剩余的兩個(gè)樣本中，CNV-GABP 檢測(cè)的CNV 個(gè)數(shù)均多于其他方法。這表明針對(duì)特征不明顯的CNV 時(shí)，該算法具備較強(qiáng)的檢測(cè)能力。5種方法在真實(shí)數(shù)據(jù)上的性能比較結(jié)果如圖6所示。

由圖6 可以看出，3 個(gè)樣本中CNV-GABP 的F1 評(píng)分值均最高，平均值為0.85；就精度而言，MFCNV 表現(xiàn)優(yōu)于CNV-GABP，其次是CNVnator；就靈敏度而言，CNV-GABP 獲得兩個(gè)最高值，其次是MFCNV。因此，CNV-GABP 在實(shí)際測(cè)序樣本數(shù)據(jù)中的整體性能表現(xiàn)優(yōu)于其他4種算法，這說(shuō)明該算法在實(shí)際應(yīng)用中是有效且可靠的。

圖6 真實(shí)數(shù)據(jù)中5種方法的性能比較Fig.6 Performance comparison of five methods in real data

針對(duì)GC 含量偏差校正過(guò)程產(chǎn)生的誤差對(duì)現(xiàn)有CNV 檢測(cè)方法造成影響的問(wèn)題，本文提出了一種遺傳算法優(yōu)化的BP神經(jīng)網(wǎng)絡(luò)拷貝數(shù)變異檢測(cè)算法。該算法無(wú)需進(jìn)行GC含量偏差校正，采用BP神經(jīng)網(wǎng)絡(luò)解決不同特征之間的聯(lián)合作用，然后使用遺傳算法優(yōu)化BP 神經(jīng)網(wǎng)絡(luò)模型獲取最優(yōu)參數(shù)，最后使用優(yōu)化的BP 神經(jīng)網(wǎng)絡(luò)進(jìn)行CNV 預(yù)測(cè)。本文基于仿真數(shù)據(jù)和1 000 基因組計(jì)劃的真實(shí)數(shù)據(jù)，將CNV-GABP與其他4種現(xiàn)有算法進(jìn)行比較。實(shí)驗(yàn)結(jié)果表明，CNV-GABP算法在靈敏度、精度、F1 評(píng)分和邊界偏差上均獲得了明顯的優(yōu)勢(shì)，尤其提高了對(duì)特征不明顯的CNV的檢測(cè)能力。

對(duì)于將來(lái)的工作，將從以下幾方面著手對(duì)CNVGABP 方法進(jìn)行改進(jìn)和擴(kuò)展：（1）將算法擴(kuò)展至腫瘤-正常配對(duì)樣本的研究，對(duì)腫瘤樣本和正常樣本分別提取特征，并將提取后的特征作為CNV-GABP算法的輸入，以進(jìn)行差異基因分析；（2）單細(xì)胞測(cè)序技術(shù)蓬勃發(fā)展，下一步可收集此類(lèi)測(cè)序數(shù)據(jù)的基因組特征，將該算法擴(kuò)展至單細(xì)胞測(cè)序數(shù)據(jù)中，以分析單個(gè)細(xì)胞的基因組突變；（3）同時(shí)可根據(jù)混合優(yōu)化算法[29-31]優(yōu)化BP神經(jīng)網(wǎng)絡(luò)的結(jié)構(gòu)，提高算法運(yùn)行效率。

4 結(jié)束語(yǔ)