俞盛健，麻懷露，陳王洋，丁曉飛，陳光，梁勇

（1.溫州醫(yī)科大學 第一臨床醫(yī)學院，浙江 溫州 325035；2.臺州學院 醫(yī)學院，浙江 臺州 318000；3.河北北方學院 研究生院，河北 張家口 075000）

膠質(zhì)瘤是頭部最常見的原發(fā)性惡性腫瘤[1]。目前，主要的腦膠質(zhì)瘤治療方式有手術(shù)、化療、放療、靶向治療等，但是膠質(zhì)瘤患者的預后仍較差，其中膠質(zhì)母細胞瘤平均生存時間不超過18個月[2]。因此，探討腦膠質(zhì)瘤發(fā)生發(fā)展的潛在分子機制來發(fā)現(xiàn)有效的診斷和治療靶點非常重要。長鏈非編碼RNA（long non-coding RNA，lncRNA）是一類長度超過200個堿基的非編碼RNA。雖然它不參與蛋白質(zhì)的編碼，但是近來許多研究發(fā)現(xiàn)lncRNA參與了體內(nèi)多種生命活動，如胚胎干細胞的多樣性、細胞周期的調(diào)節(jié)，以及癌癥的發(fā)生發(fā)展等[3]。癌組織中l(wèi)ncRNA的差異表達可以通過調(diào)控DNA的突變，如改變增強子和啟動子的活性或者染色質(zhì)狀態(tài)廣泛影響轉(zhuǎn)錄來實現(xiàn)[4]。同樣，lncRNA的差異表達也可以用來預測因組織DNA突變而導致的癌癥發(fā)生。如今越來越多的研究證實lncRNA可以作為癌癥的診斷標志物，甚至作為一種治療選擇[5]。例如在早期手術(shù)切除的乳腺癌中，lncRNA HOTAIR的過表達對于腫瘤轉(zhuǎn)移和整體生存率具有很高的預測作用[6]。lncRNA PCA3作為前列腺癌患者尿中的一類生物標志物，與血清前列腺特異性抗原相比具有更高的特異性和預測價值，現(xiàn)已成為在前列腺癌診斷中一個非常有效可靠的指標[7]。

癌癥基因圖譜（the Cancer Genome Atlas，TCGA）是開放的公共數(shù)據(jù)平臺，其包括30種癌癥11 000例患者信息，對于這些巨大數(shù)據(jù)的挖掘分析，有助于研究者全面研究各種癌癥發(fā)生的分子機制，確定新的治療靶點和與腫瘤相關(guān)的分子標志物，從而提高對多種人類惡性疾病的診斷能力和治療效果[8]。基因型組織表達（GTEx）數(shù)據(jù)庫研究來自449名生前健康的人類捐獻者的7 000多份尸檢樣本，GTEx對幾乎所有轉(zhuǎn)錄基因的基因表達模式進行了觀察，從而能夠確定基因組中影響基因表達的特定區(qū)域。在本研究中，我們通過對聯(lián)合GTEx和TCGA數(shù)據(jù)庫中的腦膠質(zhì)瘤樣本進行基因轉(zhuǎn)錄水平的差異表達分析，構(gòu)建了可用于預測患者預后的風險模型，并提取了可能與膠質(zhì)瘤惡性進展有關(guān)的lncRNA靶點。

1 材料和方法

1.1 數(shù)據(jù)收集 從GDC Data Portal（https://portal.gdc.cancer.gov/）中下載人類腦膠質(zhì)瘤的RNA-Seq數(shù)據(jù)和臨床數(shù)據(jù)，時間截至2019-01-14。RNA-Seq數(shù)據(jù)包含lncRNA和mRNA，包括700個腦膠質(zhì)瘤樣本。按照TCGA官方指導說明中的方法，使用Data Transfer Toll進行數(shù)據(jù)下載。從UCSC Xena（https://xena.ucsc.edu/）中下載GTEx對應正常腦部組織表達數(shù)據(jù)，得到1 152個正常樣本。

1.2 差異基因的篩選 利用edgeR包對腦膠質(zhì)瘤以及相對應的正常樣本進行分析，篩選出差異表達的lncRNA和mRNA。篩選標準為P＜0.05，F(xiàn)old Change＞2。

