LncRNA m6A修飾在肺鱗癌的預(yù)后及免疫療效預(yù)測中的鑒定和價值分析*

顧 睿，莫春芬，林 蘋，張 潔，李 凱△

1.四川大學(xué)華西醫(yī)院 腫瘤中心實驗?zāi)[瘤研究室(成都 610041); 2.成都醫(yī)學(xué)院 基礎(chǔ)醫(yī)學(xué)院免疫學(xué)教研室(成都 610500)

據(jù)統(tǒng)計，2018年全球新增肺癌病例209萬人，且肺癌的病死率逐年升高[1]。肺鱗癌作為肺癌主要的組織學(xué)亞型，發(fā)病率占所有肺癌的25%～30%[2]。肺鱗癌臨床預(yù)后較差，常規(guī)放化療及靶向治療療效欠佳[3]。免疫檢查點抑制劑(immune checkpoint inhibitors，ICI)的應(yīng)用對于無免疫治療禁忌證的患者，無論是單一免疫療法還是聯(lián)合免疫療法，都已成為晚期肺鱗癌的一線治療標(biāo)準(zhǔn)[4]。因此，迫切需要可靠的生物標(biāo)志物用以評估肺鱗癌患者的預(yù)后及免疫療效，從而針對患者病情制定個性化的診療方案。LncRNA是一類轉(zhuǎn)錄本長度>200 nt的非編碼RNA,通過表觀遺傳調(diào)控、轉(zhuǎn)錄水平調(diào)控等多種方式參與編碼蛋白質(zhì)基因的功能調(diào)控[5]。M6A修飾是一種可逆性的RNA表觀遺傳調(diào)控方式，即腺苷酸(A)在甲基轉(zhuǎn)移酶的催化作用下，第6位N發(fā)生甲基化修飾，通過特定的閱讀蛋白識別發(fā)生m6A修飾的堿基，激活下游的調(diào)控通路[6]。LncRNA m6A修飾參與肺鱗癌發(fā)生發(fā)展的過程[7]。然而，lncRNA m6A修飾在肺鱗癌預(yù)后及免疫功能中的作用多為單個lncRNA的功能研究，不具有獨立的預(yù)測價值。因此，探究聯(lián)合lcRNA m6A修飾在預(yù)測肺鱗癌預(yù)后價值及免疫療效價值中的意義重大。

1 資料與方法

1.1 數(shù)據(jù)提取與整理

從TCGA肺鱗癌數(shù)據(jù)庫中下載轉(zhuǎn)錄組數(shù)據(jù)、體細(xì)胞突變及臨床信息，得到501個肺鱗癌組織和49個非腫瘤組織的樣本數(shù)據(jù)，剔除6例患者生存時間不詳?shù)哪[瘤組織樣本。利用Perl編程語言對肺鱗癌組織的原始數(shù)據(jù)進行整理和注釋，得到14 086個LncRNA的表達矩陣。從TCGA數(shù)據(jù)庫中提取已被廣泛報道的23個m6A調(diào)控因子(甲基轉(zhuǎn)移酶METTL3、METTL14、METTL16、WTAP、VIRMA、ZC3H13、RBM15、RBM15B；閱讀蛋白YTHDC1、YTHDC2、YTHDF1、YTHDF2、YTHDF3、HNRNPC、FMR1、LRPPRC、HNRNPA2B1、RBMX、IGFBP1、IGFBP2、IGFBP3；去甲基轉(zhuǎn)移酶FTO、ALKBH5)的表達矩陣，通過R語言Rx64 4.1.2 “l(fā)imma”軟件包對lncRNA和m6A調(diào)控因子的表達矩陣整合鑒定出1 136個lncRNA m6A修飾，在此基礎(chǔ)上提取肺鱗癌患者樣本的生存時間和生存狀態(tài)，最終得到27個具有肺鱗癌預(yù)后意義的lncRNA m6A修飾。

1.2 預(yù)后風(fēng)險模型的構(gòu)建及評估

基于上述27個lncRNA m6A修飾，利用R語言Rx64 4.1.2 “Survival”“caret”“glmnet”“survminer”軟件包構(gòu)建Lasso回歸分析的內(nèi)置函數(shù)cvfit和lambda，根據(jù)該函數(shù)中的“λ”值對構(gòu)建預(yù)后風(fēng)險模型的lncRNA m6A修飾個數(shù)進行取值，“λ”值代表該模型的過擬合程度，當(dāng)“λ”值為13時模型的擬合程度最佳，對應(yīng)地鑒定出13個lncRNA m6A修飾以構(gòu)建預(yù)后風(fēng)險模型，即LINC02332、AL591686.1、DHRS4-AS1、AC008734.1、AC138035.1、AL390195.3、AP001189.3、WT1-AS、AL391095.4、AC130651.1、AC015922.2、HORMAD2-AS1、AC243919.2。隨后，利用風(fēng)險評分公式對上述13個lncRNA m6A修飾的相關(guān)回歸系數(shù)進行加權(quán)計算。為驗證該預(yù)后風(fēng)險模型對肺鱗癌預(yù)后的預(yù)測能力，首先根據(jù)5∶5留出法對該預(yù)后風(fēng)險模型中的肺鱗癌樣本進行隨機抽樣得到訓(xùn)練集和測試集，訓(xùn)練集和測試集的患者樣本在年齡、性別、TMN分期等變量上基本無差異。之后將訓(xùn)練集的患者樣本根據(jù)風(fēng)險評分的中位值劃分為高、低風(fēng)險組，并以訓(xùn)練集的風(fēng)險評分中位值劃分測試集的高、低風(fēng)險組。本研究利用“tidyverse” “ggplot2” “ggExtra”軟件包分析13個LncRNA m6A修飾與納入的23個m6A調(diào)控因子之間的相關(guān)性；“survival”“timeROC”“survminer”軟件包比較高、低風(fēng)險組間的總生存期(OS)、ROC曲線下面積(AUC)差異及獨立風(fēng)險因素；“org.Hs.eg.db”“clusterProfiler”“enrichplot”等軟件包對高、低風(fēng)險組間差異基因進行GO富集分析；“GSVA”“GSEABase”“reshape2”“pheatmap”軟件包進行高低風(fēng)險組差異基因免疫功能分析。

1.3 免疫療效評估

從TCGA數(shù)據(jù)庫中下載肺鱗癌的體細(xì)胞突變數(shù)據(jù)，并應(yīng)用perl編程語言和“maftools”軟件包進行數(shù)據(jù)提取和可視化整理，得到高、低風(fēng)險組每個樣本的腫瘤突變負(fù)荷(tumor mutation burden,TMB)，即每100萬個堿基中發(fā)生突變的堿基數(shù)，通過“l(fā)imma” “ggplot2”軟件包分析高、低風(fēng)險組間的TMB差異。從腫瘤免疫功能障礙和排斥(tumor immune dysfunction and exclusion,TIDE)數(shù)據(jù)庫中下載NSCLC患者的TIDE評分，利用“ggpubr” “l(fā)imma”軟件包分析高、低風(fēng)險組間的TIDE差異。通過CIBERSORT算法計算高、低風(fēng)險組樣本的免疫細(xì)胞評分、基質(zhì)細(xì)胞評分及TME總評分，分析高、低風(fēng)險組間的TME差異。

1.4 統(tǒng)計學(xué)方法

采用Perl編程語言和R語言Rx64 4.1.2進行數(shù)據(jù)處理，采用SPASS 22.0軟件進行數(shù)據(jù)統(tǒng)計分析。單因素Cox回歸分析取95%置信區(qū)間(CI)。采用雙側(cè)Log-rank檢驗分析Kaplan-Meier(K-M)曲線，Wilcoxon檢驗分析兩組變量之間的差異。檢驗水準(zhǔn)α除Pearson分析設(shè)定為0.001外，其余統(tǒng)計學(xué)分析均設(shè)定為0.05。

2 結(jié)果

2.1 構(gòu)建預(yù)后風(fēng)險模型

基于TCGA肺鱗癌數(shù)據(jù)庫，利用相關(guān)性分析得到1 136個lncRNA m6A修飾(圖1A)。利用單因素Cox回歸分析得到具有預(yù)后意義的27個lncRNA m6A修飾(圖1B)，其中AL591686.1、AC008734.1、AP001189.3、WT1-AS、AL391095.4、AC015922.2、HORMAD2-AS1、DSCR9與肺鱗癌預(yù)后呈負(fù)相關(guān)，其他19個lncRNA m6A修飾呈正相關(guān)。利用Lasso回歸分析cvfit和lambda函數(shù)，鑒定出13個聯(lián)合lncRNA m6A修飾預(yù)后標(biāo)志物用以構(gòu)建預(yù)后風(fēng)險模型，即LINC02332、AL591686.1、DHRS4-AS1、AC008734.1、AC138035.1、AL390195.3、AP001189.3、WT1-AS、AL391095.4、AC130651.1、AC015922.2、HORMAD2-AS1、AC243919.2(圖1C～D)。上述13個LncRNA m6A修飾的表達與不同的m6A調(diào)控因子的表達相關(guān)，其中AC138035.1與HNRNPA2B1、HORMAD2-AS1與IGFBP1、AC243919.2與YTHDC1的表達之間存在較強的正相關(guān)，AP001189.3與HNRNPC的表達之間存在較強的負(fù)相關(guān)(圖1E)。

圖1 預(yù)后風(fēng)險模型的構(gòu)建

2.2 預(yù)后風(fēng)險模型的評估

訓(xùn)練集和測試集基本互斥且不存在分組差異(表1)。隨著樣本風(fēng)險評分增加，患者的生存率下降(圖2A～B)。K-M曲線分析表明,在訓(xùn)練集和測試集中均顯示高風(fēng)險組患者的OS低于低風(fēng)險組(圖2C～D)。同時，單因素Cox回歸分析顯示,相較于年齡、性別、臨床分期等臨床特征，風(fēng)險評分是肺鱗癌相關(guān)的OS預(yù)后指標(biāo)(圖3A～B);多因素Cox回歸分析表明,風(fēng)險評分作為獨立的預(yù)后因子，預(yù)測價值高(P<0.001)(圖3C～D)。

表1 訓(xùn)練集與測試集TCGA數(shù)據(jù)庫肺鱗癌患者樣本的臨床信息及分組情況[n(%)]

圖2 預(yù)后風(fēng)險模型的風(fēng)險趨勢與K-M曲線

圖3 單因素、多因素Cox回歸分析預(yù)后風(fēng)險因素

2.3 GO通路富集分析及免疫功能差異分析

GO通路富集分析結(jié)果顯示，上述67個差異基因主要執(zhí)行信號受體激活因子活性、受體配體活性等功能，并參與ERK1/2級聯(lián)調(diào)控等生物學(xué)過程(圖4)，高、低風(fēng)險組的差異基因與趨化因子受體(CCR)、人類白細(xì)胞抗原(HLA)、干擾素(IFN)Ⅰ/Ⅱ型反應(yīng)、MHC Ⅰ類分子、細(xì)胞活性、炎癥刺激、T細(xì)胞共抑制/共刺激分子、免疫檢查點等13種免疫分子及功能密切相關(guān)(圖5)。

圖4 高、低風(fēng)險組差異基因的GO通路富集分析 圖5 高、低風(fēng)險組免疫相關(guān)功能差異分析

2.4 免疫治療療效分析

目前對于免疫治療療效的分析主要包括TMB、TIDE及TME評分等方法。TMB能夠衡量體細(xì)胞的非同義突變，即產(chǎn)生新抗原或肽段被自身免疫系統(tǒng)識別為非自身抗原的概率[8]。TMB瀑布圖顯示,TP53在高、低風(fēng)險組中突變頻率最高(圖6A～B)，高、低風(fēng)險組的TMB并無明顯差異(圖6C)。TIDE評分表明,高風(fēng)險組的TIDE評分高于低風(fēng)險組(圖6D)。此外，TME評分表明,高風(fēng)險組的免疫細(xì)胞評分(圖6E)、基質(zhì)細(xì)胞評分(圖6F)以及TME總評分(圖6G)均高于低風(fēng)險組。

圖6 免疫療效分析

3 討論

相較于其他類型的非小細(xì)胞肺癌，肺鱗癌體細(xì)胞突變頻率更高、免疫原性更強，患者更適合應(yīng)用ICI免疫療法改善預(yù)后[9]，因此探索能夠評估肺鱗癌患者的預(yù)后及ICI療效的生物標(biāo)志物不可或缺。研究[10]發(fā)現(xiàn),AL049840.3、AC008770.3、AL355312.3等11個lncRNA m6A修飾與胃癌的不良預(yù)后密切相關(guān)，通過聚類分析發(fā)現(xiàn),C1和C2亞組的胃癌患者免疫治療的應(yīng)用前景更好。在膀胱癌患者樣本的分析中，研究[11]發(fā)現(xiàn),PTOV1-AS2、EHMT2-AS1、KCNQ1OT1等lncRNA m6A修飾可能是參與膀胱癌免疫微環(huán)境的重要組成。本研究涉及的LncRNA m6A修飾中，DHRS4-AS1被報道通過拮抗TP53和TET1調(diào)節(jié)腫瘤干細(xì)胞在NSCLC中的干性，從而發(fā)揮抑癌作用[12]。研究[13]發(fā)現(xiàn)，AP001189.3與鐵蛋白代謝相關(guān)，并能夠影響結(jié)腸腺癌患者的結(jié)局。WT1-AS廣泛表達于實體腫瘤和血液腫瘤中，在膠質(zhì)瘤中通過調(diào)節(jié)miR-494-3p抑制惡化[14]，在非霍奇金淋巴瘤的患者中檢測到高表達WT1的特異性T細(xì)胞[15]，同時研究[16]發(fā)現(xiàn)，WT1能夠作為胸腺上皮腫瘤的免疫治療靶點，且能作為免疫抑制細(xì)胞因子IL-10(IL-10)新的轉(zhuǎn)錄激活因子，調(diào)控IL-10基因表達，從而發(fā)揮免疫調(diào)節(jié)作用[17]。因此，本研究基于上述lncRNA m6A修飾所建立的預(yù)后風(fēng)險模型可信度高，以預(yù)后模型的回歸系數(shù)為基礎(chǔ)得到風(fēng)險評分不僅能夠作為獨立的預(yù)后因子，同時能夠較準(zhǔn)確地預(yù)測患者的1、3、5年生存情況。上述結(jié)果表明，本研究所構(gòu)建的預(yù)后風(fēng)險模型在肺鱗癌預(yù)后的預(yù)測方面是高度可信且穩(wěn)定的。

ICI免疫治療作為肺鱗癌患者改善預(yù)后的關(guān)鍵突破口，其應(yīng)答效率可通過TMB、TIDE和TME評分進行評估[18]。在一項晚期NSCLC的研究[19]中發(fā)現(xiàn),高TMB患者ICI的有效作用時間為64周，比低TMB患者更長。而另一項研究[20]中，在阿特珠單抗(PD-L1抑制劑)單藥治療時，評估晚期NSCLC患者治療前的TMB發(fā)現(xiàn)，治療后高TMB患者的無進展生存期(PFS)得到明顯改善。本研究中，低風(fēng)險組和高風(fēng)險組整體突變頻率無明顯差異。有研究[21]通過對TGFB1和SOX10等影響TIDE的基因進行綜合分析，在轉(zhuǎn)移性黑色素瘤中發(fā)現(xiàn)，TIDE評分較高的患者OS更高。在肝癌的療效評估中，研究[22]發(fā)現(xiàn)，TIDE評分較高的高危組(113/186)比低危組(89/186)免疫治療應(yīng)答患者人數(shù)更多。類似地，本研究中高風(fēng)險組的TIDE得分明顯高于低風(fēng)險組，初步說明ICI在高風(fēng)險組中療效更佳。此外，本研究發(fā)現(xiàn)，高風(fēng)險組的患者TME評分更高，即高風(fēng)險組患者的ICI應(yīng)用更有效。上述結(jié)果表明，本研究所構(gòu)建的預(yù)后風(fēng)險模型預(yù)測肺鱗癌的免疫療效可靠性強。

然而，本研究仍存在一些不足。首先，用于構(gòu)建肺鱗癌患者預(yù)后風(fēng)險模型的樣本量相對不足，因其僅包含了TCGA數(shù)據(jù)庫。其次，研究結(jié)果僅依賴于生物信息學(xué)分析，后期仍需通過大量基礎(chǔ)實驗和臨床試驗驗證13個lncRNA m6A修飾在肺鱗癌中的作用，同時對接受免疫治療和化療肺鱗癌患者進行多中心大規(guī)模的驗證性臨床試驗。

本研究鑒定出LINC02332、AL591686.1、DHRS4-AS1、AC008734.1、AC138035.1、AL390195.3、AP001189.3、WT1-AS、AL391095.4、AC130651.1、AC015922.2、HORMAD2-AS1、AC243919.2共13個lncRNA m6A修飾，并以此構(gòu)建肺鱗癌的預(yù)后風(fēng)險模型，K-M曲線顯示,高風(fēng)險組非鱗癌患者OS明顯低于低風(fēng)險組；AUC顯示,此預(yù)后風(fēng)險模型用于預(yù)測肺鱗癌預(yù)后情況可信度高。隨后,通過分析TIDE、TME證實了此預(yù)后風(fēng)險模型能夠有效預(yù)測ICI療效。

綜上所述，本研究鑒定的13個lncRNA m6A修飾聯(lián)合應(yīng)用對評估肺鱗癌預(yù)后及ICI療效具有重要價值，為肺鱗癌ICI療法的應(yīng)用及敏感性化療藥物的開發(fā)提供了重要的理論依據(jù)。