陳念華,袁細(xì)國(guó)

(西安電子科技大學(xué) 計(jì)算機(jī)科學(xué)與技術(shù)學(xué)院，陜西 西安 710071)

0 引言

拷貝數(shù)變異(Copy Number Variation，CNV)是人類基因組中一種重要的結(jié)構(gòu)變異類型，長(zhǎng)度通常在1K base pairs (bp)到3Mbp之間，包括拷貝數(shù)擴(kuò)增(amplification)和拷貝數(shù)缺失(deletion)兩種類型[1, 2]。人類基因在正常情況下是二倍體，所以對(duì)于人類基因組來說，拷貝數(shù)擴(kuò)增是指基因組區(qū)域的拷貝數(shù)從正常二倍體到多倍體的變化，拷貝數(shù)缺失則是基因組區(qū)域中拷貝數(shù)減少的變異，若拷貝數(shù)缺失至單倍體，稱作雜合性缺失；若拷貝數(shù)缺失至0，則稱作純合性缺失.研究表明，CNV在人類基因組中十分常見，它會(huì)引起基因表達(dá)發(fā)生異常，與人類復(fù)雜疾病緊密關(guān)聯(lián)，例如自閉癥[3]、精神分裂癥[4]、自身免疫性疾病[5]以及癌癥[6]等疾病。

自1975年第一代DNA測(cè)序技術(shù)開創(chuàng)至今，人類已經(jīng)積累了大量的測(cè)序數(shù)據(jù)，這使得利用計(jì)算機(jī)技術(shù)對(duì)這些數(shù)據(jù)進(jìn)行分析成為可能。相比于直接用醫(yī)學(xué)手段檢測(cè)CNV，利用計(jì)算機(jī)技術(shù)檢測(cè)CNV更加便捷，成本也十分低廉。當(dāng)下檢測(cè)CNV的主要難點(diǎn)在于如何區(qū)分驅(qū)動(dòng)CNV[7]和隨機(jī)CNV。所謂驅(qū)動(dòng)CNV，是指對(duì)疾病有直接影響或者關(guān)聯(lián)較大的CNV，找出驅(qū)動(dòng)CNV對(duì)理解疾病的發(fā)病機(jī)理有很大幫助；隨機(jī)CNV則是指在基因中隨機(jī)出現(xiàn)、與疾病的發(fā)生關(guān)聯(lián)不大的CNV。在多樣本檢測(cè)[8]的前提下，CNV按照在不同樣本中發(fā)生的頻率不同可以分為復(fù)發(fā)CNV[9]和個(gè)體CNV[10]，其中復(fù)發(fā)CNV指在不同患者基因組中相同位置發(fā)生的CNV，而個(gè)體CNV在不同患者基因組中發(fā)生的位置則是隨機(jī)的.研究表明，復(fù)發(fā)CNV更有可能是驅(qū)動(dòng)CNV，即更有可能包含疾病相關(guān)基因，因此本文算法的目標(biāo)就是從多樣本數(shù)據(jù)中檢測(cè)出復(fù)發(fā)CNV。

當(dāng)前有許多檢測(cè)復(fù)發(fā)CNV的方法，例如PLA(Piecewise-constant and low-rank approximation for identification of recurrent copy number variations)[11]是將多樣本CNV檢測(cè)問題轉(zhuǎn)化為矩陣分解問題，其中原始數(shù)據(jù)矩陣被分解為低秩分量，稀疏分量和噪聲分量。 這3個(gè)成分分別對(duì)應(yīng)復(fù)發(fā)CNV，個(gè)體CNV和隨機(jī)噪聲。FLLat(A fused lasso latent feature model for analyzing multi-sample aCGH data)[12]則是使用潛變量模型對(duì)基于陣列的比較基因組雜交技術(shù)(array-based Comparative genomic hybridization，aCGH)數(shù)據(jù)進(jìn)行建模，其中每個(gè)樣本均通過固定數(shù)量特征的加權(quán)組合來建模。這些特征代表了樣本組CNV的關(guān)鍵區(qū)域，并與權(quán)重相結(jié)合，描述了每個(gè)單獨(dú)樣本中的CNV區(qū)域。SAIC(Genome-wide identification of significant aberrations in cancer genome)[13]使用置換檢驗(yàn)方法來評(píng)估每個(gè)位點(diǎn)的重要程度，以此來檢測(cè)復(fù)發(fā)CNV。