1.3 生存分析以及多因素Cox模型的建立 根據(jù)lncRNA表達量的中位值將RNA表達量分為高表達組和低表達組，結(jié)合臨床信息，通過survival包對lncRNA做Kaplan-Meier生存分析并畫出生存曲線圖。使用survival包對lncRNA做單因素Cox分析來挑選生存相關(guān)基因（survival related genes，SRGs），并且根據(jù)P值排序篩選出P值最小的前20個lncRNA，對挑選出的20個lncRNA進行多因素Cox模型的建立，其中剔除有共表達關(guān)系的lncRNA。根據(jù)風險值的中位數(shù)將患者分為高風險組和低風險組并畫出風險生存曲線。最后通過ROC曲線評估多因素Cox模型的鑒別能力。

1.4 共表達基因預測 通過使用雙側(cè)Pearson相關(guān)系數(shù)檢驗lncRNA與蛋白編碼基因的表達關(guān)系（|皮爾遜相關(guān)系數(shù)|＞0.4，P＜0.001）。符合條件的蛋白編碼基因被認為是該lncRNA的共表達基因。

1.5 富集分析 根據(jù)前述所得的差異表達基因，通過DAVID在線軟件（https://david.ncif-crf.gov/）和京都基因與基因組百科全書KEGG通路數(shù)據(jù)庫（https://www.kegg.jp/）進行功能注釋和通路富集分析。

1.6 蛋白相互作用網(wǎng)絡構(gòu)建 根據(jù)預測共表達基因中重復出現(xiàn)基因進行蛋白相互作用網(wǎng)絡的構(gòu)建。利用在線工具STRING11（https://string-db.org/）構(gòu)建蛋白質(zhì)的相互作用網(wǎng)絡，統(tǒng)計互作網(wǎng)絡各基因鄰接節(jié)點數(shù)目。利用Cytoscape對獲得的蛋白互作網(wǎng)絡進行可視化，并通過MCODE插件對蛋白相互作用網(wǎng)絡中各節(jié)點進行關(guān)聯(lián)度分析（K值＞5）。

1.7 患者標本和臨床評估 患者標本和臨床評估來源于具有腦膠質(zhì)瘤組織學診斷，且在手術(shù)前未曾接受放療和化療的患者。從臺州學院附屬臺州市立醫(yī)院和臺州市中心醫(yī)院獲取5對腦膠質(zhì)瘤癌組織及癌旁組織。所有的樣本均在手術(shù)切除后立即冷凍在液氮中。研究通過臺州學院醫(yī)學院倫理委員會批準，并獲得患者的知情同意。

1.8 熒光實時定量PCR TRIzol試劑提取總RNA，Takara試劑盒完成反轉(zhuǎn)錄和qPCR檢測，HOXA-AS2引物序列為5'-CCCGTAGGAAGAACCGATGA-3'（正）和5'-TTTAGGCCTTCGCAGACAGC-3'（反），GAPDH引物序列為5'-ACCACAGTCCATGCCATC-3'（正）和5'-TCCACCCTG TTGCTGTA-3'（反）。PCR擴增循環(huán)：95 ℃ 10 min；40循環(huán)的95 ℃ 10 s，60 ℃ 60 s。通過ΔΔCT法半定量HOXA-AS2相對水平。

1.9 統(tǒng)計學處理方法 應用R3.3.5軟件進行相關(guān)圖形制作，差異顯著性的篩選閾值為log2foldchange（log2FC）不小于2，采用Kaplan-Meier法繪制生存曲線并應用log-rank方法比較生存率，多變量因素利用Cox回歸進行生存分析，雙側(cè)Pearson相關(guān)系數(shù)檢驗lncRNA與蛋白編碼基因的表達關(guān)系。采用雙側(cè)檢驗，P＜0.05為差異有統(tǒng)計學意義。

2 結(jié)果

2.1 膠質(zhì)瘤中差異表達的lncRNA和mRNA 通過比較腦膠質(zhì)瘤樣本和正常腦組織樣本，獲得差異表達的lncRNA 424個（211個lncRNA上調(diào)，213個lncRNA下調(diào)；見圖1A），mRNA 3 827個（1 618個mRNA上調(diào)，2 209個mRNA下調(diào)；見圖1B），差異有統(tǒng)計學意義（P＜0.05）。

圖1 腦膠質(zhì)瘤和正常腦組織中差異表達RNA的火山圖

