邢 斌 周從華 張付全 張 婷 蔣躍明

（1.江蘇大學計算機科學與通信工程學院 鎮(zhèn)江 212013）（2.無錫市精神衛(wèi)生中心 無錫 214151）（3.無錫市婦幼保健院 無錫 214002）（4.無錫市第五人民醫(yī)院 無錫 214073）

1 引言

遺傳病是由致病基因所控制的疾病，這種由于遺傳物質(zhì)發(fā)生改變引起的疾病病種多，目前已發(fā)現(xiàn)的遺傳病超過3000種，而且具有發(fā)病率高、先天性、終生性、家族性等特點，對人類健康產(chǎn)生巨大影響。近年來隨著世界范圍內(nèi)人類基因組研究（Ge?nome-Wide Association Study，GWAS）開展實施，基因測序技術(shù)依靠對遺傳特征的有效挖掘，在疾病診斷、預測和治療等方面發(fā)揮著更加重要作用［2］，而GWAS是在單核苷酸多態(tài)（Single Nucleotide Poly?morphism，SNP）性的基礎(chǔ)上展開研究的，單核苷酸多態(tài)性主要是指基因組水平上的單核苷酸變異導致脫氧核糖核酸（DNA）序列多態(tài)性的現(xiàn)象，是人類最常見的一種遺傳變異，研究表明SNP對人類疾病有直接或間接的聯(lián)系，并且SNP具有穩(wěn)定性好、頻率高、采集容易等優(yōu)點，因此對SNP數(shù)據(jù)進行研究有重要意義。在此背景下，許多機器學習算法在SNP數(shù)據(jù)中得到了廣泛的應用。但是由于SNP數(shù)量過大，維數(shù)過高，數(shù)據(jù)中存在冗余和噪聲，研究中必須要考慮“維數(shù)災難”問題。因此，SNP數(shù)據(jù)分析的初始階段一般是選擇SNP中信息量最大的子集即信息SNP子集，以提高算法的性能，降低時間要求。

信息SNP子集的選擇本質(zhì)上是特征選擇，特征選擇（Feature Selection，F(xiàn)S）是在保持原始數(shù)據(jù)準確表示的同時，顯著降低特征空間的維數(shù)的過程。然而由于SNP數(shù)據(jù)之間的非相互獨立的特點，現(xiàn)有的特征選擇方法難以挖掘SNP位點之間的相關(guān)性，從而漏掉重要的遺傳信息，最終降低算法的分類效果。鑒于上述問題，本文采用基于信息論的新的相似度度量方法，并將其應用到K-means中，把SNP依據(jù)一定的關(guān)系劃分為若干簇，然后使用粒子群算法從每個簇中選出一個或多個能夠代表整個簇的信息SNP，構(gòu)造出最終的信息SNP子集。

2 相關(guān)工作

2.1 SNP選擇研究現(xiàn)狀

目前SNP選擇主要有兩種方法，一種是基于統(tǒng)計檢驗的關(guān)聯(lián)研究，分兩個步驟，第一步，通過生物實驗技術(shù)在全基因組上掃描篩選有效位點，第二步查驗SNP位點基因分型，通過關(guān)聯(lián)分析識別信息SNP；另一種是基于機器學習的SNP選擇，這種SNP選擇本質(zhì)上是特征選擇問題，可分為過濾式和封裝式兩種。過濾式設(shè)置評價指標給每個SNP打分，優(yōu)點是計算量小，但忽略了SNP之間的內(nèi)在聯(lián)系，不容易得到最優(yōu)的特征子集；而包裹式方法雖然計算量大，卻可以把學習器和評價指標結(jié)合起來，最終得到最優(yōu)解。近年來，許多人嘗試對這兩種方法進行改進，Raid Alzubi等［1］將過濾式和包裹式方法結(jié)合在一起提出一種混合特征選擇方法來檢測信息SNP并選擇最優(yōu)SNP子集，取得了不錯的效果。為解決SNP數(shù)據(jù)高維問題，Cong等［2］提出了一種基于主成分分析的基因數(shù)據(jù)降維算法，用于低維空間SNP位點的聚類。上述兩篇文獻雖然在結(jié)果上較前人有一定的進步，但是仍忽略了SNP與SNP子集之間的相關(guān)性問題，Liao［3］將信息論引入到SNP選擇中，用于衡量SNP子集之間的相關(guān)性，該方法顯著提高了標簽SNP選擇的效率和預測精度，但該方法沒有將SNP之間的相似度和聚類有效結(jié)合在一起。

