，， ，， ，，，，

肺癌是我國(guó)常見惡性腫瘤之一，其發(fā)病率和病死率在惡性腫瘤中均居首位。非小細(xì)胞肺癌(non-small cell lung cancer,NSCLC)是肺癌最常見的類型，其中肺腺癌(lung adenocarcinoma,LA)是主要的NSCLC類型，超過40%的患者發(fā)現(xiàn)時(shí)已屬晚期[1]。表皮生長(zhǎng)因子受體(epidermal growth factor receptor，EGFR)、間變淋巴瘤激酶(anaplastic lymphoma kinase,ALK)和鼠肉瘤病毒致癌基因(Kirsten rat sarcoma viral oncogene,KRAS)基因突變是LA中常見的突變類型。近年來，隨著分子靶向治療技術(shù)的不斷進(jìn)步，以吉非替尼、克唑替尼為代表的靶向藥，為EGFR、ALK基因突變陽性的LA患者帶來了新的希望。盡管如此，依然有10%～30% KRAS突變型LA患者無法從現(xiàn)有治療方案中獲益[2]。到目前為止，KRAS基因突變所產(chǎn)生的耐藥機(jī)制仍然不明確。

隨著二代測(cè)序和大數(shù)據(jù)分析方法的不斷進(jìn)步，通過全基因組轉(zhuǎn)錄組測(cè)序獲得的表達(dá)譜能鑒別腫瘤組織與正常組織的差異表達(dá)基因[3]。其中，以基因表達(dá)數(shù)據(jù)庫(gene expression omnibus，GEO)的數(shù)據(jù)量和類型最為豐富，提供了各種類型的轉(zhuǎn)錄組信息[4-5]。本研究擬采用GEO數(shù)據(jù)庫中的NSCLC基因芯片數(shù)據(jù)篩選出KRAS突變型NSCLC差異表達(dá)基因，并對(duì)差異表達(dá)基因進(jìn)行生物信息學(xué)分析，預(yù)測(cè)和篩選出KRAS突變型LA的轉(zhuǎn)錄調(diào)控網(wǎng)絡(luò)，為進(jìn)一步研究KRAS突變型LA的耐藥機(jī)制、早期診斷標(biāo)志和潛在分子靶點(diǎn)提供理論依據(jù)。

1 資料和方法

1.1 資料

1.1.1 數(shù)據(jù)來源 以“KRAS”AND“l(fā)ung cancer”和“KRAS”AND“l(fā)ung”為檢索式在GEO公共數(shù)據(jù)庫進(jìn)行檢索，查看所有“series”。具體數(shù)據(jù)篩選標(biāo)準(zhǔn)如下：①物種為人；②至少2個(gè)生物重復(fù)；③實(shí)驗(yàn)設(shè)計(jì)思路清晰、數(shù)據(jù)質(zhì)量好。在詳細(xì)查看所有數(shù)據(jù)的注釋文件，確定數(shù)據(jù)集中樣本類型及數(shù)量后進(jìn)行分析。本研究使用的肺腺癌KRAS突變數(shù)據(jù)的具體描述如表1，其中，GSE是指一個(gè)實(shí)驗(yàn)項(xiàng)目中的芯片實(shí)驗(yàn)編號(hào)，GPL是芯片平臺(tái)。本研究所使用的數(shù)據(jù)集，主要為Affymetrix芯片數(shù)據(jù)。

表1 肺腺癌數(shù)據(jù)的基本信息

1.1.2 實(shí)驗(yàn)設(shè)計(jì) 本研究所使用的數(shù)據(jù)集均為經(jīng)GEO數(shù)據(jù)庫預(yù)處理后的轉(zhuǎn)錄組表達(dá)譜矩陣。根據(jù)數(shù)據(jù)來源實(shí)驗(yàn)描述以及各樣本類型，將數(shù)據(jù)類型設(shè)定為對(duì)照與突變兩組，即：在未區(qū)分具體人種的情況下，正常的肺組織設(shè)為對(duì)照組，僅具有KRAS突變的LA組織設(shè)為突變組。進(jìn)行差異分析并進(jìn)行后續(xù)的分析研究。

1.2 方法