2.2 生存分析以及多因素Cox模型的建立 運用單因素Cox分析探究差異表達的lncRNA與膠質(zhì)瘤患者生存期的關(guān)系，挑選出統(tǒng)計學最顯著的20個lncRNA（見表1）。利用這20個lncRNA，我們進行了多因素Cox分析?？紤]到臨床應用的經(jīng)濟性，我們剔除了有共表達關(guān)系的lncRNA，只保留其中的一個，最終得到包含9個lncRNA（見表2）的生存預后模型，并且進行生存曲線分析（見圖2A）。為了測試模型的準確性，我們計算了模型的AUC值并進行ROC曲線分析（見圖2B），AUC值達0.907，說明該模型能夠相對準確地預測患者的預后。

表1 經(jīng)單因素Cox模型挑選的lncRNA

2.3 膠質(zhì)瘤中顯著差異表達的lncRNA富集分析 對挑選出的9個lncRNA中的CRNDE、LINCOO994、LBX2-AS1、HOXA-AS2進行蛋白編碼基因共表達關(guān)系的預測，以相關(guān)系數(shù)大于0.6為臨界點，并對所得基因進行篩選，挑選出重復出現(xiàn)的基因進行富集分析。GO功能和KEGG通路富集分析結(jié)果顯示（見圖3），對其富集注釋結(jié)果提示（見表3-4），這些基因功能主要富集在通道蛋白活性調(diào)控，分子黏附等方面；介導的信號轉(zhuǎn)導通路主要集中于免疫缺陷病毒感染，神經(jīng)信號轉(zhuǎn)導相關(guān)通路。

表2 經(jīng)多因素Cox模型挑選的lncRNA

2.4 蛋白相互作用網(wǎng)絡的構(gòu)建和關(guān)鍵基因的篩選 在預測的共表達基因中，構(gòu)建重復基因的蛋白相互作用網(wǎng)絡，并與差異表達基因比對，發(fā)現(xiàn)33個上調(diào)基因，35個下調(diào)基因（見圖4A）。對蛋白相互作用網(wǎng)絡中各節(jié)點的鄰接數(shù)進行統(tǒng)計（見圖4B），其中KIF2C與CDC20的鄰接節(jié)點最多，分別是16個和13個。使用MCODE插件篩選出網(wǎng)絡中的關(guān)聯(lián)度最緊密的2個基因簇（見圖4C），其中一個基因簇包括KIF1A、KIF3A、KIF3C、SPTBN2、CAPZA1、KDELR2、KIF2C、KIF21B、TMED9、KDELR3、KDELR1。對其蛋白相互作用網(wǎng)絡關(guān)鍵基因簇功能進行注釋，該基因簇功能主要富集在逆行小泡介導轉(zhuǎn)運，微管轉(zhuǎn)運，ATP結(jié)合和細胞胞質(zhì)中。另一個基因簇包括GABBR1、ADCY5、GNG12、HTR1A、ANXA1、SSTR2、GNG5，功能主要富集在細胞胞膜上。

圖2 Cox分析篩選出的差異基因以及構(gòu)建的生存模型

圖3 預測共表達mRNA的GO（A）及KEGG（B）通路富集結(jié)果

2.5 HOXA-AS2在腦膠質(zhì)瘤組織及正常腦組織中的表達 為進一步驗證生物信息學分析結(jié)果，我們收集5例臨床腦膠質(zhì)瘤樣本，抽提RNA，進行Real time PCR定量檢測。結(jié)果如圖5所示，HOXA-AS2在腦膠質(zhì)瘤組織中表達水平明顯高于正常腦組織（見圖5）。

表3 GO富集注釋結(jié)果

表4 KEGG通路富集注釋結(jié)果

3 討論

膠質(zhì)瘤是最常見的腦腫瘤之一，在早期腦腫瘤中占50%～60%[9-10]。在WHO分類中，神經(jīng)膠質(zhì)瘤被歸類為I-IV級。I級和II級星型細胞瘤和少突神經(jīng)膠質(zhì)瘤是低級別膠質(zhì)瘤，III級和IV級星型細胞瘤和少突神經(jīng)膠質(zhì)瘤是高級別膠質(zhì)瘤[11]。目前有證據(jù)表明lncRNA作為一類新的非編碼基因調(diào)節(jié)因子，可能在膠質(zhì)瘤的廣泛的生物學過程中發(fā)揮重要作用，例如，通過充當癌基因或腫瘤抑制基因，它有助于膠質(zhì)瘤的發(fā)生、進展和其他惡性表型的發(fā)展[12]。在本研究中，我們使用TCGA膠質(zhì)瘤轉(zhuǎn)錄水平數(shù)據(jù)和臨床數(shù)據(jù)，分析篩選出了424個差異表達的lncRNA，再對這些差異表達lncRNA進行單因素Cox分析，挑選出20個代表性的lncRNA，其中有研究發(fā)現(xiàn)DGCR9在胃癌組織中高表達，并且有利于胃癌細胞系A(chǔ)GS，MGC803的增殖、侵襲和糖代謝[13]。

