王 薇, 蔣俊正
(桂林電子科技大學 信息與通信學院,廣西 桂林 541004)
如今,癌癥的病發率越來越高,根據國際癌癥研究機構最新發表的《2020年全球癌癥統計報告》,2020年全球新發惡性腫瘤1 930萬例,全球惡性腫瘤死亡病例約為1 000萬例,癌癥已經成為了全球第二大死因[1]。隨著計算機技術和基因芯片技術的發展,海量的基因表達數據可以被獲取,這為致癌機理的研究和癌癥的早期診斷提供了一種全新的方法。然而,相對于可以大量獲取的基因表達數據,可獲得的用于分析基因數據特性的樣本較少,且僅有一小部分的基因與癌癥分類識別相關。因此,如何降低基因表達數據的維度,從成千上萬的基因表達數據中篩選出分類能力高的基因,對基因表達數據的研究具有重要意義[2-3]。
特征選擇是常用的降維方法之一,其主要思想是基于一定準則,篩選出符合要求的特征[4]。目前,許多學者已將特征選擇應用到基因表達數據的降維中。Chris等[5]提出了最小冗余最大相關(minimum redundancy-maximum relevance, 簡稱MRMR)方法,并用于基因選擇,該方法能有效減少冗余,并具有較好的泛化能力。Sun等[2]利用多個不同的濾波器以選擇出具有不同分類能力的候選基因子集,然后使用交叉熵算法剔除候選基因子集中的冗余基因。Li等[6]將基因選擇問題建模為一個多元優化問題,提出了一種多目標排序二元人工蜂群算法,對基因表達數據進行降維。然而,在基因表達過程中,一個基因會受到其他基因或者分子(如蛋白質)的調控作用,并且這個基因又可能調控其他基因的表達[7]。但是,現有的特征選擇方法在對基因表達數據進行降維時,大多未考慮基因之間的調控作用,僅僅依據特定的準則來選取基因,通過這種方式選擇出的基因很可能缺失基因之間調控作用的信息,不利于癌癥的致癌機理研究。因此,如何對特征選擇方法進行改進,使其在篩選基因時可以納入基因之間的調控作用,是一個值得研究的問題。
基因之間的調控作用可通過推斷基因調控網絡來實現。利用基因表達數據推斷基因調控網絡的模型主要有3類:貝葉斯網絡模型、布爾網絡模型、常微分方程模型[7-9]。由于基因在表達時受不同的基因調控,基因調控網絡的結構是非規則的。隨著信號處理理論的發展,圖信號處理(graph signal processing,簡稱GSP)理論被用來處理非規則網絡結構上的數據。此外,傳統信號處理中一些重要的理論分析概念和方法也被拓展到圖信號處理理論中,并衍生出了圖傅里葉變換(graph Fourier transform,簡稱GFT)[10]、圖采樣[11]、頻域分析[12]、圖濾波器[13]、圖濾波器組[13-15]等圖信號處理理論分析工具。在圖信號處理中,圖模型可以刻畫節點之間的相似性等關系,因此基因間的調控作用可通過圖模型來刻畫。若將基因作為圖上節點,推斷出基因調控網絡可建模為圖模型中的鄰接矩陣,那么在選擇具有高分類能力的基因時,可利用圖濾波器等理論工具構建新的特征選擇方法,使其不僅能選擇出分類能力高的基因,也能保留基因表達時的調控關系。
為了在篩選分類能力強的基因時保留基因間的調控關系,提出了一種基于C3NET(conservative causal core,簡稱C3NET)和圖濾波器的基因特征選擇算法。基因之間的調控關系可用基因調控網絡表示,因此利用C3NET算法推斷出基因表達數據集的基因調控網絡,以此保留基因間的調控關系。此外,相較于冗余基因,分類能力高的基因在不同的樣本(病人)類別中的基因表達值差異較大。依據圖信號處理理論,這些分類能力高的基因所攜帶的信息在圖頻域呈現高頻特性,可利用高通濾波器得到分類能力高的基因。因此,將基因調控網絡建模為圖上鄰接矩陣,同時將基因表達數據中的基因建模為圖上節點,將基因數據建模為圖信號。通過設計高通圖濾波器對基因數據進行濾波,并根據圖傅里葉變換提出了一種評估基因分類能力的指標Duk,篩選出濾波后分類能力高的基因。仿真實驗結果表明,本算法篩選出的基因在不同的分類器中的分類準確率優于其他對比的算法。