如前所述，現(xiàn)有多樣本CNV檢測(cè)方法更關(guān)注數(shù)據(jù)的數(shù)學(xué)特性，而忽略了數(shù)據(jù)所包含的生物特性，因此本文提出一種基于聚類的從多樣本中檢測(cè)復(fù)發(fā)CNV的新算法DBSCAN-CNV(A recurrent copy number variation detection algorithm from multi-sample based on clustering)，該算法首先提取兩種與復(fù)發(fā)CNV的發(fā)生緊密關(guān)聯(lián)的特征，分別是每個(gè)位點(diǎn)發(fā)生CNV樣本的比率和每個(gè)位點(diǎn)的幅度均值，然后根據(jù)這兩個(gè)特征進(jìn)行聚類.由于發(fā)生復(fù)發(fā)CNV的位點(diǎn)相較于正常位點(diǎn)僅占少數(shù)，在整體數(shù)據(jù)中屬于異常點(diǎn)，因此本文采用的聚類方法為DBSCAN(A density-based algorithm for discovering cluster in large spatial databases with noise)，DBSCAN的優(yōu)勢(shì)在于可以對(duì)任意形狀的簇進(jìn)行聚類，并且如果對(duì)參數(shù)恰當(dāng)?shù)卦O(shè)定，它可以將噪聲點(diǎn)剔除，這可以解決發(fā)生復(fù)發(fā)CNV位點(diǎn)在全部數(shù)據(jù)中占比低的問題。

本文分別將該算法應(yīng)用在模擬數(shù)據(jù)和真實(shí)數(shù)據(jù)上，并與三種同行方法進(jìn)行比較(PLA、FLLat、SAIC)，實(shí)驗(yàn)結(jié)果表明，本算法對(duì)于復(fù)發(fā)CNV的檢測(cè)性能確實(shí)有一定提升。

1 方法

本文算法的流程如圖1所示，該算法通過以下4個(gè)主要步驟實(shí)現(xiàn)對(duì)復(fù)發(fā)CNV的檢測(cè)：(1) 數(shù)據(jù)預(yù)處理，這一步主要是將數(shù)據(jù)中的拷貝數(shù)信息(以2為基準(zhǔn))，轉(zhuǎn)化為以0為基準(zhǔn)的數(shù)據(jù)，即將原始數(shù)據(jù)除以2然后取對(duì)數(shù)，這樣方便對(duì)CNV的類型做區(qū)分：信號(hào)值小于0代表缺失，大于0代表擴(kuò)增(若數(shù)據(jù)本身就是以0為基準(zhǔn)的，跳過該步驟)；(2) 特征提取；(3) 根據(jù)上一步得到的特征進(jìn)行聚類；(4) 根據(jù)聚類結(jié)果獲得復(fù)發(fā)CNV的檢測(cè)結(jié)果。下面是對(duì)第(2) (3)步驟的詳細(xì)說明。

圖1 DBSCAN-CNV的主要步驟

1.1 特征提取

由于測(cè)序錯(cuò)誤、正常細(xì)胞污染等噪聲的存在，原始數(shù)據(jù)往往呈現(xiàn)為雜亂無章的信號(hào)序列，因此本文采用循環(huán)二元分割算法(Circular binary segmentation, CBS)[14]對(duì)每個(gè)樣本進(jìn)行分段平滑.如圖2所示，分段平滑后會(huì)將原始單個(gè)樣本數(shù)據(jù)分為多個(gè)連續(xù)區(qū)域，區(qū)域內(nèi)部的信號(hào)值是相同的。

圖2 使用CBS對(duì)原始數(shù)據(jù)進(jìn)行分段平滑，曲線是原始數(shù)據(jù)，直線是平滑過后的數(shù)據(jù)

1.1.1 每個(gè)位點(diǎn)發(fā)生CNV的樣本比率。在對(duì)每個(gè)樣本進(jìn)行分段平滑之后，根據(jù)分段區(qū)域內(nèi)的信號(hào)值越大，則代表該區(qū)域的拷貝數(shù)越大的原理，選定合適的閾值，判斷每個(gè)樣本在每個(gè)位點(diǎn)處是否發(fā)生CNV。對(duì)于拷貝數(shù)擴(kuò)增，選定正閾值，分段內(nèi)信號(hào)值若大于該閾值則認(rèn)為該分段內(nèi)的所有位點(diǎn)均發(fā)生拷貝數(shù)擴(kuò)增；相應(yīng)的，對(duì)于拷貝數(shù)缺失，選定負(fù)閾值，分段內(nèi)信號(hào)值若小于該閾值則認(rèn)為該分段內(nèi)的所有位點(diǎn)均發(fā)生拷貝數(shù)缺失。

