肺癌關鍵基因篩選及預測模型的研究

潘以紅, 范雪萌, 朱 平

(江南大學 理學院， 無錫 214122)

肺癌是世界上最常見的惡性腫瘤之一，已成為城市人口惡性腫瘤死亡原因的第一位[1]。肺癌是一種異質性疾病，基因突變、細胞環境和生活環境的變化都可能影響腫瘤的形成、生長和轉移[2]。在過去的幾十年里，人們應用了手術切除、化療和放射治療在內的多模式治療策略，但肺癌患者的預后仍然不令人滿意[3]。由于人們對肺癌潛在發展的分子機制缺乏精確了解，不能夠對晚期肺癌進行有效治療。因此，迫切需要新的生物標志物來用于診斷、預后和藥物反應，以開發有效的治療方法來改善肺癌的臨床結果。

近年來，許多研究表明有多個基因和多種細胞通路參與了肺癌的發生和發展。例如，HMGN5(High mobility group nucleosome 5)在肺癌組織中過表達，其表達程度與肺癌組織的分化程度和TNM分期有關[4]。ROR2(Receptor tyrosine kinase-like orphan receptor2)或Wnt5a的過表達和肺癌預后不良顯著相關，對ROR2和Wnt5a表達水平的檢測有助于肺癌的預后[5]。肺癌組織中HMGA2(High-mobility group AT-hook 2)、Tiam1(T-lymphoma invasion and metastasis1)、Notch1(Notch homolog 1)的表達量顯著升高，高表達的HMGA2能夠促進癌細胞的上皮間質轉化、Tiam1能夠促進細胞侵襲、Notch1能夠促進細胞增殖[6]。在肺癌的不同TNM分期(IA、IB、IIA、IIB、IIIA、IIIB、IV)中，細胞周期通路(Cell cycle)、孕酮介導的卵母細胞成熟通路(Progesterone-mediated oocyte maturation)和卵母細胞減數分裂通路(Oocyte meiosis)都同時被激活，這些通路可能是肺癌診斷和治療的潛在標記[7]，等等。盡管這些研究已經確定了肺癌的一些預測基因和通路，但還不足以為肺癌的治療提供更有針對性的方案，因此需要更多的預測生物標志物。

某些基因表達水平的改變會導致含有這些基因的通路功能異常，同時也會對其他通路功能產生影響，促使癌癥的發生和發展。串話基因是指同時出現在兩條或多條通路中的基因，這些基因將不同的通路關聯起來，可以將某條通路的異常傳遞給其他通路。研究發現EPAS-1 / HIF-2α與PXR信號通路之間的串話基因是胃癌細胞多藥耐藥的調節因子[8]。SKP2被證實可以作為串話基因，通過影響不同的信號通路，參與癌癥的發展[9-10]。

微陣列是用于分析基因表達的高通量平臺，是醫學腫瘤學領域中很有前景的技術。微陣列具有較廣的臨床應用，在從癌癥的分子診斷到分子分類、從患者分層到預后預測、從新的藥物靶點發現到腫瘤反應的預測等方面都有較好的表現[11-13]。微陣列技術與生物信息學分析工具相結合能夠全面分析肺癌的發生和發展過程中基因的表達變化情況。

GEO(Gene expression omnibus)是一個包含微陣列數據和下一代測序數據的公共數據儲存庫。為探討肺癌潛在預測和治療的新生物標志物，從GEO數據庫篩選和下載了肺癌的微陣列數據，將肺癌細胞的基因表達數據與正常細胞的基因表達數據進行比較，使用logistic回歸和倍數法相結合篩選出差異表達基因。通過對差異表達基因的蛋白質交互(protein-protein interaction, PPI)網絡中的關鍵基因篩選以及關鍵基因的通路富集分析進一步了解肺癌在分子水平上的發展情況。最后利用失調的通路中串話基因建立肺癌的預測模型，對潛在的候選生物標志物進行探索，以用于肺癌的診斷、預后以及作為藥物靶點。

1 數據來源

本文從GEO數據庫下載了基因芯片GSE19188、GSE40791。GSE19188數據集中包含91個肺癌患者樣本和65個正常樣本[14]。GSE40791數據集中包含100個肺癌患者樣本和94個正常樣本[15]。

