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

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

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

1 資料與方法

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

從TCGA肺鱗癌數(shù)據(jù)庫中下載轉(zhuǎn)錄組數(shù)據(jù)、體細(xì)胞突變及臨床信息，得到501個(gè)肺鱗癌組織和49個(gè)非腫瘤組織的樣本數(shù)據(jù)，剔除6例患者生存時(shí)間不詳?shù)哪[瘤組織樣本。利用Perl編程語言對(duì)肺鱗癌組織的原始數(shù)據(jù)進(jìn)行整理和注釋，得到14 086個(gè)LncRNA的表達(dá)矩陣。從TCGA數(shù)據(jù)庫中提取已被廣泛報(bào)道的23個(gè)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)的表達(dá)矩陣，通過R語言Rx64 4.1.2 “l(fā)imma”軟件包對(duì)lncRNA和m6A調(diào)控因子的表達(dá)矩陣整合鑒定出1 136個(gè)lncRNA m6A修飾，在此基礎(chǔ)上提取肺鱗癌患者樣本的生存時(shí)間和生存狀態(tài)，最終得到27個(gè)具有肺鱗癌預(yù)后意義的lncRNA m6A修飾。

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

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

1.3 免疫療效評(píng)估

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

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

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

2 結(jié)果

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

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

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

2.2 預(yù)后風(fēng)險(xiǎn)模型的評(píng)估

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

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

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

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

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

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

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

2.4 免疫治療療效分析

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

圖6 免疫療效分析

3 討論

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

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

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

本研究鑒定出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個(gè)lncRNA m6A修飾，并以此構(gòu)建肺鱗癌的預(yù)后風(fēng)險(xiǎn)模型，K-M曲線顯示,高風(fēng)險(xiǎn)組非鱗癌患者OS明顯低于低風(fēng)險(xiǎn)組；AUC顯示,此預(yù)后風(fēng)險(xiǎn)模型用于預(yù)測(cè)肺鱗癌預(yù)后情況可信度高。隨后,通過分析TIDE、TME證實(shí)了此預(yù)后風(fēng)險(xiǎn)模型能夠有效預(yù)測(cè)ICI療效。

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