經(jīng)過以上操作可以得到每個(gè)樣本在各位點(diǎn)發(fā)生CNV的情況，據(jù)此可以在每個(gè)位點(diǎn)計(jì)算發(fā)生CNV的樣本占總樣本的比率，即

freq(i)=count(i)/S，

(1)

其中count(i)指在第i個(gè)位點(diǎn)處發(fā)生CNV的樣本數(shù)，S指總樣本數(shù)。因?yàn)閺?fù)發(fā)CNV正是指那些在不同樣本間發(fā)生頻率較高的相同CNV，因此每個(gè)位點(diǎn)發(fā)生CNV的樣本比率是檢測(cè)復(fù)發(fā)CNV的重要特征。

1.1.2 每個(gè)位點(diǎn)的幅度均值。在經(jīng)過數(shù)據(jù)預(yù)處理之后，數(shù)據(jù)都是以0為基準(zhǔn)的(0代表拷貝數(shù)為2)，不論是大于0還是小于0都代表拷貝數(shù)發(fā)生了變異.對(duì)每個(gè)位點(diǎn)處各個(gè)樣本的信號(hào)值取絕對(duì)值，然后再取均值，便得到每個(gè)位點(diǎn)的幅度均值，其代表了每個(gè)位點(diǎn)的拷貝數(shù)均值與正常拷貝數(shù)偏離的程度，即

(2)

其中data[j,i]表示在第j個(gè)樣本，第i個(gè)位點(diǎn)處的信號(hào)值，S指總樣本數(shù).幅度均值越大，說明該位點(diǎn)的拷貝數(shù)偏離正常值越多，因此幅度均值也是檢測(cè)復(fù)發(fā)CNV的重要特征。

1.2 DBSCAN聚類

經(jīng)過上述操作，我們得到每個(gè)位點(diǎn)發(fā)生CNV的比率以及每個(gè)位點(diǎn)的幅度均值這兩個(gè)特征，接下來需要根據(jù)這兩個(gè)特征對(duì)所有位點(diǎn)進(jìn)行聚類。本文采用的聚類算法DBSCAN[15]，是一種基于密度的聚類方法，通過每個(gè)點(diǎn)Eps半徑之內(nèi)點(diǎn)的個(gè)數(shù)來衡量每個(gè)點(diǎn)的密度，如圖3所示，可以對(duì)任意形狀的數(shù)據(jù)進(jìn)行檢測(cè)。

圖3 基于中心的密度，點(diǎn)A的密度是9(包含A本身)

基于密度的聚類方法將數(shù)據(jù)集內(nèi)的點(diǎn)分為核心點(diǎn)、邊界點(diǎn)和噪聲點(diǎn)三類.核心點(diǎn)是在基于密度的簇內(nèi)部的點(diǎn).點(diǎn)的鄰域由距離函數(shù)和距離參數(shù)Eps決定。如果在一個(gè)點(diǎn)的半徑為Eps的鄰域內(nèi)，包含的點(diǎn)的個(gè)數(shù)超過閾值MinPts，則這個(gè)點(diǎn)是一個(gè)核心點(diǎn)；若某個(gè)點(diǎn)落在某個(gè)核心點(diǎn)的鄰域內(nèi)，并且該點(diǎn)是非核心點(diǎn)，則這個(gè)點(diǎn)為邊界點(diǎn)；噪聲點(diǎn)是既非核心點(diǎn)也非邊界點(diǎn)的任何點(diǎn)。圖4是這三種點(diǎn)的圖示。

圖4 若MinPts=7,則A是核心點(diǎn)，B是邊界點(diǎn)，C是噪聲點(diǎn)

本文距離函數(shù)使用歐氏距離，兩點(diǎn)之間的距離由公式定義

(3)