2 方法

2.1 Logistic回歸篩選差異表達基因

本文使用的是GEO數據庫中經過標準化后的基因表達的矩陣序列數據。首先進行數據整合，剔除沒有相對應基因名的數據；如果多個探針號對應一個基因名，則取多個基因表達值的均值。

本文將數據集中樣本性狀信息設為0-1的M矩陣，0代表正常人樣本，1代表肺癌患者樣本。將樣本性狀信息矩陣M與樣本中基因表達值矩陣N進行逐列分析，來篩選出肺癌樣本和正常樣本中差異表達的基因。因變量為是否患癌，自變量為基因。由于樣本性狀信息矩陣M為0-1型因變量矩陣，因此采用針對 0-1 型因變量的 Logistic回歸模型[16]。

考慮具有n個獨立變量的向量X=(X1,X2,…,Xn)，設條件概率P(Y=1|X)=p為根據觀測量相對于某事件X發生的概率。Logistic回歸模型可以表示為

(1)

事件發生與不發生的概率之比為

(2)

(3)

(4)

可對Y進行線性回歸。

本文是逐列比較，概率方程為：

(5)

式(5)中Xi為第i個基因的表達值，設p為患癌概率，取分界點0.5為是否患癌的分界值。但不知道患癌與否的概率p的準確值，為了方便做回歸運算取區間中值，即p=0.75對應Y=1，p=0.25對應Y=2。

本文使用matlab R2012b實現回歸分析，計算出p值。回歸分析過程中，當進行一次假設時，得到的結果的可靠性較高，但當連續進行檢驗的時候，其結果的可靠性就會大大降低，從而出現假顯著性問題，因此本文采用FDR錯誤控制法對p值進行校正。

本文以校正p值小于0.05及|logFC|>1.5為標準，篩選出數據集中的差異表達基因。

2.2 肺癌特異蛋白互作網絡的構建

STRING數據庫是一個搜尋已知蛋白質之間和預測蛋白質之間相互作用的系統。這種相互作用既包括蛋白質之間直接的物理相互作用，也包括蛋白質之間間接的功能相關性[17]。Cytoscape是一個用于圖形化顯示蛋白質互作網絡數據并進行分析和編輯的軟件平臺。它能夠為網絡添加豐富的注釋信息，并且可以利用自身以及第三方開發的大量功能插件，針對網絡問題進行深入分析[18]。本文使用STRING在線工具，以得分>0.4為標準，選取差異表達基因的互作網絡數據。利用Cytoscape軟件繪制肺癌的特異PPI網絡，并對特異PPI網絡進行拓撲特征分析，包括節點的度分布、接近中心性、聚類系數以及平均最短路徑等。最后使用MCODE插件獲得網絡中節點大于10的重要子模塊，通過DAVID在線工具[19]對最大的子網絡中的差異基因進行GO(Gene Ontology)功能富集分析。Gene Ontology可分為分子功能、生物過程和細胞組成3個部分。蛋白質或者基因可以通過ID對應或者序列注釋的方法找到與之對應的GO號，即功能類別或者細胞定位。

2.3 特異PPI網絡中關鍵基因的篩選

根據特異PPI網絡中節點(基因)的度和基因表達的差異得分篩選出差異表達基因中的關鍵基因。首先定義每個基因的表達區間I:(AVG-AD,AVG+AD)，其中AVG為基因在正常樣本中的表達值均值，AD為基因在正常樣本中表達值的平均差。平均差含義明確，能夠較全面、客觀地反映變量數列標志值平均變動程度[20]。在癌癥樣本中如果某個基因的表達值超出了區間I，我們就把這個超出的數值用于計算差異得分。差異得分計算公式如下:

(6)

W=degree×score

(7)

其中，di表示癌癥樣本i中某基因表達值，若di大于AVG+AD，則d值為AVG+AD;若di小于AVG-AD，那么d值就為AVG-AD。我們用log2函數將節點的度值標準化。最后，計算得出W值，并將差異表達基因排序。

