宮曉娜，孔 薇

（上海海事大學(xué) 信息工程學(xué)院，上海 201306）

阿爾茨海默?。ˋlzheimer disease，AD）是由巴伐利亞的神經(jīng)病理學(xué)家阿爾茨海默[1]（Alois Alzheimer）于1907年首先發(fā)現(xiàn)，并以其名字而命名，這是一種漸進(jìn)性的神經(jīng)變性疾病，這種疾病表現(xiàn)為全面的認(rèn)知障礙，包括：記憶、定位、判斷和推理。約60%-70%老年癡呆患者為AD型[2]：AD患病率隨著年齡增加而明顯上升，60-64歲人群患病率為1%，而85歲以上老年人中則高達(dá)40%。AD給社會帶來沉重的經(jīng)濟(jì)及社會負(fù)擔(dān)。

近些年來隨著生物芯片技術(shù)的快速發(fā)展，研究人員可以同時測定不同樣本中成千上萬的基因表達(dá)水平，為我們進(jìn)行相關(guān)研究提供數(shù)據(jù)基礎(chǔ)，但是每個樣本都包含成千上萬個基因，而我們采集的樣本只有有限多個（最多幾百個，對我們的研究造成很大的困難），即維數(shù)災(zāi)難（少量的樣本對應(yīng)著巨大數(shù)量的特征）。

針對基因表達(dá)數(shù)據(jù)的這個特點(diǎn)，采取了以下兩種相應(yīng)的解決方法[3]：

1）進(jìn)行特異性基因選擇，在分類的時候只使用那些與分類緊密相關(guān)的基因，這樣可以篩除與分類無關(guān)的基因，減小維數(shù)、噪聲和冗余，從而降低計算復(fù)雜度并提高分類的準(zhǔn)確度和可靠度，從而減少用于實(shí)際臨床的診斷費(fèi)用。

2）構(gòu)造分類器時采用交叉檢驗(yàn)法。由于樣本數(shù)量太少，把樣本分為訓(xùn)練集和測試集幾乎是不可能的，從而使用交叉檢驗(yàn)法對樣本進(jìn)行重復(fù)利用。

1 特異性基因提取

目前，針對差異表達(dá)基因的基因排序選擇算法有很多方法[4-5]，包括傳統(tǒng)的倍數(shù)法、t檢驗(yàn)和方差分析、Wilcoxon非參數(shù)發(fā)，Bonferroni修正法等，還有專門針對芯片數(shù)據(jù)特點(diǎn)的SAM[6]（significance analysis of microarray）方法等，這些方法或多或少存在某些缺陷[7]。常用的特征選取的方法包括主成分分析方法、因子分析方法、獨(dú)立成分分析方法、小波變換方法和離散余弦變換方法等。相比t檢驗(yàn)的方法，小波包變換-SAM對基因表達(dá)譜可以更有效的將不相關(guān)的冗余基因剔除。目前的聚類方法用于特征基因選取只能把一個基因歸類到一類中，這不符合生物學(xué)特性，即一個基因可以參與多個信號傳導(dǎo)過程；而ICA方法是一種雙向聚類方法，能夠把基因歸入到不同類別中，從而較好的反應(yīng)基因在不同信號傳導(dǎo)通路中的作用[8]。本文提出的先小波包變換-SAM-FastICA算法可以充分利用每種算法的優(yōu)點(diǎn)，克服它們的不足。實(shí)驗(yàn)表明這種方法比單純的對基因表達(dá)譜數(shù)據(jù)進(jìn)行ICA分析具有更高的準(zhǔn)確度。

1.1 小波-SAM原理

由于微陣列數(shù)據(jù)自身含有很大的噪聲，如何有效的消除基因表達(dá)譜數(shù)據(jù)的各種噪聲，對于特異性表達(dá)基因的篩選和后續(xù)的生物學(xué)分析等有著重要的影響。

小波變換具有良好的時頻局部化特性，因而能有效的從信號中提取資訊，通過伸縮和平移等運(yùn)算功能對函數(shù)或信號進(jìn)行多尺度細(xì)化分析（Multiscale Analysis）。小波去噪是將信號映射到小波域，根據(jù)噪聲和噪聲的小波系數(shù)在不同尺度上具有不同的性質(zhì)和機(jī)理，對含噪信號的小波系數(shù)進(jìn)行處理。實(shí)質(zhì)是減少剔除噪聲產(chǎn)生的小波系數(shù)，最大限度的保留真實(shí)信號的系數(shù)。