給定核心點(diǎn)、邊界點(diǎn)和噪聲點(diǎn)的定義后，DBSCAN算法可以非形式地描述如下：任意兩個(gè)相互距離在Eps之內(nèi)的核心點(diǎn)將放在同一個(gè)簇內(nèi)。落在某個(gè)核心點(diǎn)鄰域內(nèi)的邊界點(diǎn)和該核心點(diǎn)放在同一個(gè)簇內(nèi).噪聲點(diǎn)不屬于任何一個(gè)簇。下面是DBSCAN算法的詳細(xì)描述：(1) 將所有點(diǎn)標(biāo)記為核心點(diǎn)、邊界點(diǎn)或噪聲點(diǎn)；(2) 刪除噪聲點(diǎn)；(3) 為距離在Eps之內(nèi)的所有核心點(diǎn)之間賦予一條邊；(4)每組連通的核心點(diǎn)形成一個(gè)簇；(5) 將每個(gè)邊界點(diǎn)指派到一個(gè)與之關(guān)聯(lián)的核心點(diǎn)的簇中。

如前文所述，雖然DBSCAN的實(shí)現(xiàn)十分簡(jiǎn)單，但是檢測(cè)結(jié)果十分依賴半徑Eps的設(shè)定。如果設(shè)定的半徑足夠大，則所有點(diǎn)的密度都等于數(shù)據(jù)集中所有點(diǎn)的個(gè)數(shù)；類似地，如果半徑太小，則所有點(diǎn)的密度都是1(僅包含該點(diǎn)本身)。因此，可以通過觀察點(diǎn)到它的第k個(gè)最近鄰的距離(稱為k-距離)來選取合適的Eps。對(duì)于屬于某個(gè)點(diǎn)的簇，如果k不大于簇的大小的話，則k-距離將很小。然而對(duì)于不在簇中的點(diǎn)(如噪聲點(diǎn))，k-距離將相對(duì)較大。因此，如果我們對(duì)于某個(gè)k，計(jì)算所有點(diǎn)的k-距離，以遞增次序?qū)⑺鼈兣判颍缓罄L制排序后的值，則我們會(huì)看到k-距離的急劇變化，如圖5所示。

圖5 k距離的變化趨勢(shì)

選取k-距離發(fā)生急劇變化的點(diǎn)對(duì)應(yīng)的k-距離作為Eps是一個(gè)比較合適的值。如果我們選取該距離為Eps參數(shù)，而k的值作為MinPts參數(shù)，則k-距離小于Eps的點(diǎn)將被標(biāo)記為核心點(diǎn)，而其他點(diǎn)將被標(biāo)記為噪聲或邊界點(diǎn).由[15]可知，k=4對(duì)于大多數(shù)數(shù)據(jù)集都是一個(gè)合適的參數(shù)設(shè)定，因此本文算法默認(rèn)設(shè)k=4。Eps默認(rèn)取排序后的k-距離數(shù)組中第turn個(gè)位置的值，turn定義為

turn=P·turnPercent，

(4)

其中P是k-距離數(shù)組的長(zhǎng)度，turnPercent是k-距離導(dǎo)數(shù)驟增的點(diǎn)與k-距離數(shù)組長(zhǎng)度的比值，經(jīng)過實(shí)驗(yàn)?zāi)J(rèn)取0.9625。

因?yàn)榘l(fā)生復(fù)發(fā)CNV的位點(diǎn)在所有位點(diǎn)中所占比率很低，并且其特征與正常位點(diǎn)有顯著差異，因此DBSCAN的聚類結(jié)果中最大的簇代表非復(fù)發(fā)CNV位點(diǎn)，而剩下的簇代表發(fā)生在不同位置處的復(fù)發(fā)CNV位點(diǎn)。由于本文的目標(biāo)是檢測(cè)出發(fā)生復(fù)發(fā)CNV的位點(diǎn)，所以本文將DBSCAN聚類結(jié)果中除了最大簇以外的簇都作為檢測(cè)結(jié)果，噪聲點(diǎn)也視作檢測(cè)結(jié)果，至此對(duì)復(fù)發(fā)CNV的檢測(cè)全部完成。

2 實(shí)驗(yàn)結(jié)果