在PPI網絡中，差異得分和度較大的基因和肺癌關系較為密切，即為我們挑選出的關鍵基因。本文中，我們分別挑選了兩個數據集中W值排名前100和后150的基因，并進行重疊之后選出共同的基因作為關鍵基因。

2.4 關鍵通路的富集分析及串話基因的選取

使用DAVID在線工具對篩選出的關鍵基因進行KEGG通路分析。KEGG通路數據庫是一個手工畫的代謝通路的集合，包含以下幾方面的分子間相互作用和反應網絡：新陳代謝、遺傳信息加工、環境信息加工、細胞過程、生物體系統、人類疾病和藥物開發。本文選取P<0.05的通路作為最后結果，根據通路得分對樣本進行分層聚類分析。每個樣本中的通路得分計算公式如下：

(8)

通過斯皮爾曼等級相關系數，得到通路間相關程度，公式如下：

(9)

相關系數在(0，1)之間，則通路正相關;相關系數在(-1, 0)之間，則通路負相關。最后在顯著相關的通路中，挑選出差異表達基因中的串話基因。

2.5 預測模型的建立及驗證

將2.4中得到的串話基因作為預測基因，采用支持向量機(support vector machine,SVM)建立肺癌預測模型。將所得串話基因在數據集GSE19188的156個樣本中表達值作為訓練集，采取五折交叉驗證確定核函數和最優參數C和g，用數據集GSE40791進行測試。最后使用生物信息學中預測(分類)的評價標準：總體準確率(accuracy, ACC)、敏感性(sensitivity, SN)、特異性(specificity, SP)和馬修斯相關系數(Matthew’s correlation coefficient, MCC)來衡量模型的分類結果。ACC、SN、SP、MCC的定義如下：

(10)

(11)

(12)

(13)

其中，TP表示正確預測到的正例個數；TN表示正確預測到的反例個數；FP表示把反例預測成為正例的個數；FN表示把正例預測為反例的個數。

圖1 兩個數據集的重疊分析結果Figure 1 Overlap analysis of two data sets

3 結果及分析

3.1 差異表達基因

從數據集GSE19188和GSE40791中分別篩選出差異基因1020個和1770個。對兩個數據集中差異表達基因做重疊分析篩選出808個相同基因(圖1)，其中，有805個基因在兩個數據集中都呈現了相同的表達趨勢，分別為214個上調基因和591個下調基因。

3.2 肺癌特異PPI網絡分析

本文從STRING 數據庫中得到了805個差異表達基因中的632個基因之間6312對互作關系，利用Cytoscape軟件構建出差異表達基因的特異PPI網絡，即此網絡包含632個節點，6312條邊(圖2)。

特異PPI網絡的節點拓撲特征參數如圖3所示，基因數目以1-632的數字編碼表示。網絡的節點度分布遵循了冪律分布(圖3-A)，大多數節點的度較小，少數節點的度較大，節點度分布集中在0-20及80-100的區域，這種差異可能是由于網絡中存在高密度連接的模塊造成的。節點連接的接近中心性分布在0.2～0.4之間(圖3-B)，表明了特異PPI網絡中節點良好的連通性，差異表達基因之間的功能關系較為密切。網絡高連通蛋白的聚集系數大于0.7(圖3-C)，表明聚集程度較高。平均最短路徑的峰值集中在2.5～4區域(圖3-D)，因此，特異PPI網絡拓撲參數均呈現出一定的聚集區域，表明在網絡中存在多個功能模塊。

圖2 特異蛋白交互網絡Figure 2 Protein-protein interaction network

A：度分布；B：接近中心性；C：聚類系數；D：平均最短路徑

圖3特異PPI網絡的節點拓撲特征參數
Figure 3 Topological properties of protein-protein interaction network

3.3 特異PPI網絡中重要模塊的功能性分析

通過MCODE插件得到了特異PPI網絡中節點大于10的4個重要子模塊，如表1所示。

表1 4個節點數大于10的子模塊Table 1 Four modules with more than 10 nodes

圖4特異蛋白網絡的子網絡
Figure 4 The largest subnetwork of protein-protein interaction network