HOTAIRM1是HOXA基因簇的反義鏈的轉(zhuǎn)錄產(chǎn)物，最近有研究顯示其與急性髓性白血病、結(jié)直腸癌、多形性膠質(zhì)母細胞瘤有關(guān)[14-16]。然而，HOTAIRM1在胃癌中又通過靶向miR-17-5p抑制胃癌細胞系的增殖和遷移并誘導其凋亡[17]。我們還通過多因素Cox分析構(gòu)建了與腦膠質(zhì)瘤患者生存相關(guān)的生存模型。通過繪制出ROC曲線和計算AUC值來確定該模型的準確性。最終AUC值達0.907，說明該模型能夠相對準確地預測患者的預后，并且，我們對篩選出的CRNDE、LINCOO994、LBX2-AS1、HOXA-AS2這4個lncRNA進行了mRNA共表達的預測，利用這些mRNA的GO功能和KEGG通路富集結(jié)果來預測lncRNA可能參與的功能和通路。我們發(fā)現(xiàn)，這些基因主要在控制通道蛋白活性，分子黏附等功能富集。其中，有10個基因富集在控制鈣通道活性的功能里。鈣黏蛋白是鈣依賴性細胞-細胞黏附分子，在神經(jīng)系統(tǒng)發(fā)育成熟過程中起至關(guān)重要的作用，根據(jù)其結(jié)構(gòu)和功能，可分為12個家族，包括經(jīng)典鈣黏蛋白、橋粒鈣黏蛋白、T-鈣黏蛋白、7D-鈣黏蛋白、原鈣黏蛋白、CDH15和23、脂肪、Dachsous、Flamingo、Celsr、鈣調(diào)素和Ret[18]。在癌癥中，E-鈣黏蛋白作為侵襲抑制因子，常常被下調(diào)，而作為侵襲起始因子的N-鈣黏蛋白則被上調(diào)[19]。原鈣黏蛋白PCDH10甲基化抑制了腫瘤細胞的生長、遷移、侵襲和克隆形成[20]?；蜻€在免疫缺陷病毒感染，神經(jīng)信號轉(zhuǎn)導相關(guān)通路上富集。研究表明，HIV-1胞膜糖蛋白120可誘導膠質(zhì)瘤細胞分泌基質(zhì)金屬蛋白酶[21]。另外，已經(jīng)有報道稱CRNDE在非小細胞肺癌、肝癌、宮頸癌和腦膠質(zhì)瘤中起到促癌作用[22-25]。對中國膠質(zhì)瘤基因組圖譜數(shù)據(jù)庫的RNA測序數(shù)據(jù)進行分析，發(fā)現(xiàn)LBX2-AS1與膠質(zhì)瘤患者的預后和惡性進展有顯著關(guān)系[26]。HOXA-AS2不僅在急性髓系白血病細胞對阿霉素的耐藥性中起重要作用，也參與了肝癌的惡性進展過程[27-28]。

圖4 預測共表達基因中重復基因的蛋白相互作用網(wǎng)絡及關(guān)鍵基因簇分析

圖5 HOXA-AS2在腦膠質(zhì)瘤中的表達情況

綜上所述，本研究利用生物信息學方法分析TCGA數(shù)據(jù)，篩選潛在的膠質(zhì)瘤中表達差異顯著的基因，通過Cox分析構(gòu)建生存模型，并篩選出了可能為膠質(zhì)瘤分子致病機制和分子靶向治療提供參考的幾個lncRNA和mRNA，包括CRNDE、LINCOO994、LBX2-AS1和HOXA-AS2。并在有限例數(shù)臨床患者膠質(zhì)瘤樣本HOXA-AS2表達水平和病理進程相關(guān)性的分析中，進一步證實HOXA-AS2可能為診斷預警分子，但其潛在與膠質(zhì)瘤分子機制的關(guān)系仍需臨床及基礎(chǔ)實驗驗證。