為了評(píng)估DBSCAN-CNV算法對(duì)復(fù)發(fā)CNV的檢測(cè)性能，本文將DBSCAN-CNV應(yīng)用在模擬數(shù)據(jù)上，并將PLA、FLLat、SAIC也應(yīng)用在這些數(shù)據(jù)上進(jìn)行比較.除此之外，本文還將DBSCAN-CNV應(yīng)用在真實(shí)數(shù)據(jù)上，看是否可以檢測(cè)出已被驗(yàn)證過的疾病相關(guān)基因，以此來驗(yàn)證該算法的可用性.接下來是對(duì)這些實(shí)驗(yàn)的詳細(xì)說明。

2.1 模擬數(shù)據(jù)

本文實(shí)驗(yàn)分別生成了兩種模擬數(shù)據(jù)，一種是根據(jù)文獻(xiàn)[16]的描述生成的包含六種場(chǎng)景的高頻率復(fù)發(fā)CNV數(shù)據(jù)，另一種則是本團(tuán)隊(duì)設(shè)計(jì)的相對(duì)低頻的復(fù)發(fā)CNV數(shù)據(jù)。

2.1.1 高頻率模擬數(shù)據(jù)。在文獻(xiàn)[16]里，作者詳細(xì)地定義了6種不同的復(fù)發(fā) CNV場(chǎng)景。本文根據(jù)其描述，在每種場(chǎng)景下生成50組數(shù)據(jù)，每組數(shù)據(jù)是50×2000的矩陣，其中50代表樣本數(shù)，2000代表位點(diǎn)數(shù)，即每一行數(shù)據(jù)都代表一個(gè)樣本.在生成數(shù)據(jù)時(shí)，將未發(fā)生CNV的位點(diǎn)的信號(hào)值設(shè)為0；復(fù)發(fā)CNV區(qū)域位于750-1250位點(diǎn)之間，其模式參考圖6。

圖6 在 Rueda and Diaz-Uriarte (2010)里定義的六種常見復(fù)發(fā)CNV的模式.每個(gè)場(chǎng)景的縱軸代表樣本，橫軸代表位點(diǎn)

將擴(kuò)增區(qū)域和缺失區(qū)域位點(diǎn)的信號(hào)值分別設(shè)為1和-1。每個(gè)樣本還需要在非復(fù)發(fā)CNV區(qū)域隨機(jī)選取一個(gè)位置，添加一個(gè)長(zhǎng)度為200的個(gè)體CNV，其值從{-2，-1，1，2}中隨機(jī)選取。最后再向整個(gè)數(shù)據(jù)添加噪聲水平為1的高斯噪聲，圖7是場(chǎng)景1和場(chǎng)景2模擬數(shù)據(jù)的生成過程示例。

圖7 場(chǎng)景1-2模擬數(shù)據(jù)的生成過程

本文用靈敏度(sensitivity)和準(zhǔn)確率(precision)來評(píng)估各個(gè)方法對(duì)模擬數(shù)據(jù)的檢測(cè)性能，其中靈敏度和準(zhǔn)確率的含義由公式定義

(5)

(6)

圖8是四種方法的檢測(cè)結(jié)果圖示，橫軸為準(zhǔn)確率，縱軸為靈敏度，圖中曲線是F1-score等高線，F(xiàn)1-score是準(zhǔn)確率和靈敏度的調(diào)和平均值，其定義為

圖8 不同方法對(duì)場(chǎng)景1-6模擬數(shù)據(jù)的復(fù)發(fā)CNV檢測(cè)結(jié)果的性能比較

(7)

F1-score值越大，說明算法性能越好，對(duì)應(yīng)到圖中就是越靠近右上方的點(diǎn)，F(xiàn)1-score值越大。

從圖中可以看出，除了在場(chǎng)景4里FLLat和DBSCAN-CNV的F1-score相近，在剩下的五種場(chǎng)景里DBSCAN-CNV的F1-score值都明顯比另外三種方法要大。比如在場(chǎng)景3里，雖然PLA、FLLat和DBSCAN-CNV的靈敏度幾乎都達(dá)到了1，但是PLA的準(zhǔn)確率只有0.657，F(xiàn)LLat的準(zhǔn)確率是0.801，而DBSCAN-CNV的準(zhǔn)確率卻達(dá)到了0.98；又比如在復(fù)發(fā)CNV模式比較復(fù)雜的場(chǎng)景6里，另外三種方法中F1-score最高的FLLat也只有0.810，而DBSCAN-CNV的F1-score卻達(dá)到了0.908，其中靈敏度是0.8424，準(zhǔn)確率是0.96。