由美國斯坦福大學(xué)開發(fā)的SAM軟件作為插件被安裝到Excel軟件中使用，該軟件專門用于篩選差異表達(dá)基因，其算法是基于傳統(tǒng)的t檢驗(yàn)和方差分析，用統(tǒng)計量衡量基因表達(dá)與反應(yīng)變量之間關(guān)系的強(qiáng)弱。通過數(shù)據(jù)的重復(fù)排列或抽樣來計算FDR，調(diào)節(jié)FDR篩選特異性表達(dá)基因。這種分析方法允許一定的假陽性率，適用于發(fā)現(xiàn)性的實(shí)驗(yàn)方法，有利于發(fā)現(xiàn)低拷貝或差異小的基因[9]。

1.2 ICA（Independent Component Analysis）算法原理

ICA算法的實(shí)質(zhì)是在假設(shè)源信號統(tǒng)計獨(dú)立的基礎(chǔ)上，在不知道源信號及混合矩陣任何信息的情況下，試圖將一組隨機(jī)變量表示成統(tǒng)計獨(dú)立的變量的線性組合。

設(shè) X=（x1，…，xn）T為 n 維隨機(jī)觀測向量，由 n 個未知源信號 S=（s1，s2，…，sn）T線性組合而成。用矩陣形式來定義 ICA 線性模型

其中 si稱為獨(dú)立成分 （Independent Component，IC），A=（a1，a2，…，an）∈Rn×n為一滿秩矩陣，稱為混合矩陣。 由此可以看出，觀測數(shù)據(jù)X是由獨(dú)立信號源S經(jīng)A線性加權(quán)得到的。進(jìn)行ICA處理的目的就是找到混合信號X的一個線性變換矩陣W，使得輸出盡可能的獨(dú)立，即

ICA實(shí)際上是一種尋優(yōu)過程，即如何使分離出的獨(dú)立成分最大限度地逼近各源信號?？梢酝ㄟ^改變P中的系數(shù)來觀察Y=PS的分布如何變化。因此，ICA包括兩個主要方面：目標(biāo)函數(shù)和尋優(yōu)算法。

1）極大化非高斯性的ICA目標(biāo)函數(shù)

由ICA的估計原理2：極大非高斯性：在y的方差為常數(shù)的約束下，求線性組合非高斯的局部極大值。每個局部極大值給出一個獨(dú)立成分。

在實(shí)際應(yīng)用中，可以使用峭度來度量非高斯性。y的峭度kurt（y）可以定義為：

2）采用峭度的梯度算法

對混合量進(jìn)行白化預(yù)處理，意味著將X線性變換成一個隨機(jī)變量:

在實(shí)踐中為了極大化峭度的絕對值，可以從某個向量W開始，依據(jù)可用的樣本值 Z（1），…，Z（T），計算出使 Y=WTZ 的峭度絕對值增大的最快的方向，然后將向量W轉(zhuǎn)到該方向。利用梯度的原理，WTZ的峭度的絕對值的梯度可以用下式計算得到：

對于白化過的數(shù)據(jù)，有E{（WTZ）2}=‖W‖2。因此，可以得到下面的梯度算法：

1.3 FastICA算法原理

上一節(jié)介紹了以峭度絕對值作為非高斯度量，導(dǎo)出了極大化非高斯的一種梯度優(yōu)化方法。該梯度的算法收斂慢，且依賴于合理的學(xué)習(xí)速度序列的選擇，如果學(xué)習(xí)速度選擇不當(dāng)，收斂性可能會被破壞。針對這個問題，1997年芬蘭赫爾辛基大學(xué)的 A.Hyv?rinen和 E.Oja等人提出的快速固定點(diǎn)算法（FastICA）[10]。

在梯度算法的一個穩(wěn)定（收斂）點(diǎn)出，梯度必須指向W的方向，也就是說梯度必須等于一個常數(shù)標(biāo)量與W的乘積。只有在這種情況下，將梯度與W相加才不改變其方向，且算法在此處收斂。令公式（5）中峭度的梯度與W相等，可以得到：

由該公式直接隱含著一個不動點(diǎn)算法，可以首先計算右邊的項(xiàng)，然后將其賦給W作為新值：

最后收斂的W以WTZ的線性組合形式可以給其中一個獨(dú)立成分。

該算法能夠快速且可靠地收斂。與ICA相比，該算法沒有學(xué)習(xí)速度或其它可調(diào)節(jié)的參數(shù)，因此，F(xiàn)astICA更易用且更可靠。

1.4 小波-SAM-FastICA在基因表達(dá)譜中的應(yīng)用