給定一個無向圖G=(V,E,A),圖信號f=[f1,f2,…,fN]T是由圖上每個節點的信號值組成的N維向量。與傳統信號處理類似,圖信號處理理論中也有相應的圖傅里葉變換和圖濾波器。對于圖信號f,其圖傅里葉變換GFT定義為
(1)
其中,U=[u1,u2,…,uN]為L的特征向量矩陣。同時,逆-圖傅里葉變換(inverse graph Fourier transform,簡稱iGFT)定義為

(2)
圖濾波器有多種形式,常用的圖濾波器為圖拉普拉斯矩陣多項式的形式,即
(3)
其中,K為濾波器的階數。由于Lk與L的特征向量相同,所以式(3)可寫為
(4)

一個基因表達數據集S∈RM×Q包含了M個樣本(病人)的Q個基因數據,每個樣本(病人)的類別信息用y∈RM×1表示,若數據集S包含了P個類別,則y的取值范圍為[1,2,…,P]。
將樣本(病人)的每個基因建模為圖上節點,得到數據集S的基因圖模型GS=(VS,ES,AS),其中:VS={1,2,…,Q}為Q個節點(基因)集合,AS∈RQ×Q為圖GS的鄰接矩陣,ES為圖GS上節點(基因)邊的集合。目前,構建圖模型的算法有很多種,如最近鄰算法[16]、基于圖學習的方法[17]等。圖1為利用最近鄰算法將Gastric數據集建模為圖的示意圖。
對于基因圖模型GS=(VS,ES,AS),將每個樣本(病人)的基因數據建模為圖信號,即xi=S(i,∶),i=1,2,…,M,其中,S(i,∶)為基因表達數據集S的第i行,則可得到M個圖信號。圖2為將Gastric數據集的第一個樣本(病人)x1建模為圖信號的示意圖。

圖1 Gastric數據集圖模型構建示意圖

圖2 Gastric數據集圖信號建模示意圖
基因表達數據具有高維度、小樣本、冗余多的特點,如何從基因表達數據中篩選出分類能力較高的基因,對癌癥的早期診斷具有重要意義。為此,提出了一種基于C3NET和圖濾波器的基因特征選擇算法。首先,利用C3NET算法推斷基因表達數據集的基因調控網絡;其次,將得到的基因調控網絡作為基因圖GS的鄰接矩陣,計算其圖拉普拉斯矩陣和圖傅里葉變換,并定義評估基因分類能力的指標Duk;然后通過圖濾波器對基因表達數據進行濾波;最后根據Duk的大小篩選出分類能力高的基因。
C3NET是一種基于統計推斷構建基因調控網絡的算法,其基本思想是當2個基因之間的互信息至少是其中一個基因與其他所有基因的互信息的最大值時,2個基因才能相互連接[9]。C3NET算法主要分為2步。
第1步,計算基因之間互信息。對于2個隨機變量X、Y,其互信息計算方式為
(5)
其中:ΩX、ΩY為隨機變量X、Y的樣本空間;p(x,y)為X,Y的聯合概率分布;p(x)、p(y)分別為X、Y的邊緣分布。實際上,互信息需要通過合適的估計器從數據中估計出來,使其接近總體的理論值,因此,使用參數高斯估計器估計互信息的值,
(6)