2.1.2 低頻率模擬數(shù)據(jù)。低頻率數(shù)據(jù)按照以下幾個(gè)步驟生成：(1) 設(shè)定復(fù)發(fā)CNV的區(qū)域、變異樣本所占比率以及拷貝數(shù)擴(kuò)增(或缺失)的大小；(2) 添加個(gè)體CNV噪聲；(3) 添加腫瘤純度噪聲；(4) 添加高斯噪聲.接下來是對(duì)這些步驟的詳細(xì)介紹。

在低頻率的模擬數(shù)據(jù)里，本文按照?qǐng)D9的模式進(jìn)行復(fù)發(fā)CNV的模擬。每組數(shù)據(jù)是100×2000的數(shù)組，其中100代表樣本數(shù)，2000代表探針數(shù)，即每一行數(shù)據(jù)代表一個(gè)樣本。在生成數(shù)據(jù)時(shí)，首先使所有樣本的所有位點(diǎn)都為2(代表正常二倍體)，然后將拷貝數(shù)擴(kuò)增變異區(qū)域設(shè)定在100-149、500-529、900-919位點(diǎn)之間，設(shè)定發(fā)生CNV的樣本分別占總樣本數(shù)的0.2、0.25、0.20，拷貝數(shù)分別擴(kuò)增到6、4、5；將拷貝數(shù)缺失變異區(qū)域設(shè)定在1100-1149、1500-1529、1900-1919位點(diǎn)之間，設(shè)定發(fā)生CNV的樣本分別占總樣本數(shù)的0.15、0.2、0.25，拷貝數(shù)分別缺失至0、0/1、1。

圖9 低頻率模擬數(shù)據(jù)的復(fù)發(fā)CNV區(qū)間分布

設(shè)置復(fù)發(fā)CNV的變異區(qū)間后，在每個(gè)樣本的任意位置添加大小為100的個(gè)體CNV，其拷貝數(shù)在{0、1、3、4}中任意選取.在真實(shí)數(shù)據(jù)中，通常得到的信號(hào)數(shù)據(jù)并不是全部來自腫瘤細(xì)胞的，而是在正常細(xì)胞和腫瘤細(xì)胞混雜的情況下進(jìn)行信號(hào)的測(cè)量，這就導(dǎo)致了信號(hào)數(shù)據(jù)是被正常細(xì)胞“污染”過的數(shù)據(jù).腫瘤純度指腫瘤細(xì)胞占所有細(xì)胞的比率，腫瘤純度越高，說明數(shù)據(jù)受到正常細(xì)胞的“污染”越少。本文為了模擬這種在真實(shí)數(shù)據(jù)中不可避免的噪聲，每個(gè)樣本數(shù)據(jù)在0.3-0.7之間選擇一個(gè)腫瘤純度，然后與正常的拷貝數(shù)進(jìn)行加權(quán)平均，這樣便得到添加腫瘤純度噪聲之后的數(shù)據(jù)，最后對(duì)每個(gè)樣本添加指定水平的高斯噪聲。

在上述過程中，腫瘤純度有0.3-0.7五種選擇，高斯噪聲有0.2和0.4兩種選擇，通過對(duì)這兩種參數(shù)選擇的不同組合，共生成十組數(shù)據(jù)，每組有50個(gè)100×2000的模擬數(shù)據(jù)。

為了可量化地比較DBSCAN-CNV、PLA、FLLat、SAIC這四種方法在這十組數(shù)據(jù)上的檢測(cè)性能，這里依舊使用靈敏度(sensitivity)和準(zhǔn)確率(precision)以及它們的調(diào)和平均值F1-score作為衡量標(biāo)準(zhǔn)。

圖10是四種方法的檢測(cè)結(jié)果展示。從圖中可以看出，DBSCAN-CNV在大多數(shù)情況下的F1-score的值都是最大的，例如當(dāng)數(shù)據(jù)的noise(高斯噪聲水平)=0.4，purity(腫瘤純度)=0.3時(shí)，信號(hào)數(shù)據(jù)的各類噪聲水平是最大的，相應(yīng)的檢測(cè)難度也是最大的，DBSCAN-CNV的檢測(cè)結(jié)果依舊有0.683的靈敏度和0.815的準(zhǔn)確率，而FLLat的靈敏度只有0.21，準(zhǔn)確率只有0.56，SAIC的靈敏度雖然有0.79，高于DBSCAN-CNV，但準(zhǔn)確率只有0.27，所以F1-score依舊遠(yuǎn)小于DBSCAN-CNV，PLA在這組數(shù)據(jù)的靈敏度和準(zhǔn)確率為0。

