基于N6-甲基腺苷修飾與免疫細(xì)胞浸潤(rùn)的結(jié)直腸癌預(yù)后模型構(gòu)建及驗(yàn)證

楊菁菁, 郭懷娟, 茅靜賢, 王佳欣, 王 穎, 嚴(yán)雪冰, 潘曉萍

(揚(yáng)州大學(xué)附屬醫(yī)院, 1. 腫瘤科, 2. 健康管理中心, 江蘇 揚(yáng)州, 225000)

結(jié)直腸癌(CRC)是人類常見(jiàn)的消化道惡性腫瘤,在所有腫瘤中發(fā)病率排名第3位,死亡率排名第2位[1]。盡管近年來(lái)風(fēng)險(xiǎn)篩查、外科技術(shù)和靶向治療方面取得了巨大進(jìn)步,但仍有相當(dāng)一部分CRC患者預(yù)后不佳。CRC的發(fā)生發(fā)展是一個(gè)多步驟、多基因的復(fù)雜生物學(xué)過(guò)程,而如何選擇合適的分子標(biāo)志物幫助判斷CRC預(yù)后尤為重要。近年來(lái)研究[2]表明, N6-甲基腺苷(m6A)修飾不僅參與CRC發(fā)生發(fā)展,而且在調(diào)控抗癌藥物敏感性過(guò)程中也發(fā)揮著重要作用。例如,甲基轉(zhuǎn)移酶樣蛋白14(METTL14)作為催化m6A修飾主要甲基化酶,不僅參與調(diào)控CRC細(xì)胞的生長(zhǎng)和轉(zhuǎn)移,而且其表達(dá)水平可能影響免疫及靶向治療的敏感性[3]。在臨床實(shí)踐中, CRC組織特異性免疫細(xì)胞浸潤(rùn)也是一個(gè)至關(guān)重要的預(yù)后因素[4]。例如, CD3+T細(xì)胞和細(xì)胞毒性CD8+T細(xì)胞浸潤(rùn)與CRC患者的術(shù)后復(fù)發(fā)及總體生存密切相關(guān)[5]。盡管m6A修飾及免疫細(xì)胞浸潤(rùn)均與CRC預(yù)后相關(guān),但是兩者相關(guān)分子能否共同用于預(yù)后判斷尚不清楚。本研究擬利用生物信息學(xué)方法篩選m6A修飾及免疫細(xì)胞浸潤(rùn)相關(guān)分子,構(gòu)建預(yù)后模型并進(jìn)行驗(yàn)證,為CRC患者的精準(zhǔn)預(yù)后評(píng)估提供新的策略,現(xiàn)報(bào)告如下。

1 資料與方法

1.1 CRC患者相關(guān)數(shù)據(jù)來(lái)源