因此,對于一個基因表達數據集S∈RM×Q,利用式(6)計算任意2個基因之間的互信息,可得到一個互信息矩陣I∈RQ×Q,其元素I(i,j)表示第i和第j個基因的互信息值。
第2步,根據得到的互信息矩陣I,確定每個基因與其相連的鄰居,以得到基因調控網絡。給定一個全連接矩陣C∈RQ×Q,其所有元素C(i,j)均為1,給定閾值α,若I(i,j)<α,則令C(i,j)=0,將第i個基因的鄰居集合定義為
NS(i)={j∶C(i,j)=1且j≠i},
(7)
則基因表達數據集S的基因調控網絡W∈RQ×Q可通過式(8)得到:
W(i,jc(i))=W(jc(i),i)=
(8)
其中,i=1,2,…,Q。
C3NET具體算法流程為
算法1C3NET算法
輸入:全1矩陣C∈RQ×Q,全零矩陣W∈RQ×Q,閾值α;
輸出:基因調控網絡W;
根據式(7)計算任意2個基因之間的互信息,得到互信息矩陣I∈RQ×Q;
fori=1:Q
forj=1:Q
ifI(i,j)<α
C(i,j)=0
end
end
end;
fori=1:Q
Ns(i)={j:C(i,j)=1且j≠i}
ifNS(i)≠?
else
jc(i)=?
end
W(i,jc(i))=W(jc(i),i)=1
end;
return 基因調控網絡W
首先,將C3NET算法得到的基因調控網絡W建模為基因表達數據集S的基因圖模型GS的鄰接矩陣,即令AS=W。其次,計算拉普拉斯矩陣LS,并對LS進行特征分解,即LS=USΛSUST。然后,利用式(1)計算每個圖信號的圖傅里葉變換,得到所有病人的圖傅里葉變換,并將其拼成矩陣

為了確定基因的分類能力,定義了計算基因分類能力的公式:
(9)


(10)
在圖信號處理理論中,若圖上節點的圖信號與其鄰居節點上的信號差異較大,則認為這是一個高頻信號。與傳統信號處理類似,若要保留圖信號的高頻信息,則需設計高頻圖濾波器對圖信號進行濾波。由式(4)可知,圖濾波器的頻率響應為
(11)


(12)
因此,高通圖濾波器的濾波器系數c可以通過求解如下優化問題得到:
(13)
其中:h=[λ0,λ1,…,λK]T;c=[c0,c1,…,cK]T。
在基因表達數據集中,具有較高分類能力的基因在不同樣本(病人)類別中基因表達值差異較大,而分類能力較低的基因表達值在不同的類別中可能不會表現出明顯的差異。因此,認為具有較高分類能力的基因是高頻圖信號,利用高通圖濾波器過濾出分類能力高的基因。
本算法主要分為3步:第1步,利用C3NET算法推斷基因表達數據集的基因調控網絡;第2步,利用高通圖濾波器對基因表達數據集進行濾波;第3步,根據Duk來篩選分類能力高的基因,基于C3NET和圖濾波器的基因選擇算法流程如下所示。
算法2基于C3NET和圖濾波器的基因選擇算法
輸入:基因表達數據集S∈RM×Q,圖濾波器的階數K,篩選出的基因個數d;
輸出:篩選后的基因表達數據集Y∈RM×d;
根據算法1得到基因表達數據集S∈RM×Q的基因調控網絡W∈RQ×Q;
令AS=W;
計算歸一化圖拉普拉斯矩陣Lnor=I-D-1/2AD-1/2,并對其進行特征分解Lnor=UΛUT;

利用式(9)計算Duk,k=1,2,…,Q,篩選出最大的d個Duk值,并記錄其序號k,組成集合index;
whilei≤ddo

end while
使用從公開數據庫[18]中下載的2個基因表達數據集驗證算法的有效性。數據集的具體信息如表1所示。