圖10 不同方法對(duì)低頻率模擬數(shù)據(jù)的復(fù)發(fā)CNV檢測(cè)結(jié)果的性能比較，標(biāo)題中的noise值代表高斯噪聲水平，purity值代表腫瘤純度

相較另外三種方法，DBSCAN-CNV在噪聲水平高的數(shù)據(jù)里有著明顯優(yōu)勢(shì)，而在噪聲水平較低的數(shù)據(jù)里，雖然不能保證所有結(jié)果都是最好的，但是也都有著很高的靈敏度和準(zhǔn)確度。所以綜合來看，DBSCAN-CNV的性能穩(wěn)定，表現(xiàn)突出，是四種方法里最優(yōu)的。

2.2 真實(shí)數(shù)據(jù)

為了證明DBSCAN-CNV在真實(shí)數(shù)據(jù)上的可用性，本文實(shí)驗(yàn)將DBSCAN-CNV分別應(yīng)用在乳腺癌真實(shí)數(shù)據(jù)和肺癌真實(shí)數(shù)據(jù)上，并將檢測(cè)結(jié)果與現(xiàn)有研究已發(fā)現(xiàn)的疾病相關(guān)基因進(jìn)行比對(duì)，結(jié)果證明該方法可以檢測(cè)出正確的疾病相關(guān)基因。接下來是對(duì)這兩種數(shù)據(jù)檢測(cè)的詳細(xì)介紹。

2.2.1 乳腺癌真實(shí)數(shù)據(jù)。該數(shù)據(jù)包含了112個(gè)乳腺癌樣本的數(shù)據(jù)，每個(gè)樣本都有23條染色體上的不同數(shù)據(jù).在進(jìn)行實(shí)驗(yàn)時(shí)，首先把不同染色體的數(shù)據(jù)分割開來，形成23個(gè)不同的信號(hào)數(shù)據(jù)矩陣，然后對(duì)前22個(gè)數(shù)據(jù)進(jìn)行檢測(cè)(僅在常染色體上進(jìn)行檢測(cè))。由于真實(shí)數(shù)據(jù)更加雜亂無章，在檢測(cè)過程中聚類的數(shù)目會(huì)比較多，因此設(shè)定閾值T=0.1，如果某些類中包含點(diǎn)的個(gè)數(shù)加起來不超過閾值T，則認(rèn)為這些類代表的是發(fā)生復(fù)發(fā)CNV的位點(diǎn)。

以已報(bào)道的文獻(xiàn)為標(biāo)準(zhǔn)， DBSCAN-CNV在9號(hào)染色體上檢測(cè)出的與乳腺癌相關(guān)的基因如表1所示。例如DBC1基因，文獻(xiàn)[17]中說明了缺失DBC1對(duì)于乳腺癌的影響。又例如文獻(xiàn)[18]中說明了MTAP基因?qū)θ巳橄侔┘?xì)胞侵襲和遷移的影響。由此可知，DBSCAN-CNV具有從真實(shí)乳腺癌患者數(shù)據(jù)中檢測(cè)出乳腺癌相關(guān)基因的能力。

表1 DBSCAN-CNV在9號(hào)染色體上檢測(cè)出的與乳腺癌相關(guān)的基因

2.2.2 肺癌真實(shí)數(shù)據(jù)。 本文實(shí)驗(yàn)使用的肺癌真實(shí)數(shù)據(jù)[19]中共包含371個(gè)肺癌患者的樣本數(shù)據(jù)，每個(gè)樣本都包含了23個(gè)染色體上的所有數(shù)據(jù)。與對(duì)乳腺癌真實(shí)數(shù)據(jù)的處理類似，首先將數(shù)據(jù)按照不同染色體分割為23組數(shù)據(jù)，然后對(duì)不同染色體的數(shù)據(jù)分別做檢測(cè)。