1.2.1 數(shù)據(jù)處理及差異表達(dá)分析 在Bioconductor(http://www.bioconductor.org/)網(wǎng)站下載生物信息分析的R語言程序包[10]。利用Bioconductor中的Impute程序包對(duì)已獲取表達(dá)譜進(jìn)行歸一化處理，使各個(gè)樣本的數(shù)據(jù)歸一化[11]。同時(shí)利用Bioconductor中的注釋包對(duì)數(shù)據(jù)進(jìn)行注釋，將探針對(duì)應(yīng)到基因上，如果多個(gè)探針對(duì)應(yīng)一個(gè)基因，那么用所有探針表達(dá)值的平均值代表這個(gè)基因的表達(dá)值。得到3組數(shù)據(jù)集的表達(dá)譜用于后續(xù)分析。隨后利用Bionconductor中的Limma程序包對(duì)表達(dá)譜篩選KRAS突變LA樣本與正常樣本間的差異基因。首先使用P<0.05初步篩選差異基因。而后利用Bonferroni算法來校正假發(fā)現(xiàn)率(FDR)的原始P值并計(jì)算倍數(shù)變化(FC)[12]。本研究中，差異基因的篩選標(biāo)準(zhǔn)是FDR<0.05，Log2FC>1。3組數(shù)據(jù)集的差異結(jié)果保存用于后續(xù)分析。

1.2.2 篩選顯著性基因 通過上述分析得到的3組差異基因，隨后利用R語言中的RRA算法篩選3組差異基因中的顯著性基因。RRA方法主要是對(duì)不同數(shù)據(jù)集獲得的基因進(jìn)行排序。如果一個(gè)基因在所有實(shí)驗(yàn)中表達(dá)倍數(shù)越高，且P值越小，則其越顯著[13]。

1.2.3 構(gòu)建顯著性基因的蛋白質(zhì)互作網(wǎng)絡(luò) 利用STRINGv9.1(http://www.string-db.org/)開源數(shù)據(jù)庫構(gòu)建蛋白質(zhì)互作(protein-protein interaction,PPI)網(wǎng)絡(luò)[14]。為更好地理解重要基因的相互作用，使用STRING數(shù)據(jù)庫預(yù)測(cè)其編碼蛋白的PPI網(wǎng)絡(luò)，本研究使用的可靠性閾值為>0.4。Cytoscape軟件(http://cytoscape.org/)用于構(gòu)建PPI網(wǎng)絡(luò)[15]。在PPI網(wǎng)絡(luò)中，每個(gè)節(jié)點(diǎn)代表基因、蛋白質(zhì)或其他分子，節(jié)點(diǎn)之間的連線代表生物分子之間的相互作用。PPI網(wǎng)絡(luò)可用于識(shí)別由EGFR突變型LA中顯著基因編碼的蛋白質(zhì)之間的相互作用和通路關(guān)系。核心節(jié)點(diǎn)，即具有重要生物功能的蛋白質(zhì)通過計(jì)算蛋白質(zhì)之間的連線數(shù)以及每個(gè)節(jié)點(diǎn)與特定節(jié)點(diǎn)有無直接連接來確定，并給出排名。本研究中利用R語言設(shè)計(jì)腳本篩選核心節(jié)點(diǎn)，結(jié)果輸出前10名基因。

1.2.4 富集分析 DAVID(The Database for Annotation，Visualization and Integrated Discovery，http://david.abcc.ncifcrf.gov/)[16]，是生物信息學(xué)常用的聚類分析數(shù)據(jù)庫，其中包含GO(Gene Ontology)[17]和KEGG pathway(The Kyoto Encyclopedia of Genes and Genomes，http://www.genome.jp/kegg/pathway.html)[18]分析。本研究對(duì)出現(xiàn)在PPI網(wǎng)絡(luò)中的56個(gè)顯著性基因進(jìn)行聚類，聚類結(jié)果以P<0.05為閾值進(jìn)行篩選，最后利用R語言對(duì)結(jié)果進(jìn)行可視化展示。

2 結(jié)果

2.1 差異表達(dá)分析及顯著性基因篩選 3個(gè)數(shù)據(jù)集GSE31210、GSE32863、GSE75037分別篩選出差異表達(dá)基因2 625個(gè)、902個(gè)、2 581個(gè)。其中上調(diào)表達(dá)的基因分別為1 115個(gè)、351個(gè)、1 191個(gè)，下調(diào)表達(dá)的基因分別為1 510個(gè)、551個(gè)、1 390個(gè)(圖1A、圖1B、圖1C)。利用RRA算法對(duì)3組差異結(jié)果綜合統(tǒng)計(jì)得到顯著性基因共120個(gè)，其中56個(gè)基因表達(dá)上調(diào)，64個(gè)基因表達(dá)下調(diào)，部分顯著性基因結(jié)果，見圖1D。

注：FDR<0.05，Log2FC>1圖1 各數(shù)據(jù)集差異表達(dá)分析及顯著性基因

2.2 顯著性基因蛋白質(zhì)相互作用調(diào)控網(wǎng)絡(luò)構(gòu)建 用STRING在線軟件對(duì)顯著性基因進(jìn)行PPI網(wǎng)絡(luò)構(gòu)建，其中有56個(gè)蛋白出現(xiàn)在PPI網(wǎng)絡(luò)中，顯示出與其他蛋白的相互作用關(guān)系，基因互作網(wǎng)絡(luò)見圖2A。而后對(duì)互作網(wǎng)絡(luò)進(jìn)一步分析，依據(jù)基因連線數(shù)目篩選核心基因，統(tǒng)計(jì)結(jié)果圖形化展示見圖2B。

圖2 蛋白質(zhì)互作網(wǎng)絡(luò)及核心基因

2.3 蛋白質(zhì)互作網(wǎng)絡(luò)中出現(xiàn)的顯著性基因的富集分析 對(duì)120個(gè)顯著性基因進(jìn)行富集分析，GO分析結(jié)果顯示，這些基因的富集情況分為3種類型，包括生物學(xué)進(jìn)程、細(xì)胞成分、分子功能，圖3。其中，在生物學(xué)進(jìn)程方面，這些基因大多與血管生成、免疫反應(yīng)、細(xì)胞間信號(hào)轉(zhuǎn)導(dǎo)以及細(xì)胞粘附有關(guān)。在細(xì)胞成分方面，大多與藥物吸收、藥物結(jié)合以及耐藥有關(guān)。在分子功能方面，大多與細(xì)胞對(duì)腫瘤壞死因子的反應(yīng)、T細(xì)胞遷移及細(xì)胞對(duì)白介素1的反應(yīng)有關(guān)。

注：橫坐標(biāo)代表GOID，縱坐標(biāo)代表富集在各GO的基因數(shù)目圖3 GO富集分析結(jié)果

KEGG分析結(jié)果顯示這些顯著性基因主要富集在免疫相關(guān)通路、各種感染性疾病以及細(xì)胞周期、凋亡、轉(zhuǎn)錄調(diào)控等細(xì)胞增殖相關(guān)通路和藥物代謝相關(guān)通路上，圖4。

注：分析結(jié)果橫坐標(biāo)代表富集在各通路的基因數(shù)目，縱坐標(biāo)代表通路名稱，P<0.05圖4 KEGG富集分析結(jié)果

3 討論

KRAS突變型LA有難治、預(yù)后差、耐藥性強(qiáng)和容易復(fù)發(fā)等特點(diǎn)，是臨床治療最為棘手的肺癌類型。盡管目前針對(duì)KRAS突變LA有抑制KRAS基因(如法尼基蛋白轉(zhuǎn)移酶抑制劑、香葉基轉(zhuǎn)移酶、RAS轉(zhuǎn)換酶1等)、改變KRAS膜定位(如BKM120、GDC0941和XL147等)等治療方法[19]，但并沒有解決臨床的治療難題。本研究以KRAS突變相關(guān)的全轉(zhuǎn)錄組數(shù)據(jù)庫為研究對(duì)象，通過對(duì)KRAS突變型LA的轉(zhuǎn)錄組學(xué)數(shù)據(jù)進(jìn)行生物信息學(xué)分析發(fā)現(xiàn)，KRAS突變產(chǎn)生的影響不僅限于改變細(xì)胞周期使腫瘤細(xì)胞無限增殖或使EGFR和ALK突變陽性患者出現(xiàn)靶向藥耐藥，實(shí)際上對(duì)整個(gè)基因組均具有不同程度的影響，且這種作用具有指向性。突變影響的基因表達(dá)包括對(duì)細(xì)胞增殖具有顯著影響的細(xì)胞周期以及腫瘤預(yù)后密切相關(guān)的免疫分子等，尤其以IL-6、MMP9、CDH5、TIMP1、TOP2A、TEK、CD36、CLDN5、LCN2、SPP1等基因最為顯著[20-29]，且這些基因與KRAS突變的關(guān)系已有研究證實(shí)。KRAS基因突變后對(duì)LA的預(yù)后可能與其表觀調(diào)控作用有關(guān)，如Caetano等[30]發(fā)現(xiàn)血清IL-6水平與LA存活率低及預(yù)后較差相關(guān)；Xu等[31]發(fā)現(xiàn)MMP9在LA的發(fā)生發(fā)展中發(fā)揮重要作用；Zhang等[29]發(fā)現(xiàn)SSP1上調(diào)PD-1使得LA細(xì)胞發(fā)生免疫逃逸等。以上均說明，采用轉(zhuǎn)錄組數(shù)據(jù)庫分析基因突變的方法學(xué)具有實(shí)用性，且可與DNA測(cè)序等相結(jié)合研究基因突變的病理生理學(xué)意義，指導(dǎo)基礎(chǔ)和臨床研究。

在基因本體方面，本研究發(fā)現(xiàn)，在KRAS突變影響的基因中，大量基因與腫瘤微環(huán)境有關(guān)，例如血管生成、細(xì)胞粘附、免疫反應(yīng)等相關(guān)基因。包括TEK、IL-6、MMP9等在內(nèi)的基因已經(jīng)有明確的文獻(xiàn)證實(shí)其與腫瘤的轉(zhuǎn)移、預(yù)后不良有關(guān)[26,32]。同時(shí)，還有大量基因與藥物反應(yīng)、藥物結(jié)合以及藥物代謝有關(guān)，這提示，KRAS突變所影響的基因中包含與藥物代謝相關(guān)的基因，證實(shí)了KRAS突變極易發(fā)生耐藥的特點(diǎn)，IL-6、CDH3、TOP2A為代表的KRAS突變影響與藥物代謝相關(guān)的分子，有望為未來研究KRAS突變LA的耐藥機(jī)制提供思路[33-34]。

在信號(hào)通路方面，本研究發(fā)現(xiàn)，KRAS突變影響的信號(hào)通路不僅限于KRAS突變后RAS通路的異?；罨?，還包括對(duì)原本具有抗腫瘤作用的免疫相關(guān)通路的抑制，例如對(duì)IL-17信號(hào)通路、白細(xì)胞內(nèi)皮遷移等信號(hào)通路的抑制作用。這2個(gè)信號(hào)通路的抑制，可以改變腫瘤微環(huán)境，出現(xiàn)免疫耐受而導(dǎo)致預(yù)后不良。提示，KRAS突變對(duì)抗腫瘤免疫的抑制，可能是KRAS預(yù)后不良的重要機(jī)制之一。Akbay等[35]發(fā)現(xiàn)IL-17可以通過促炎反應(yīng)促進(jìn)肺癌的發(fā)展，同時(shí)這種促炎反應(yīng)對(duì)PD-1的阻斷有抑制作用，導(dǎo)致肺癌細(xì)胞免疫耐受，這與我們的預(yù)測(cè)結(jié)果是一致的。另一方面，KRAS突變也可以使細(xì)胞周期相關(guān)通路異常活化，增強(qiáng)腫瘤細(xì)胞的增殖能力，同時(shí)抑制腫瘤細(xì)胞的凋亡信號(hào)通路，導(dǎo)致腫瘤細(xì)胞過度增殖。此外，KRAS突變還可以影響與人體藥物敏感性有關(guān)的細(xì)胞色素P450藥物代謝信號(hào)通路，發(fā)生耐藥。

由此可見，KRAS突變會(huì)產(chǎn)生全轉(zhuǎn)錄組范圍的影響，這種組學(xué)影響稱為“組學(xué)漣漪”。正是由于“組學(xué)漣漪”的出現(xiàn)，才使得KRAS突變成為L(zhǎng)A的發(fā)生和發(fā)展中的“核心基因”。這種“組學(xué)漣漪”一方面通過促進(jìn)細(xì)胞周期轉(zhuǎn)換，抑制凋亡而增強(qiáng)腫瘤細(xì)胞的生存能力，另一方面，抑制免疫系統(tǒng)的抗腫瘤作用，使腫瘤微環(huán)境有利于腫瘤細(xì)胞的生存。正是由于這一系列的組學(xué)改變，最終使得KRAS突變型LA出現(xiàn)耐藥和預(yù)后不良。

本研究表明，從全轉(zhuǎn)錄組范圍內(nèi)進(jìn)行研究，可以更加深入了解KRAS突變型LA的深層次病理機(jī)制，明晰KRAS突變?cè)贚A中的表觀調(diào)控機(jī)制，為KRAS突變型LA的臨床治療和醫(yī)學(xué)研究提供新的思路。