基于最大最小爬山算法的肺癌預后模型

(山東科技大學 數學與系統科學學院，山東 青島 266590)

肺癌是發病率和死亡率增長最快、對人類健康和生命威脅最大的惡性腫瘤之一，世界衛生組織國際癌癥研究機構發布最新報告稱肺癌死亡人數最多，占預計癌癥死亡總人數的18.4%。另外，肺癌還具有預后差的特點，影響其預后的因素主要包括患者個體相關因素、腫瘤相關因素和治療相關因素[1]。目前，臨床醫學主要根據手術病理分期判斷預后，由于考慮影響腫瘤發生的因素減少，其預測效果較差[2]，因此，建立一個適用于臨床醫學且考慮多因素的肺癌預后模型具有重要意義。

早期，國內外在疾病預測方面多采用統計學中的COX回歸方法構建模型。隨著數據挖掘技術被應用到醫學研究領域，眾多學者采用機器學習的方法進行疾病研究。劉雅琴等[3]使用logistic回歸、決策樹和人工神經網絡方法研究預后模型的預測效果，是國內機器學習領域研究腫瘤預測的有效嘗試。Kim等[4]利用支持向量機預測了乳腺癌患者術后5年生存情況。Chen等[5]對4個醫療機構的非小細胞肺癌患者，使用人工神經網絡建立了患者生存狀況風險模型。牟冬梅等[6]構建了妊娠高血壓綜合征危險因素決策樹預測模型。宋一鳴[7]基于SEER數據庫使用決策樹、神經網絡、支持向量機、Logistic回歸、深度神經網絡等分別建立了肺癌患者預后的相關研究模型。

復發、轉移、風險評估及生存情況評價是腫瘤預后的主要研究內容[8]，本研究針對患者術后5年后生存情況進行研究。選取SEER數據庫[9]中部分肺癌患者的數據，根據相關研究提取影響患者生存情況的預后因素，通過貝葉斯網絡方法利用訓練集構建肺癌預后模型，其中采用最大最小爬山算法建立模型，并采用貝葉斯估計進行模型參數學習與概率推理，最后將本研究模型與Logistic回歸、人工神經網絡、決策樹及支持向量機方法在測試集上進行分類實驗比較，驗證所建立模型的有效性。

1 數據來源及變量選擇

數據選自美國國立癌癥研究所“監測、流行病學和結果數據庫”(SEER數據庫)[9]中2008年至2014年期間被確診為肺癌的患者，其中包括5年內直接因癌細胞致死和隨訪期滿5年仍生存的患者。刪除數據缺失嚴重、記錄錯誤及因非肺癌致死的患者記錄，最終共計879位患者數據。

表1 肺癌患者變量信息

根據腫瘤信息，參考文獻[10,11]和其他相關研究[2-3,7]中提及的與患者生存相關的預后因素，從數據庫中導出包含這些因素的16個信息變量，具體如表1所示，其中后四項為連續型變量，其余為離散型變量。

5年后生存情況是預后效果的重要評價指標，所以選擇患者術后5年生存情況為結果變量(生存時間以月為單位)。生存時間60個月及以上患者生存情況為“生存”(記為1)，低于60個月的患者生存情況為“死亡”(記為0)。

2 肺癌預后模型的建立

2.1 特征選擇

為提高模型的預測準確性，對上述16個信息變量進行特征選擇。首先，利用SPSS進行卡方檢驗，在p<0.05下通過檢驗的變量有12個，分別為：婚姻狀況、組織學分級、腫瘤分期、轉移程度、擴散程度、淋巴結累積程度、手術類型、是否放療、確診年齡、腫瘤大小、淋巴結受檢數量及淋巴結陽性數量。然后，在卡方檢驗基礎上利用SPSS進行Logistic回歸分析，在p<0.05下最終篩選出的特征變量有6個，分別為：組織學分級、腫瘤分期、確診年齡、腫瘤大小、淋巴結受檢數量及淋巴結陽性數量。篩選結果如表2所示。

2.2 數據離散化

(1)

表2 Logistic回歸分析篩選變量結果

2.3 模型建立方法

在疾病生存預測方面，傳統的統計模型難以計算后驗概率，不能直觀地表示變量之間的關系，本研究利用貝葉斯網絡方法建立肺癌預后模型。

貝葉斯網絡是一個帶參數的有向無環圖，用二元組〈G,Θ〉表示，其中G=(V,E)表示節點關系的有向無環圖，稱為貝葉斯網絡結構，節點集合V={X1,X2,…,Xn}表示隨機變量，有向邊集合E={eij|Xi→Xj,i,j=1,2,…,n}表示變量之間的依賴關系；Θ={Θ1,Θ2,…,Θn}表示節點Xi的條件概率，稱為貝葉斯網絡參數，節點Xi的參數Θi表示其自身和父節點集Pa(Xi)的條件概率分布，即Θi=P(Xi|Pa(Xi))。另外，任意給定的貝葉斯網絡都滿足馬爾科夫條件，即?Xi∈V，Xi獨立于除其父節點集合Pa(Xi)外的所有非子孫節點，因此，變量集V=(X1,X2,…,Xn}聯合概率分布可分解為：