本次實(shí)驗(yàn)選用的是GEO數(shù)據(jù)庫Series GSE5281中海馬區(qū)（HIP）的23組AD樣本數(shù)據(jù)集（13組control無病樣本，10組affected AD樣本，54 675個基因表達(dá)數(shù)據(jù)）。由于基因表達(dá)譜數(shù)據(jù)具有高噪聲、高維性及數(shù)據(jù)分布不均衡等特點(diǎn)，因此，在對數(shù)據(jù)集進(jìn)行FastICA特征提取之前，首先要對數(shù)據(jù)集進(jìn)行預(yù)處理，即將大量的無關(guān)基因通過一定的算法篩選掉。這里采用小波包變換-SAM（Significance Analysis of Microarrays）的預(yù)處理機(jī)制，預(yù)處理后基因表達(dá)數(shù)據(jù)中的噪聲大大地降低，并且能夠更好地描述樣本特征，有利于后期特異性基因的提取。

1.4.1 對HIP數(shù)據(jù)集預(yù)處理

首先對HIP數(shù)據(jù)集進(jìn)行小波-SAM預(yù)處理，圖1為不同預(yù)處理的數(shù)據(jù)聚類效果。

圖1 不同預(yù)處理的數(shù)據(jù)聚類效果Fig.1 Data clustering effect of different pretreatment

（a），（b）分別為不同預(yù)處理后的數(shù)據(jù)聚類效果。從圖1可以看出經(jīng)過小波包-SAM篩選后的基因，經(jīng)過聚類后，可以準(zhǔn)確地將control無病樣本和affected AD樣本分開，證實(shí)了小波包-SAM預(yù)處理的有效性。

1.4.2 對基因表達(dá)譜提取特征基因

由于混合矩陣A可以反映特征樣本在樣本基因表達(dá)譜中的活躍情況，文中主要研究混合矩陣A的權(quán)重來分析基因表達(dá)譜。

圖中白色和黑色分別代表正值和負(fù)值，而每個正方形的大小則對應(yīng)于每個樣本中成分的數(shù)量。原始數(shù)據(jù)集由13個無病樣本（前 13行）和 10個AD患病 樣本（最后 10行）組成。 根據(jù)正負(fù)值，從圖 2（a）、2（b）中可以看出，圖 2（b）第 8、9、10、11、12、13、18、21、23 列可以將無病樣本和患病樣本區(qū)分開，而圖2（a）基本沒有將無病樣本和患病樣本區(qū)分開。

圖2 未經(jīng)預(yù)處理和經(jīng)小波變換的對比圖Fig.2 Camparison chart of without pretrea tment and after the wavelet packet transform

根據(jù)混合矩陣A的權(quán)重圖，提取FastICA分析后的相對應(yīng)的第 8、9、10、11、12、13、18、21、23 個特征樣本。

圖3 Histogram of the characteristic sampleFig.3 特征樣本的柱狀圖

由于過高表達(dá)和過低表達(dá)的基因很有可能是阿爾茨海默癥致病基因或者相關(guān)基因。這里特征樣本的柱狀圖，將基因表達(dá)值大于5或者基因表達(dá)值小于-5的基因提取出來，得到172個特征基因。

2 基于學(xué)習(xí)向量量化（LVQ）神經(jīng)網(wǎng)絡(luò)的分類

學(xué)習(xí)向量量化（LVQ）神經(jīng)網(wǎng)絡(luò)[10]是一種有監(jiān)督的訓(xùn)練競爭層的方法。學(xué)習(xí)向量量化網(wǎng)絡(luò)能夠?qū)θ我廨斎胂蛄窟M(jìn)行分類，不管它們是不是線性可分，這點(diǎn)比感知器神經(jīng)網(wǎng)絡(luò)要優(yōu)越得多。

實(shí)驗(yàn)分析：

將FastICA提取的23個樣本（訓(xùn)練樣本6個，測試樣本17個）的172個基因進(jìn)行LVQ神經(jīng)網(wǎng)絡(luò)分類，在MATLAB里多次運(yùn)行，得到以下結(jié)果：

從表格中可以看出，經(jīng)過小波包變換-SAM-FastICA提取的特征基因比小波包-SAM-ICA提取的特征基因更能高精確地將樣本分類，從而驗(yàn)證了算法的有效性。

表1 對提取的基因進(jìn)行LVQ神經(jīng)網(wǎng)絡(luò)分類準(zhǔn)確率結(jié)果Tab.1 LVQ neural network classification accuracy results of the extracted gene

3 基因集合富集分析

基因功能富集分析[11]（Gene Set Enrichment Analysis，GSEA）又稱功能聚類分析，利用超幾何分布型統(tǒng)計原理來檢驗(yàn)一組基因（共表達(dá)或差異表達(dá)）中某個功能類的顯著性，通過離散分布的顯著性分析、富集度分析和假陽性分析，得出與實(shí)驗(yàn)?zāi)康挠酗@著關(guān)聯(lián)的、低假陽性率的及靶向性的基因功能類別，該功能類別即是導(dǎo)致樣本性狀差異的最重要的功能差別，而其所屬基因是需進(jìn)一步驗(yàn)證的重要目標(biāo)基因，其功能特征將闡明樣本性狀變化的內(nèi)在生物學(xué)意義。