2.2 K-means聚類算法

K-means算法是一種迭代求解的無監(jiān)督聚類算法，通過選取K個對象作為初始聚類中心計算每個樣本點與聚類中心的距離即相似度，將每個樣本點分配到最近的聚類中心，最終使得相同簇內(nèi)的樣本點相似度盡可能高，不同簇間的相似度盡可能低。

2.2.1 歐氏距離

歐氏距離是衡量兩個同一樣本集樣本對象差異性，距離越大樣本差異度越大。歐氏距離公式如下：

式中，xiμ為樣本對象xi的第μ個屬性。

2.2.2 簇內(nèi)誤差平方和

在歐氏距離基礎(chǔ)上，進行變形，將其中一個樣本對象換為簇中心得到簇內(nèi)誤差平方和，簇內(nèi)誤差平方和用來度量聚類效果的好壞。

3 基于改進的K-means算法的SNP選擇方法

3.1 傳統(tǒng)K-means在SNP選擇上的缺陷

傳統(tǒng)的K-means聚類使用SNP之間的歐氏距離來度量SNP之間的相似度，這并不能挖掘出SNP位點之間生物學上的聯(lián)系性，鑒于上述問題［1］將互信息引入到聚類中，用作SNP相似度度量，雖然起到很好的效果，但仍沒解決傳統(tǒng)K-means算法中另一缺陷，只考慮兩個SNP位點之間的相似度，然而在實際SNP數(shù)據(jù)集中，單個SNP往往和某個SNP子集有著強關(guān)聯(lián)性，鑒于上述兩個問題，提出一種新的基于信息論的SNP選擇算法K-MIGS。

3.2 基于信息論的相似度度量

3.2.1 信息熵與互信息

信息熵用來衡量信息的不確定性，單個SNP位點的信息熵可表示為

假設(shè)樣本S中某一SNPX有T個可能取值，其中p(xi)表示位點X中第i類樣本出現(xiàn)的頻率。Tx表示位點X所有屬性的個數(shù)。聯(lián)合信息熵：

則SNP位點X和Y的互信息可表示為

我們使用兩個位點的互信息來衡量它們之間的相似度，互信息越大，表示相似度越大，每個特征與初始簇中心的距離度量公式表示如下：

接下來我們試圖利用MI(S1;SNPX)和MI(S2;SNPX)來構(gòu)建SNPx與SNP子集之間的相似度，有下式：

假設(shè)n個SNP位點上存在m種，那么SNP子集內(nèi)聯(lián)合熵可表示為

但是如果S1、S2兩個子集SNP個數(shù)是不相等的，那么I(S1;SNPX)和I(S2;SNPX)將沒有可比性，因為單個SNP和大的SNP子集會得到大的互信息值，為了解決這個問題，我們結(jié)合了互信息和信息熵提出了一種單個SNP和SNP子集的相似度度量方式MIGS，公式如下：

其中，H(S)、H(SNPX)分別表示SNP子集S和單個SNP的信息熵，MI(S;SNPX)表示它們之間的互信息，考慮到SNP位點之間存在強相關(guān)性的特性，我們在原K-means聚類的距離度量歐氏距離中引入互信息的概念，則第一輪迭代計算中，每個特征與初始簇中心的距離度量公式表示如下：

其中MId(xi，uj)表示特征xi與初始簇中心uj的相似度，表示特征與初始簇中心的歐氏距離，MIGS(xi，Ci)表示特征x i與初始簇i之間的相似度。

3.2.2 簇中心的更新

在K-MIGS中，初始的簇中心的選擇與傳統(tǒng)K-means一致，不同的是在簇中心更新時，傳統(tǒng)K-means采用均值向量的方式更新簇中心，這并不適合用在K-MIGS算法中，因此我們采用新的簇更新方式，選取與均值向量最近的一個SNP作為簇中心。

3.3 算法K-MIGS步驟