第一個模塊中包含85個差異表達基因，是特異網絡中最大的子網絡(圖4)。其他模塊分別包含31、31和34個差異表達基因。本文對最大的子網絡做GO富集分析，結果如表2所示。可以看出，子網絡中所含的差異基因，參與多個生物學過程包括細胞分裂、胚胎形成，蛋白結合等，與細胞周期有密切的相關性。

表2 最大子網絡的GO富集分析Table 2 GO enrichment analysis of the largest subnetwork

3.4 關鍵基因的通路分析

對差異表達基因的差異得分和進行log2轉化后的度值做乘積運算，得到W值，并根據W值從兩個數據集挑選出209個關鍵基因。這些基因的表達值與正常樣本相比存在顯著偏差，并且處于特異PPI網絡中的中心位置，與肺癌的發生和發展都有重要關系。

對兩個數據集中W絕對值排名前150個關鍵基因中的123個相同關鍵基因進行KEGG通路富集分析，一共得到P<0.05的11條通路(表3)。其中，上調通路7條，下調通路4條。關鍵基因主要在細胞周期、癌癥、信號以及代謝等通路上富集，這表明了肺癌的發生和發展的潛在機制。細胞周期的紊亂可能導致癌癥的發生，而信號轉導和代謝的異常可能與癌細胞的轉移、疾病的復發和凋亡有關。

表3 關鍵基因的KEGG通路富集分析Table 3 Pathway enrichment of the key genes

隨機選取數據集GSE19188中61個癌癥(Cancer)樣本和50個正常(Normal)樣本，分別計算每個樣本的通路得分得到層次聚類結果(圖5)，表明11條通路的得分可以清楚地對正常樣本和癌癥樣本進行區分。上調通路的變化趨勢顯著，下調通路在某些癌癥樣本中下調趨勢不顯著，可能是由于在癌癥樣本中，這些通路對癌癥的發展還存在一定的抑制作用。

通路間相關系數如表4所示。標注*表明兩條通路相關性不顯著。在相關性顯著的通路中同時存在的基因，我們選取其作為串話基因，這些基因的差異表達暗示其所在通路出現異常。一共得到了18個串話基因，其中，CCNE2在5條通路中出現；MAD2L1、CDK1在4條通路中出現；CHEK1、CDC20、PTTG1、CCNB1、CCNA2、FOS及IL6在3條通路中出現；BUB1B、MCM2、PPARG、MMP1、MMP9、BIRC5、CD36和POLE2在2條通路中出現。

圖5 11條通路的層次聚類分析Figure 5 Hierarchical clustering analysis of 11 pathways

3.5 預測模型建立及評價

輸入數據集GSE19188中18個串話基因的表達值，分別試用不同的核函數來訓練模型，利用數據集GSE40791進行測試，分別統計其ACC、SN、SP與MCC的值(表5)。

表5 不同核函數的預測性能Table 5 Prediction performance of different kernel functions

從表5可以看出，高斯核函數的特異性達到98.94%，說明其能夠很好地分辨出正常樣本；線性核函數和sigmoid核函數的特異性達到99%，說明其對癌癥樣本的識別率高。總體來看，高斯核函數的分類效果最優，準確率達到97%。表明18個串話基因對正常樣本和癌癥樣本具有很高的區分能力，可作為肺癌的預測標記物和化療的靶點。

4討論

肺癌是發病率和死亡率增長最快、對人群健康和生命威脅最大的惡性腫瘤之一。本文從肺癌樣本中得到了805個呈現相同表達趨勢的差異表達基因，其中有214個上調基因和591個下調基因。最終經篩選得到209個與肺癌相關的關鍵基因，其中部分基因已被研究證明與肺癌的發生發展密切相關。比如：KIAA0101是本文特異PPI網絡中一個樞紐基因，度為98。研究發現在肺癌等多種人類惡性腫瘤中都出現了KIAA0101基因過表達的情況，表明KIAA0101表達水平會影響腫瘤的發展[22-23]。TOP2A(DNA Topoisomerase II Alpha)的W值較大，它的表達水平在一定程度上反映了腫瘤細胞的增殖狀態，TOP2A在腫瘤組織中過表達預示著一些腫瘤患者預后不良，同時也會導致惡性腫瘤的淋巴結轉移和遠處轉移[24]。TTK(TTK Protein Kinase)被發現是癌癥治療的新靶點，可以作為一個獨立的預后生物標志物[25]。這些研究結果表明本文得到的還未被證實與肺癌關鍵基因有一定的后續研究價值。