(2)

貝葉斯網絡模型用有向無環圖表示變量之間的依賴和獨立關系，用條件概率分布刻畫變量對其父節點的依賴關系，因此，建立貝葉斯網絡模型包括兩部分：①確定變量間關系，找到網絡結構，即結構學習；②確定每個節點的條件概率表，即參數學習。

2.3.1 結構學習方法

利用最大最小爬山(Max-Min hill-climbing, MMHC)算法對貝葉斯網絡結構進行學習。該算法是Tsamardinos等[12]于2006年提出的一種經典的貝葉斯網絡結構學習算法，結合了依賴分析和評分搜索等方法，分為兩個階段進行學習：第一階段利用MMPC(max-min parents and children)算法確定出每個節點的候選父子節點集，構建出貝葉斯網絡結構的無向框架；第二階段利用貪婪爬山算法對已經得到的網絡結構的框架進行搜索評分，找出使評分函數最大的網絡結構。

MMPC算法是從給定數據集中利用最大-最小啟發式策略確定目標變量T的候選父子節點(candidate parents and children，CPC)集，分為兩個階段。第一階段通過定義一個關聯度函數來確定其他變量與目標變量T在給定CPC下的條件依賴程度，函數值越大表示變量間的條件依賴關系越強；當函數值為零時，表示變量間沒有依賴關系，也就是條件獨立。最大最小啟發式策略每次選擇與目標變量T在給定CPC條件下最小關聯度最大的那個變量進入CPC，當除了CPC中變量所有其他變量與目標變量T在給定CPC條件下都條件獨立時，第一階段停止。第二階段檢驗候選父子節點集CPC中的變量，移去不該有變量，即對于CPC中的變量X，如果存在CPC的子集S使得Assoc(X,T|S)，則將變量X從CPC中移去。

變量X與T在給定變量集Z下的關聯度函數定義為：

(3)

(4)

其中，S表示變量集Z的子集。MMPC算法如下所示：

算法1：MMPC算法

輸入：目標變量T，數據集D

輸出：目標變量T的候選父子節點集CPC

第一階段：

1：令CPC≠φ；

2：WhileCPC不再變化 do

3： 〈F,assocF〉=MaxMinHeuristic(T,CPC)

4： ifassocF≠0 then

5：CPC=CPC∪F

6： end if

7：end

第二階段：

8：for 任意X∈CPC

9： if 存在S?CPC，使Assoc(X,T|S)=0，即Ind(X,T|S)then

10：CPC=CPC{X}

11： end if

12：end for

13：返回CPC

子程序MaxMinHeuristic(T,CPC)

輸入：目標變量T，CPC子集

輸出：以CPC為條件集，與T的最小關聯度最大的變量

14：assocF=maxX∈VMinAssoc(X,T|CPC)

15：F=arg maxX∈VMinAssoc(X,T|CPC)

16：返回 〈F,assocF〉

MMHC算法第二階段利用貪婪爬山搜索在結構空間中搜索評分最高的網絡結構，評分函數采用BDeu評分。該階段的貪婪爬山搜索從空圖開始，每一步搜索的過程是：首先在不產生有向環的情況下，對當前所得模型分別執行一次加邊、減邊、轉邊操作得到一系列候選模型，并計算出每個候選模型的評分；然后將最大評分的候選模型與當前模型比較，若最大評分的候選模型評分大，則將其作為下一個當前模型繼續搜索，否則停止搜索并返回當前模型[13]。

在MMHC算法中，貪婪爬山搜索將每個節點的搜索空間限制在其候選父子節點集上，即僅考慮當Y∈CPCX時添加邊Y→X，此約束顯著降低了搜索空間的復雜性，提高了算法的效率。MMHC算法如下：

算法2：MMHC算法

輸入：數據集D

輸出：有向無環圖

1：for 所有變量X∈Vdo

2：CPCX=MMPC(X,D)

3：end for

4：從空圖出發執行貪婪爬山搜索的3個搜索算子加邊、減邊和轉邊。

當且僅當Y∈CPCX時，添加有向邊Y→X。

5：返回最高得分的有向無環圖

2.3.2 參數學習方法

參數學習在統計學中主要有最大似然估計和貝葉斯估計兩種基本方法，本研究采用貝葉斯估計[13]對貝葉斯網絡參數進行學習。