從癌癥基因組圖譜(TCGA, https://portal.gdc.cancer.gov/)下載612例結(jié)直腸癌組織的mRNA表達(dá)譜分析結(jié)果和臨床數(shù)據(jù)。從基因表達(dá)綜合數(shù)據(jù)庫(kù)(GEO, https://www.ncbi.nlm.nih.gov/geo/)中下載歸一化基因表達(dá)數(shù)據(jù)和臨床數(shù)據(jù),數(shù)據(jù)集登錄號(hào)為GSE14333(n=226)、GSE17536(n=177)、GSE39582(n=585)和GSE38832(n=122)。

1.2 研究流程

首先,利用m6A調(diào)控因子的mRNA表達(dá)譜,基于共識(shí)聚類方法對(duì)TCGA數(shù)據(jù)庫(kù)中612例CRC患者的預(yù)后進(jìn)行分層,在同一隊(duì)列中,基于單樣本基因組富集分析(ssGSEA)檢測(cè)24種免疫細(xì)胞的豐度,并對(duì)患者預(yù)后進(jìn)行分層; 其次,分析隊(duì)列中m6A調(diào)節(jié)因子與免疫細(xì)胞的相關(guān)性,選擇612例CRC組織和51例正常組織之間的差異表達(dá)基因,基于加權(quán)基因共表達(dá)網(wǎng)絡(luò)分析(WGCNA)篩選與m6A修飾及免疫細(xì)胞浸潤(rùn)相關(guān)的基因; 最后,基于Lasso回歸分析建立預(yù)后模型,一方面利用TCGA數(shù)據(jù)庫(kù)分析預(yù)后模型在不同亞組中的預(yù)測(cè)功效,另一方面選取4個(gè)GEO數(shù)據(jù)庫(kù)隊(duì)列作為外部驗(yàn)證隊(duì)列,進(jìn)一步明確預(yù)后模型的臨床價(jià)值。見(jiàn)圖1。

1.3 CRC組織免疫細(xì)胞浸潤(rùn)分析

利用ssGSEA計(jì)算24種免疫細(xì)胞的豐度并進(jìn)一步聚類,解析樣本中不同類別。根據(jù)上一步得到的樣本類別,利用生存分析解析不同類別在預(yù)后方面是否具有顯著差異。根據(jù)不同免疫細(xì)胞豐度,結(jié)合樣本的生存情況,利用單因素Cox回歸分析解析免疫細(xì)胞類型對(duì)CRC預(yù)后的影響。

1.4 WGCNA篩選m6A修飾及免疫細(xì)胞浸潤(rùn)相關(guān)預(yù)后基因

1.5 統(tǒng)計(jì)學(xué)分析

采用t檢驗(yàn)對(duì)各組間m6A調(diào)控因子表達(dá)水平進(jìn)行比較; 采用Kaplan-Meier模型構(gòu)建生存曲線,并使用Log-rank檢驗(yàn)比較組間生存差異; 采用單變量Cox回歸分析評(píng)估免疫細(xì)胞比例與患者預(yù)后的相關(guān)性; 采用Lasso回歸分析篩選出與預(yù)后相關(guān)的基因; 根據(jù)基因表達(dá)水平和Lasso回歸系數(shù)建立預(yù)后模型; 采用受試者工作特征(ROC)曲線對(duì)所建預(yù)后模型的性能進(jìn)行評(píng)估;P<0.05為差異有統(tǒng)計(jì)學(xué)意義。

2 結(jié) 果

2.1 m6A分子表型及免疫細(xì)胞浸潤(rùn)在CRC患者預(yù)后中的意義

首先,對(duì)612個(gè)CRC患者組織樣本進(jìn)行聚類并確定最佳聚類K值為4(圖2A)。基于總生存期(OS)將612例CRC患者分為亞組1和亞組2(P=0.015, 圖2B), 發(fā)現(xiàn)14個(gè)m6A調(diào)節(jié)因子[如RNA結(jié)合序列蛋白(RBM15)、RBM15B和KIAA1429]在2個(gè)亞組間的表達(dá)差異顯著(圖2C)。同時(shí),基于ssGSEA方法計(jì)算上述患者腫瘤組織中24個(gè)免疫細(xì)胞的豐度,進(jìn)而依據(jù)免疫微環(huán)境(TME)將患者劃分為TME1組(n=180)和TME2組(n=432)(圖3A); 生存分析表明, TME1組患者的OS率低于TME2組,差異有統(tǒng)計(jì)學(xué)意義(P=0.049, 圖3B); 單變量Cox回歸分析表明,免疫細(xì)胞浸潤(rùn)比例與OS無(wú)相關(guān)性(圖3C); 進(jìn)一步研究發(fā)現(xiàn)10種m6A調(diào)節(jié)因子(如RBM15、RBM15B和WTAP)在TME1組與TME2組間的表達(dá)存在顯著差異(圖3D); 相關(guān)性分析表明,大部分m6A調(diào)節(jié)因子與部分免疫細(xì)胞的浸潤(rùn)存在顯著相關(guān)性(圖3E), 例如CRC組織中Alk B同源蛋白5(ALKBH5)表達(dá)水平與所有免疫細(xì)胞浸潤(rùn)比例均呈正相關(guān)。

2.2 篩選m6A修飾和免疫細(xì)胞浸潤(rùn)相關(guān)預(yù)后基因

從TCGA數(shù)據(jù)庫(kù)中篩選出612個(gè)CRC樣本和51個(gè)正常組織樣本的差異表達(dá)基因,納入標(biāo)準(zhǔn)為: ① |logFC|>1; ②P<0.01。與CRC組織樣本相比,對(duì)照樣本中分別有1 209個(gè)基因和922個(gè)基因顯著上調(diào)和下調(diào)(圖4A)。基于2 131個(gè)基因構(gòu)建了加權(quán)基因共表達(dá)網(wǎng)絡(luò),最佳β值定為6。對(duì)基因進(jìn)行聚類,每個(gè)模塊的基因數(shù)不少于30個(gè)。經(jīng)過(guò)動(dòng)態(tài)樹(shù)切割后,對(duì)緊密相鄰的模塊進(jìn)行聚類,最終得到9個(gè)模塊。基于本研究發(fā)現(xiàn)m6A分子表型與OS的相關(guān)性高于免疫細(xì)胞浸潤(rùn),因此選擇與m6A調(diào)節(jié)因子相關(guān)的模塊進(jìn)行下一步研究。根據(jù)納入標(biāo)準(zhǔn)(MM>0.4; GS>0.2)篩選出25個(gè)核心基因(圖4B)。同時(shí),相關(guān)性研究發(fā)現(xiàn)部分核心基因表達(dá)水平(如ATOH1)與免疫細(xì)胞浸潤(rùn)相關(guān)(圖4C、4D)。

2.3 基于m6A修飾及免疫細(xì)胞浸潤(rùn)相關(guān)基因構(gòu)建預(yù)后模型

基于Lasso回歸分析從25個(gè)核心基因中篩選出4個(gè)預(yù)后基因: 內(nèi)凝集蛋白1(ITLN1)、淋巴細(xì)胞抗原6復(fù)合位點(diǎn)G6D(LY6G6D)、無(wú)調(diào)性同族體1(ATOH1)和基質(zhì)金屬蛋白酶(MMP28)。基于4個(gè)預(yù)后基因在GEPIA數(shù)據(jù)庫(kù)中的表達(dá)差異分析,基因LY6G6D與MMP28在結(jié)腸癌與直腸癌組織中的表達(dá)較癌旁組織有差異性(圖5A、5B)。基于4個(gè)基因的Lasso分析構(gòu)建了預(yù)后模型: 風(fēng)險(xiǎn)評(píng)分=EXPITLN1×(-0.050 23)+EXPLY6G6D×(-0.022 89)+EXPATOH1×(-0.141 97)+EXPMMP28×(-0.045 74)。使用公式計(jì)算訓(xùn)練隊(duì)列每個(gè)樣本的風(fēng)險(xiǎn)評(píng)分,并根據(jù)中位值將樣本分為高風(fēng)險(xiǎn)組和低風(fēng)險(xiǎn)組。每個(gè)樣本風(fēng)險(xiǎn)評(píng)分分布及對(duì)應(yīng)患者生存信息見(jiàn)圖5C、5D。生存分析表明,高風(fēng)險(xiǎn)組的OS低于低風(fēng)險(xiǎn)組(圖5E)。ROC曲線分析表明,預(yù)后模型曲線下面積(AUC)為0.634, 優(yōu)于年齡和TNM分期等其他臨床因素(圖5F)。

2.4 預(yù)后模型在訓(xùn)練隊(duì)列不同亞組中的預(yù)測(cè)效能

為了進(jìn)一步探討預(yù)后模型的臨床價(jià)值,將TCGA數(shù)據(jù)庫(kù)隊(duì)列根據(jù)年齡、性別及腫瘤位置分為不同亞組進(jìn)行分析。結(jié)果表明,在>60歲患者中,高風(fēng)險(xiǎn)組的OS顯著低于低風(fēng)險(xiǎn)組的OS(P=0.003); 在≤60歲患者中,高風(fēng)險(xiǎn)組與低風(fēng)險(xiǎn)組OS無(wú)顯著差異(P=0.100)。在女性患者及男性患者中,高風(fēng)險(xiǎn)組OS率均顯著低于低風(fēng)險(xiǎn)組(P女性=0.014,P男性=0.047)。結(jié)腸癌高風(fēng)險(xiǎn)組OS率顯著低于低風(fēng)險(xiǎn)組,但直腸癌高風(fēng)險(xiǎn)組與低風(fēng)險(xiǎn)組OS率無(wú)顯著差異(P結(jié)腸癌=0.003,P直腸癌=0.210)。見(jiàn)圖6。

2.5 預(yù)后模型在驗(yàn)證隊(duì)列中的預(yù)測(cè)效能

從GEO數(shù)據(jù)庫(kù)中下載了4個(gè)CRC隊(duì)列作為驗(yàn)證隊(duì)列以進(jìn)一步明確預(yù)后模型的臨床價(jià)值。在GSE17536、GSE39582和GSE38832隊(duì)列中,高風(fēng)險(xiǎn)組OS顯著低于低風(fēng)險(xiǎn)組(PGSE17536=0.013,PGSE39582<0.001,PGSE38832=0.045)。在GSE14333隊(duì)列中,高風(fēng)險(xiǎn)組的無(wú)病生存率顯著低于低風(fēng)險(xiǎn)組(P=0.016)。見(jiàn)圖7。

3 討 論

研究[6-7]表明,分子標(biāo)志物能夠精準(zhǔn)指導(dǎo)CRC患者的早期診斷、預(yù)后評(píng)估及治療藥物選擇。m6A修飾是CRC發(fā)生發(fā)展中的重要分子事件之一,其調(diào)控因子(包括書(shū)寫(xiě)蛋白、閱讀蛋白及識(shí)別蛋白)在CRC組織中表達(dá)失調(diào),且可作為患者診斷及預(yù)后的潛在標(biāo)志物[8]。本課題組團(tuán)隊(duì)前期基于m6A修飾及上皮細(xì)胞間質(zhì)轉(zhuǎn)化(EMT)篩選預(yù)后分子并構(gòu)建相關(guān)臨床模型,該模型可在多個(gè)CRC隊(duì)列中識(shí)別預(yù)后較差的高風(fēng)險(xiǎn)人群[9]。本研究基于m6A修飾及免疫細(xì)胞浸潤(rùn)篩選與CRC預(yù)后相關(guān)的分子標(biāo)志物(ITLN1、LY6G6D、ATOH1和MMP28), 并構(gòu)建相應(yīng)的預(yù)后模型,該模型不僅在不同臨床特征的亞組中能識(shí)別預(yù)后不良的高風(fēng)險(xiǎn)人群,而且在多個(gè)外部驗(yàn)證隊(duì)列中表現(xiàn)出較好的臨床結(jié)局預(yù)測(cè)能力。

本研究首先從TCGA數(shù)據(jù)庫(kù)中提取612例CRC患者臨床資料及匹配的腫瘤組織轉(zhuǎn)錄組信息,利用聚類分析將患者分為不同預(yù)后差異的2個(gè)亞組。進(jìn)一步分析表明,包括RBM15、ALKBH5和脂肪量與肥胖相關(guān)蛋白(FTO)在內(nèi)的14個(gè)調(diào)控因子在高風(fēng)險(xiǎn)組中較低風(fēng)險(xiǎn)組顯著上調(diào)或下調(diào),提示上述調(diào)控因子可能參與CRC發(fā)生發(fā)展。在機(jī)制研究中,敲除RBM15可影響下游髓樣分化因子88的m6A甲基化修飾水平,進(jìn)而抑制CRC細(xì)胞的增殖及轉(zhuǎn)移能力[10]。低表達(dá)ALKBH5與CRC患者不良預(yù)后相關(guān),其可能通過(guò)m6A修飾途徑抑制PHD鋅指蛋白20的mRNA穩(wěn)定性而發(fā)揮抑癌效應(yīng)[11]。FTO通過(guò)介導(dǎo)鋅指蛋白687的m6A修飾進(jìn)而激活Wnt信號(hào)通路,最終促進(jìn)CRC的增殖、轉(zhuǎn)移及血管生成[12]。在臨床研究中, m6A修飾識(shí)別蛋白如YTH結(jié)構(gòu)域家族蛋白2和3高表達(dá)水平提示CRC患者不良預(yù)后,而胰島素樣生長(zhǎng)因子2 mRNA結(jié)合蛋白1高表達(dá)水平則提示相反結(jié)局[13]。值得注意的是, m6A修飾書(shū)寫(xiě)蛋白甲基轉(zhuǎn)移酶樣蛋白3(METTL3)及METTL14的表達(dá)水平在高風(fēng)險(xiǎn)組與低風(fēng)險(xiǎn)組間無(wú)顯著差異,可能與既往部分研究[14-15]結(jié)果相悖,因而需更多的回顧性研究來(lái)確證。

本研究同時(shí)發(fā)現(xiàn)多個(gè)免疫細(xì)胞浸潤(rùn)比例能夠共同顯著區(qū)分TCGA數(shù)據(jù)庫(kù)中CRC患者的OS, 提示其在CRC預(yù)后評(píng)估中存在潛在應(yīng)用價(jià)值。在既往文獻(xiàn)中, CD3+和CD8+T細(xì)胞浸潤(rùn)比例與CRC患者的良好預(yù)后呈正相關(guān),其機(jī)制可能與Toll樣受體介導(dǎo)的抗腫瘤免疫系統(tǒng)激活相關(guān)[16]。在微衛(wèi)星穩(wěn)定型伴肺轉(zhuǎn)移的CRC患者中,基于上述2種T細(xì)胞的免疫評(píng)分可以預(yù)測(cè)患者術(shù)后5年的OS狀態(tài)[17]。自然殺傷細(xì)胞浸潤(rùn)雖然與CRC患者OS無(wú)顯著相關(guān)性,但是與無(wú)疾病間隔時(shí)間及無(wú)復(fù)發(fā)生存時(shí)間呈正相關(guān)[18]。值得注意的是,單因素Cox回歸分析表明24種免疫細(xì)胞與預(yù)后均無(wú)顯著相關(guān)性,這與既往部分研究結(jié)果相似。例如,一項(xiàng)生物信息學(xué)研究[19]表明,在GEO數(shù)據(jù)庫(kù)CRC隊(duì)列中大多數(shù)免疫細(xì)胞浸潤(rùn)比例與CRC患者的OS無(wú)顯著相關(guān)性。因此,一方面,免疫細(xì)胞浸潤(rùn)在CRC中的預(yù)后價(jià)值仍需更多的回顧性研究來(lái)確證; 另一方面,筆者研究結(jié)果提示多種免疫細(xì)胞類型組合比單一免疫細(xì)胞類型更有利于患者預(yù)后的精準(zhǔn)評(píng)估。最后,相關(guān)性研究發(fā)現(xiàn)m6A調(diào)控因子與多數(shù)免疫細(xì)胞浸潤(rùn)比例相關(guān),提示這些因子在抗腫瘤免疫中可能扮演著重要角色。例如,最新研究[20]發(fā)現(xiàn)ALKBH5可抑制細(xì)胞核因子κB, 進(jìn)而減少趨化因子配體5表達(dá)水平并促進(jìn)CD8+T細(xì)胞浸潤(rùn),最終抑制CRC細(xì)胞轉(zhuǎn)移。METTL3可以促進(jìn)髓源性抑制細(xì)胞的積聚,進(jìn)而抑制CD4+和CD8+T細(xì)胞的激活和增殖,最終介導(dǎo)CRC細(xì)胞對(duì)免疫檢查點(diǎn)抑制劑的抵抗[21]。

為了進(jìn)一步探究m6A修飾和免疫細(xì)胞浸潤(rùn)在CRC中的預(yù)后價(jià)值,筆者團(tuán)隊(duì)基于WGCNA篩選出與兩者相關(guān)的4個(gè)預(yù)后基因(ITLN1、LY6G6D、ATOH1和MMP28)用于構(gòu)建模型。在既往文獻(xiàn)中, ITLN1高表達(dá)提示CRC患者更有利的預(yù)后,相關(guān)機(jī)制研究[22]表明其可能通過(guò)抑制表皮生長(zhǎng)因子受體、細(xì)胞外調(diào)節(jié)蛋白激酶及蛋白激酶B發(fā)揮促癌效應(yīng)。LY6G6D在CRC組織中選擇性表達(dá)并且可以作為T細(xì)胞免疫治療的潛在靶點(diǎn)[23]。生物信息學(xué)研究[24]表明ATOH1低表達(dá)為影響CRC患者總體生存的獨(dú)立不良預(yù)后因素。MMP28作為CD36的下游靶基因可通過(guò)介導(dǎo)E鈣黏蛋白的裂解促進(jìn)CRC細(xì)胞的轉(zhuǎn)移[25]。在進(jìn)一步驗(yàn)證工作中,一方面基于上述4個(gè)預(yù)后基因的臨床模型可以在不同臨床特征亞組中識(shí)別潛在預(yù)后不良人群,另一方面該臨床模型在4個(gè)GEO驗(yàn)證隊(duì)列中仍然可有效區(qū)分患者的臨床結(jié)局。因此,該臨床模型有望成為幫助臨床醫(yī)師精準(zhǔn)評(píng)估患者預(yù)后的潛在輔助工具。

本研究存在如下局限性有待改進(jìn): 首先,本研究所用驗(yàn)證隊(duì)列均來(lái)源于GEO數(shù)據(jù)庫(kù),缺少本單位及合作單位獨(dú)立回顧性隊(duì)列進(jìn)一步驗(yàn)證。其次, 4個(gè)預(yù)后基因如何通過(guò)調(diào)控m6A修飾及抗腫瘤免疫反應(yīng)影響CRC發(fā)生發(fā)展仍不清楚,有待于后續(xù)基于生命組學(xué)的機(jī)制研究來(lái)闡明。再次,在<60歲及直腸癌亞組分析中,盡管低分險(xiǎn)組OS有優(yōu)于高風(fēng)險(xiǎn)組OS的趨勢(shì),但是統(tǒng)計(jì)學(xué)上無(wú)顯著差異。因而,預(yù)后模型在上述人群中的臨床價(jià)值仍需更多的回顧性研究來(lái)確證。最后,輔助化療是CRC患者術(shù)后的重要治療手段之一,本預(yù)后模型能否識(shí)別輔助化療不敏感的風(fēng)險(xiǎn)人群也是未來(lái)值得研究的方向。

綜上所述,本研究基于m6A修飾及免疫細(xì)胞浸潤(rùn)構(gòu)建了包含4個(gè)預(yù)后基因的風(fēng)險(xiǎn)模型,該模型能夠在訓(xùn)練隊(duì)列及驗(yàn)證隊(duì)列中有效識(shí)別預(yù)后不良的潛在風(fēng)險(xiǎn)人群。上述研究成果不僅推動(dòng)相關(guān)基礎(chǔ)研究領(lǐng)域成果的臨床轉(zhuǎn)化,也為基于分子標(biāo)志物的CRC精準(zhǔn)診療提供了新的思路。