根據關鍵基因所富集的11條細胞通路的差異得分，可以清楚地分辨出癌癥樣本和正常樣本。11條通路中一共有18個串話基因，其中卵母細胞減數分裂通路(Oocyte meiosis)和細胞周期通路(Cell cycle)共同含有5個相同基因：CCNE2、MAD2L1、CDK1、CDC20和PTTG1；HTLV-I感染通路(HTLV-I infection)和細胞周期通路(Cell cycle)共同含有5個相同基因：MAD2L1、CDC20、CHEK1、PTTG1和BUB1B。這3條通路與其他通路相關性都相對顯著，因此可能參與了肺癌發生發展過程，這與文獻[7]中Cell cycle和Oocyte meiosis這兩條通路在肺癌發展的所有時期都被激活的結論是一致的。CCNE2、MAD2L1、CDK1、CDC20、PTTG1和CHEK1這6個串話基因的W值都較大，其中CCNE2、MAD2L1、CDK1和CDC20已被研究證實與肺癌相關。

CCNE2(cyclin E2)在5條通路中出現，它屬于高度保守的細胞周期家族。CCNE1(Cyclins E1)是CCNE2的重要副產物，兩者統稱為統稱為cyclin E。cyclin E被視為S期標志物，在G1/S期決定和限速中起中心調控作用。Cyclin E過表達主要由基因擴增所引起，基因擴增形成大量突變的中心體，從而促進腫瘤的惡性增殖。在許多惡性腫瘤如肺癌、乳腺癌及胃癌等腫瘤中，cyclin E常過度表達并被認為是一種疾病進展和預后的指標。研究發現miR-30d-5p能夠通過下調CCNE2的表達，誘導細胞周期G1/S期阻滯，顯著抑制胞肺癌細胞增殖和運動能力[26]。MAD2L1(Mitotic Arrest Deficient 2 Like 1)在4條通路中出現，它是紡錘體組裝檢查點的關鍵部分。紡錘體組裝檢查點能夠在有絲分裂過程中調控紡錘體微管附著染色體的著絲粒，MAD2L1表達與卵巢漿液性腫瘤的化療抵抗相關；MAD2L1能通過P53直接或間接調控TOP2A和BIRC5的生物學功能，與肺癌的化療效果有很大相關性[27]。CDK1(Cyclin Dependent Kinase 1)也在4條通路中出現，它是周期蛋白依賴性蛋白激酶家族中的一員。CDK1起重要作用的 G2/M 檢驗點主要檢測 DNA 的紡錘體中微管組裝及著絲點位置的正確連接，若 CDK1調節機制發生差錯，將直接導致細胞分化障礙、細胞周期進程發生紊亂，誘導細胞惡性增殖和異常轉化，最終形成惡性腫瘤[28]。CDC20(Cell Division Cycle 20)在3條通路中出現，它是一種調節蛋白，能調節蛋白-蛋白的相互作用。CDC20作為一個促癌因子, 參與了很多腫瘤的發生和發展過程；它能夠調節細胞周期和細胞凋亡過程，因此CDC20 抑制劑是腫瘤治療的一種新策略[29]。

其他串話基因，如IL6(Interleukin 6)，在4條通路中出現。它是一種多功能細胞因子，參與多種細胞的生長、分化和功能調節，在免疫和炎癥反應中具有重要作用，是機體重要的免疫-神經-內分泌調節因子[30]。一項有關吸煙的正常人與吸煙的肺癌患者的基因對比研究中指出IL6是一個重要的差異表達基因[31]。IL6能夠促進肺癌CD133 +細胞的生長和上皮間質轉化[32]。因此IL6可能是與肺癌相關的關鍵基因。

以上討論表明本文從失調的通路間挑選出的18個串話基因具有一定的可信度，可作為肺癌治療的預測因子和潛在靶點。