結(jié)合章節(jié)3.1和3.2，則算法K-MIGS的整體步驟如算法1所示。

算法1：K-MIGS算法

輸入：數(shù)據(jù)集D={x1，x2，…，xm}聚類簇數(shù)kMAX_N

1）從數(shù)據(jù)集D={x1，x2，…，xm}中隨機選擇k個樣本作為初始均值向量(μ1，μ2，…，μk)

2）for j=1 to k do

3） forj=1 j=1 to mdo

4）若簇CiSNP數(shù)為1，根據(jù)式（6）、（7）計算與μj的相似度距離，否則根據(jù)式（10）、（11）計算SNP子集S與SNPi的相似度距離MId(S1，SNPi)

5） end for

6）end for

7）for i=1 to m do

8）將xi劃入與其相似度距離最小的簇

9）end for

10）repeat

11）for i=1 to k do

12）計算新的均值向量，并通過式（1）計算x( )x∈Ci與μ'i的距離

13）選擇最小的d(x，μ'i)作為簇Ci新的簇中心

14）end for

15）until迭代次數(shù)達到閾值或簇中心不再更新

其中，2～6行根據(jù)式（6）、（9）、（10）分兩種情況計算SNP之間的相似度距離，分為當前簇SNP數(shù)為一和SNP數(shù)不為一；7～9行將所選SNP劃分到相似度距離最小的簇，11～13行根據(jù)式（1）將簇中心更新為與均值向量歐氏距離最小的SNP。

3.4 K-MIGS算法在SNP選擇中的應用

本節(jié)將結(jié)合K-MIGS算法和粒子群算法對每個簇中的SNP進行選擇。

粒 子 群 算 法（Particle Swarm Optimization，PSO）是一種群體智能算法，該算法模仿鳥群群體合作來最快獲得食物的基本原理。在粒子群算法中，每一個候選解都可以看作是“鳥群中的個體”，即一個粒子。在粒子運動過程中每個粒子根據(jù)自身經(jīng)驗和相鄰粒子經(jīng)驗，通過對速度進行改變來調(diào)整自己的位置。粒子通過最佳粒子流在問題空間內(nèi)運動，迭代一定次數(shù)或者達到給定的最小誤差得到最優(yōu)解。

將粒子群算法應用到SNP選擇中包括兩個主要部分，種群初始化和粒子更新。

我們?yōu)镾NP選擇問題設(shè)計了種群：

其中n代表種群P的大小，m代表SNP的數(shù)量，pij=1表示選擇了第i個粒子的第j個SNP。

通常，在種群的初始化過程中每個SNP的選擇為信息SNP的概率由確定，其中k代表選擇的SNP的數(shù)量，m代表SNP的總數(shù)。在每次迭代中依據(jù)pbest和gbest的參數(shù)對粒子進行更新，每個粒子根據(jù)以下方程式更新：

其中w是慣性權(quán)重，a1和a2是加速學習因子，r1和r2為0～1的隨機數(shù)，和為更新前和更新后的粒子i速度矢量的第j維分量，為更新前粒子i位置矢量的第j維分量，pbestij=(pi1，pi2，…，pij)表示粒子i個體經(jīng)歷過的最好位置，gbestj=(g1，g2，…，gj)表示種群所經(jīng)歷過的最好位置。

將實驗所需原始數(shù)據(jù)進行預處理，并使用K-MIGS算法進行聚類，最后使用粒子群算法從每個簇中選擇得到最終的信息SNP。

4 實驗

4.1 實驗環(huán)境和數(shù)據(jù)

實驗環(huán)境：編譯工具Python 3.7.3，操作系統(tǒng)Windows10 64位，處 理 器Intel（R）Core（TM）i7-6700HQ，CPU@2.60GHz 2.59GHz，GPU Nvidia GeForce GTX 1060 3G，運行內(nèi)存16G，硬盤1.25T。

實驗數(shù)據(jù)集：實驗數(shù)據(jù)來自于無錫精神衛(wèi)生中心，包括兩種SNP數(shù)據(jù)集EN1000（9445SNPS）、E144（2514825SNPS），另外還包含樣本的基因型，患病與否的標記，具體描述如表1所示。

表1 數(shù)據(jù)集描述

4.2 實驗評價指標