設一個貝葉斯網絡有n個節點V={X1,X2,…,Xn}，其中節點Xi有ri種取值，其父節點π(Xi)的取法有qi種組合。若Xi無父節點，則qi=1。該貝葉斯網絡的參數為：

θijk=P(Xi=k|π(Xi)=j)(i=1,2,…,n;j=1,2,…,qi;k=1,2,…,ri)。

(5)

用θ表示所有θijk組成的參數向量。設D={D1,D2,…,Dm}是一組關于貝葉斯網絡的獨立同分布的完整數據，則θ的似然函數為：

(6)

其中Nijk表示數據集D中滿足Xi=k和π(Xi)=j的樣本數量。假設參數θ的先驗概率分布服從狄利克雷分布Dir(αij1,αij2,…,αijri)(i=1,2,…n;j=1,2,…,qi)，則：

(7)

(8)

從而，p(θ|D)～Dir(Nij1+αij1,Nij2+αij2,…,Nijri+αijri)(i=1,2,…,n;j=1,2,…,qi)，因此，參數θ的貝葉斯估計為[13]：

(9)

2.4 模型建立與結果分析

將最終保留的879條完整觀測記錄的數據集按照7∶3的比例分為訓練集和測試集，其中訓練集樣本為615個，測試集樣本為264個。訓練集用來構建預后模型，測試集用來預測性能，對預后模型進行評價。

實驗環境基本配置為CPU 2.53 GHz、RAM 2.00 GB，操作系統為Windows 7，在MATLAB 7.0上利用貝葉斯網絡工具箱Full BNT1.0.4。對貝葉斯網絡結構的學習，利用MATLAB編程，最終得到的肺癌預后模型如圖1所示，其中，7個節點為表2所示的6個特征變量及1個結果變量，節點之間的連線表明變量間的相互影響關系。實驗結果顯示，腫瘤大小和組織學分級通過影響腫瘤分期間接地影響患者的生存情況；而確診時的年齡、腫瘤分期、淋巴結受檢數量以及淋巴結陽性數量直接影響患者的生存情況，這一結論符合醫學實際。

1-確診時年齡；2-腫瘤大小；3-組織學分級；4-腫瘤分期；5-淋巴結受檢數量；6-淋巴結陽性數量；7-生存情況。

進一步，對圖1得到的預后模型進行貝葉斯網絡參數學習與推理，利用測試集實現對患者生存情況的預測，從而評價該模型的性能。貝葉斯網絡的參數學習與推理過程均利用貝葉斯網絡工具箱FullBNT-1.0.4實現。最終實驗結果顯示在264個測試集樣本中，預測正確的達202例，預測準確率為76.52%，表明由MMHC算法構建的肺癌預后模型對肺癌患者5年后生存情況的預測準確性良好，可以用于對肺癌患者生存情況的預測。

3 對比試驗

在疾病預測方面，目前常用的有Logistic回歸、人工神經網絡、決策樹及支持向量機等機器學習方法[7]。為了進一步研究MMHC算法構建的貝葉斯網絡預后模型的優良性，以預測準確率為標準，將本模型與Logistic回歸、人工神經網絡、決策樹及支持向量機等方法在測試集上進行分類實驗比較。具體地在WEKA[14]上選擇上述四種方法對應的Logistic、J48、Multilayer Perceptron及SMO四個算法，采用十折交叉驗證的方法對測試集數據進行分類，與本算法在預測準確率及其他性能指標方面作比較，結果如表3所示。

由表3可知，提出的預后模型在預測準確率、精確度和ROC曲線下面積的結果均好于其他方法，說明在本研究的肺癌數據上貝葉斯網絡模型是最優的。傳統的疾病預后模型以統計學中COX回歸、Logistic回歸為主，但統計學方法通常要求變量之間滿足獨立性等條件，無法處理變量間共線性的問題，因此存在局限性。貝葉斯網絡模型是一種概率圖模型，通過有向邊和條件概率形象地刻畫出變量間的依賴關系，能夠進行有效地概率推理且預測準確率高，可以應用于疾病預測。

表3 不同算法的預測準確率及性能指標

4 結論

利用貝葉斯網絡方法建立肺癌預后模型，對患者術后5年生存情況進行研究。首先對變量進行特征選擇，最終選擇影響患者生存情況的6項預后因素；然后利用MMHC算法在訓練集上建立肺癌預后模型，在測試集上對患者進行5年后生存情況預測。實驗結果顯示，利用MMHC算法建立的肺癌預后模型的預測準確率達76.52%，高于目前常用的Logistic回歸、人工神經網絡、決策樹及支持向量機方法。但是，本研究未對數據集中所有的變量進行研究，只是根據腫瘤信息文獻提取了與生存預測相關的16個變量，故研究的模型變量具有一定的主觀性與局限性。在未來的研究中，可以對更多的變量進行系統的研究，提高模型的準確性。