以已報(bào)道的文獻(xiàn)為標(biāo)準(zhǔn)，表2匯總了DBSCAN-CNV在14號(hào)染色體上檢測(cè)出的與肺癌相關(guān)的基因。例如PAX9[20]曾被多項(xiàng)文獻(xiàn)報(bào)道其對(duì)肺癌的影響，而FOXA1[21]則與抑制肺癌抗腫瘤免疫力有關(guān)。由此可知，DBSCAN-CNV可以從真實(shí)肺癌數(shù)據(jù)中檢測(cè)出肺癌相關(guān)基因。

表2 DBSCAN-CNV在14號(hào)染色體上檢測(cè)出的與肺癌相關(guān)的基因

由DBSCAN-CNV對(duì)以上兩種真實(shí)數(shù)據(jù)的檢測(cè)結(jié)果可知，該方法可以對(duì)真實(shí)數(shù)據(jù)做出有價(jià)值的分析，這對(duì)于疾病的研究是十分重要的。

3 結(jié)論

復(fù)發(fā)CNV對(duì)人類復(fù)雜疾病的發(fā)生發(fā)展有著重要影響，因此研究復(fù)發(fā)CNV對(duì)于診斷治療這些疾病有很大意義。本文提出了一種基于聚類的可以從多樣本數(shù)據(jù)中檢測(cè)出復(fù)發(fā)CNV的算法DBSCAN-CNV，該算法首先將原始信號(hào)數(shù)據(jù)進(jìn)行平滑處理，然后提取各位點(diǎn)發(fā)生單樣本CNV的比率以及各位點(diǎn)的幅度均值，以這兩個(gè)特征作為聚類的特征；在聚類這一步，本文采用了DBSCAN聚類算法，該算法雖然實(shí)現(xiàn)簡(jiǎn)單但適用于本文的數(shù)據(jù)；最后根據(jù)聚類結(jié)果判定哪些位點(diǎn)發(fā)生了復(fù)發(fā)CNV。

在本文實(shí)驗(yàn)中，首先將DBSCAN-CNV分別應(yīng)用到高頻和低頻兩種模擬數(shù)據(jù)上，其中高頻數(shù)據(jù)共有6種不同場(chǎng)景，低頻數(shù)據(jù)共有10組不同的參數(shù)選擇，同時(shí)也將PLA、SAIC、FLLat這三種同行算法應(yīng)用在這些數(shù)據(jù)上，檢測(cè)結(jié)果以靈敏度和準(zhǔn)確率作為衡量指標(biāo)。實(shí)驗(yàn)結(jié)果表明，DBSCAN -CNV的性能顯著優(yōu)于另外三種方法。然后又將DBSCAN-CNV分別應(yīng)用在乳腺癌和肺癌真實(shí)數(shù)據(jù)集上，檢測(cè)結(jié)果中發(fā)現(xiàn)了現(xiàn)有文獻(xiàn)報(bào)道過的疾病相關(guān)基因，這表明該算法對(duì)于真實(shí)數(shù)據(jù)也具有可用性.綜上所述，DBSCAN-CNV對(duì)于復(fù)發(fā)CNV的檢測(cè)性能有著顯著提升。

針對(duì)計(jì)算復(fù)雜度，本文所提方法DBSCAN-CNV的時(shí)間復(fù)雜度近似為O(NlogN)，在實(shí)驗(yàn)過程中，與其他方法相比，本文方法的運(yùn)行時(shí)間較短。

在將來研究工作中，仍然存在不足以及可以進(jìn)一步擴(kuò)展的工作：(1) 模擬數(shù)據(jù)假設(shè)在不同樣本間發(fā)生復(fù)發(fā)CNV的位置完全相同，事實(shí)上它們的位置可能會(huì)有細(xì)小的差異，這可能對(duì)實(shí)驗(yàn)結(jié)果產(chǎn)生一定影響。(2) 本文僅在乳腺癌和肺癌的真實(shí)數(shù)據(jù)集上進(jìn)行了實(shí)驗(yàn)，然而隨著測(cè)序技術(shù)的發(fā)展，已經(jīng)積累了海量的疾病變異數(shù)據(jù)。因此下一步可以在其他疾病的真實(shí)數(shù)據(jù)集上進(jìn)行實(shí)驗(yàn)，以期發(fā)現(xiàn)更多與疾病相關(guān)的CNV，這將是十分有意義的研究。