利用基因功能富集分析的功能特點(diǎn)，本文對提取的特征基因進(jìn)行基因功能富集分析從而分析提取的特征基因在阿爾茨海默癥數(shù)據(jù)集中的效果。

將之前進(jìn)行FastICA分析篩選后的172個特征基因，經(jīng)過基因富集性分析（172個基因里有139個基因匹配）分析得到它們在阿爾茨海默癥數(shù)據(jù)集中的聚類情況及其它們的基因表達(dá)情況。

圖4 ALZHEIMERS_DISEASE_DNBlue-Pink O'Gram in the Space of the Analyzed GeneSetFig.4 阿爾茨海默癥下調(diào)基因集的藍(lán)-粉聚集圖

從圖中可以發(fā)現(xiàn)提取的特征基因RTN3、ENC1等45個基因在阿爾茨海默癥數(shù)據(jù)集的無病樣本和有病樣本中明顯的聚類，在生物學(xué)方面驗(yàn)證了文中提出的算法提取特征基因的有效性。

4 結(jié)束語

通過以上實(shí)驗(yàn)和結(jié)果可以看出，小波包變換-SAM對基因表達(dá)譜數(shù)據(jù)具有明顯的降噪作用，通過FastICA算法提取的特異性基因具有較高的分類準(zhǔn)確性。提取的特異性基因通過基因富集性分析能夠在沒有先驗(yàn)知識的情況下也能在表達(dá)譜整體層次上對數(shù)條基因進(jìn)行分析，從而從數(shù)理統(tǒng)計上把阿爾茨海默癥的基因表達(dá)譜數(shù)據(jù)與生物學(xué)意義很好地銜接起來，這將會為疾病的診斷、治療等方面提供重要參考。

[1]樂奕勤.阿爾茨海默氏病的發(fā)現(xiàn)者[J].Digest of科技文摘，2000（7）:132.LE Yi-qin.The discover of Alzheimer’s disease[J].Digest of Science and Technology，2000（7）:132.

[2]Atlas S W.Magnetic Resonance Imaging of the brain and spine[M].3rd ed.Philadelphia:Lippincott Williams and wilkins，2002.

[3]劉如云，蔡立君，易葉青.基于G-ICA的組織樣本分類算法[J].計算機(jī)工程與應(yīng)用，2010，46（31）:124-126，156.LIU Ru-yun，CAI Li-jun，YI Ye-qing.Classification algorithm of the tissue samples based on G-ICA [J]. Computer Engineering and Applications，2010，46（31）:124-126，156.

[4]Chen Y，Kamat V，Doughertyer R， et al.Ratio statistics of gene expression levels and applications to microarray data analysis[J].Bioinformatics，2002，18（9）:1207-1215.

[5]Krajweski P，Bocianowski J.Statistical methods for microarray assays[J].J Appl Genet，2002，43（3）:269-278.

[6]Tusher V G，Tibshirani R，Chu G，Significance analysis of microarrays applied to the ionizing radiation response[J].PNAS， 2001，98（9）:5116-5121.

[7]賀憲民，武建虎，賀佳，等.小樣本情況下差異表達(dá)基因鑒別的參數(shù)統(tǒng)計分析[J].中國衛(wèi)生統(tǒng)計，2005，22（3）:141-145.HE Xian-min，WU Jian-hu，HE Jia，et al.Parametric statistical analysis of differentially expressed genes identified in the case of small sample[J].China’s Health Statistics，2005，22（3）:141-145.

[8]WEI Kong.Study DNA microarray gene expression data of Alzheimer’s disease by independent component analysis[C]//International Joint Conference on Bioinformatics，Systems Biology and Intelligent Computing，Shanghai:August，2009.[9]黃得雙.基因表達(dá)譜數(shù)據(jù)挖掘方法研究[M].北京：科學(xué)出版社，2009.

[10]羅亮，史曉紅，徐進(jìn).LVQ神經(jīng)網(wǎng)絡(luò)方法預(yù)測蛋白質(zhì)結(jié)構(gòu)的二硫鍵[J].系統(tǒng)仿真學(xué)報，2007，19（9）:2077-2079.LUO Liang， SHIXiao-hong， XU Jin.Predictprotein structure of the disulfide bond based on LVQ neural network method[J].System Simulation Journal，2007，19（9）:2077-2079.

[11]Subramanian A，Tamayo P，Mootha V K.Gene set enrichiment analysis:A knowledge-based approach for interpreting genome wide expression profiles[J].PNAS，2005（102）:15545-15550.