1）SNP選擇評價指標

本文使用SNP選擇常用評價指標，信息SNP子集對非信息SNP子集重構(gòu)準確度來評價SNP選擇效果的好壞。

其中pi為信息SNP子集對非信息SNP位點i的重構(gòu)準確度，||O是非信息SNP的數(shù)量。重構(gòu)準確度越高，SNP選擇效果越好。

2）分類預測評價指標

本文使用分類實驗常用的評價指標預測的準確率（Accuracy）及F-Measure對分類結(jié)果進行評價，分類結(jié)果的“混淆矩陣”如表2所示。

表2 預測類別與實際類別的“混淆矩陣”

根據(jù)上表分類結(jié)果評價指標可由式（14～17）計算：

4.3 數(shù)據(jù)預處理

SNP數(shù)據(jù)預處理分為兩個子階段，數(shù)據(jù)編碼和數(shù)據(jù)更新，原始SNP具有三種基因型，純合子基因型（AA），雜合子基因型（Aa）以及純合變異基因型（aa），這種基因型數(shù)據(jù)不利于后續(xù)聚類操作，需要對數(shù)據(jù)進行編碼，分別將AA、Aa、aa編碼為0、1、2；數(shù)據(jù)編碼后對數(shù)據(jù)進行更新，更新又包括缺失數(shù)據(jù)填充及不符合標準數(shù)據(jù)的刪除。依據(jù)上述原則處理后，數(shù)據(jù)集G1000和E144分別剩余9298和214935條有效的SNP。

4.4 實驗結(jié)果與分析

4.4.1 聚類實驗及分析

此部分實驗包括四種算法比較實驗，分別為K-means、特征加權(quán)K-means、模糊聚類算法FCM和K-MIGS，主要評價指標為最后得到不同簇下信息SNP子集的重構(gòu)準確度。在數(shù)據(jù)集G1000和E144上，分別用四種算法進行聚類實驗，使用粒子群算法對每個簇進行信息SNPs提取，最后計算信息SNP子集對非信息SNP子集的重構(gòu)度。實驗結(jié)果如圖1和圖2所示。

圖1 G1000上算法選出的信息SNP對非信息SNP的重構(gòu)度圖

圖2 E144上算法選出的信息SNP對非信息SNP的重構(gòu)度圖

由圖1、2可看出，使用K-MIGS/粒子群算法最終提取的信息SNP在兩個數(shù)據(jù)集上對非信息SNP具有更好的重構(gòu)度，并且在聚類簇數(shù)為8時重構(gòu)度達到最大值，因此后續(xù)分類實驗采用簇數(shù)為8。

4.4.2 分類實驗及分析

分類實驗的目的是更進一步檢測K-MIGS/粒子群算法所選擇的信息SNP子集包含信息的重要程度，此實驗中，采用K-means/粒子群算法、K-MIGS/粒子群算法、特征加權(quán)K-means/粒子群算法（FW-K-means/粒子群）、ReliefF和MCMR算法進行信息SNP子集的篩選，使用SVM、DT和神經(jīng)網(wǎng)絡（Neural Networks，NN）為分類器，主要評價指標為分類準確率Acc和F1-measure。實驗結(jié)果如表3所示。

表3 不同的分類器進行SNP子集評價的結(jié)果

5 結(jié)語

針對SNP本身具有的高維少樣本特性以及自身存在的遺傳規(guī)律，本文提出了一種基于K-means的算法K-MIGS并將其應用到SNP選擇中。對數(shù)據(jù)進行預處理后，使用提出的方法對SNP數(shù)據(jù)進行聚類，最后使用粒子群算法構(gòu)造最終的SNP子集。在聚類效果評估和分類實驗中均表明，該方法很大地提升了SNP選擇的有效性。本文提出的方法相比其他SNP選擇方法優(yōu)勢在于同時考慮了單個SNP與SNP子集的相似度對聚類結(jié)果的影響。本文的后續(xù)工作主要有兩點：一是繼續(xù)優(yōu)化算法以減少單個SNP與SNP子集的相互關(guān)聯(lián)時引入的額外時間復雜度；二是繼續(xù)優(yōu)化粒子群算法，使其在每個簇中選擇出信息量更大的SNP。