表1 基因數據集信息
在仿真實驗中,分別用MATLAB自帶的5個分類器native Bayes(NB)、SVM、random forest(RF)、KNN和discriminant analysis(DA)評估篩選出的基因子集的分類準確率。首先,為了驗證圖濾波器在基因選擇時的有效性,分別對2種情況進行仿真:1)根據Duk大小,選擇數據集中Duk最大的20個值所對應的基因;2)先使用高通圖濾波器對數據集進行濾波,然后再篩選出Duk最大的20個值所對應的基因,其中,Gastric數據集中使用的高通圖濾波器的階數為K=5,Brain tumor數據集使用的高通圖濾波器的階數為K=10,實驗結果如圖3所示。

圖3 經過濾波和未經濾波分類準確率對比
由圖3可知,經過濾波后,根據Duk選擇基因的分類準確率明顯高于直接使用Duk,這是因為經過濾波后的圖信號融合了其K階鄰居的信息,也就是說,經過濾波的基因,不僅包含其自身的數據信息,也包含了與其相關的基因的調控信息。如果一個基因已經發生了癌變,那么受其影響的基因也存在異常表達的可能性,而濾波的操作將這些信息聚合到已經發生癌變的基因上,使得發生癌變的基因與正常表達的基因的差異更加明顯,區分度更高。因此,經過濾波后的基因的分類能力更高。
為了驗證本算法的有效性,將本算法分別與NNG-GS[19]、MDG-GS[19]、CDG-GS[19]、LLE[20]、PCA[21]算法進行比較,表2、3分別為從Gastric數據集和Brain tumor數據集中篩選d=20個基因的實驗結果。對于Gastric數據集,本算法的參數α=0.7,高通圖濾波器的階數K=5,對于Brain tumor數據集,本算法的參數α=0.7,高通圖濾波器的階數K=10,NNG-GS、MDG-GS、CDG-GS、LLE四種算法在2個數據集中的參數k均設置為k=3。
根據表2、3的仿真實驗結果,對于Gastric數據集,本算法雖然在DA分類器中的分類準確率與NNG-GS算法相同,但在其他4個分類器中的分類準確率均是最高的。對于Brain tumor數據集,除了KNN分類器,本算法在其他4個分類器中的分類準確率也是最高的。此外,雖然在KNN分類器中的分類準確率較低,但本算法在5個分類器中的平均分類準確率為94.71%,不僅高于CDG-GS算法在KNN分類器中的分類準確率,而且在所有對比算法中平均分類準確率也是最高的。由于基因在表達時受不同基因調控,本算法利用高通圖濾波器對基因數據集進行濾波,濾波后的基因獲得了與其相關的基因調控信息,使得發生癌變的基因與正常表達的基因的差異更加明顯,從圖3的仿真對比結果也可看出,經過高通圖濾波器濾波后的基因在不同分類器中的分類準確率更高。因此,相較于對比算法,高通濾波后的基因融合了更多基因信息,由此篩選出的基因具有更高的分類能力,可以達到較好的分類效果。

表2 Gastric數據集的分類準確率 %

表3 Brain tumor數據集的分類準確率 %
基因表達數據具有高維度、小樣本、冗余多的特點,為了篩選出分類能力高的基因,提出了一種基于C3NET算法和圖濾波器的基因特征選擇算法。由于基因調控網絡有助于分析基因間的調控作用,采用了C3NET算法推斷基因表達數據集的基因調控網絡。同時,將圖信號處理應用到基因的特征選擇問題中,將基因調控網絡建模為圖上鄰接矩陣,基因表達數據建模為圖信號,利用高通圖濾波器對基因數據進行濾波,使得濾波后的基因融合了更多調控信息,提高了基因的分類能力,并提出了一種計算基因分類能力的方法。仿真實驗結果表明,與對比算法相比,本算法篩選出的基因在不同分類器中分類準確率更高,區分能力更強。然而,本算法在KNN分類器中的分類準確率有待提高,接下來將考慮設計更合理的基因調控網絡以提高算法的性能。