基于細胞焦亡相關(guān)基因表達譜的子宮內(nèi)膜癌預(yù)后模型構(gòu)建*

常 晴,劉 佳,曲愛林,劉京康

(山東大學齊魯醫(yī)院a.檢驗科;b.基礎(chǔ)醫(yī)學研究中心;c.婦產(chǎn)科,濟南 250012)

子宮內(nèi)膜癌(uterine corpus endometrial carcinoma,UCEC)是婦科最常見的惡性腫瘤之一,近年來其發(fā)病呈年輕化趨勢,發(fā)病率和死亡率逐年上升[1]。早期診斷人群可通過手術(shù)取得較好的預(yù)后,5年生存率可達94.9%[2]。然而,晚期患者即使經(jīng)放射治療(放療)或化學治療(化療),術(shù)后復(fù)發(fā)風險仍很高[3-4]。經(jīng)陰道超聲結(jié)合子宮內(nèi)膜組織活檢是診斷UCEC的主要方法,探索新的預(yù)后有關(guān)生物標志物并與之聯(lián)合,將有助于提高UCEC的早期診斷率及制定合理的個體化治療方案。

近年來,細胞焦亡在腫瘤發(fā)生發(fā)展過程中的作用備受關(guān)注。焦亡是由炎癥介導的一種程序性細胞死亡方式,其特征是細胞腫脹、胞膜破裂、染色質(zhì)斷裂和促炎因子級聯(lián)反應(yīng)[5]。大量研究表明,惡性腫瘤中異常表達的焦亡基因可作為診斷、預(yù)后和靶向治療的分子標志物[6-7]。Gao等[8]發(fā)現(xiàn),Gasdermin D(GSDMD)是細胞焦亡的關(guān)鍵執(zhí)行蛋白,與正常組織相比,GSDMD在肺癌中高表達,并且其表達水平與患者預(yù)后呈負相關(guān)。Liu等[9]研究發(fā)現(xiàn),焦亡相關(guān)的炎性小體NOD樣受體熱蛋白結(jié)構(gòu)域相關(guān)蛋白3(NOD-like receptor thermal protein domain associated protein 3,NLRP3)在子宮內(nèi)膜癌中的表達顯著高于正常組織,其通過調(diào)節(jié)caspase-1表達,增加HEC-1A細胞的增殖、侵襲和遷移能力。Wang等[10]報道,焦亡在腫瘤微環(huán)境中可觸發(fā)強烈的免疫反應(yīng),參與抗腫瘤免疫的調(diào)控。此外,細胞焦亡與化療藥物的應(yīng)用密切相關(guān),如5-氟尿嘧啶(5-FU)活化腫瘤細胞內(nèi)的caspase-3,激活GSDME-NT的生成,從而誘導胃癌細胞焦亡[11]。以上研究表明,細胞焦亡在癌癥進展中具有促進腫瘤發(fā)生和激活抗腫瘤免疫應(yīng)答的雙重作用。同時,其對于化療、免疫治療等精準醫(yī)療方案的實施可能具有一定指示作用?；诖?借助生物信息學方法,通過分析公共數(shù)據(jù)庫,挖掘UCEC焦亡相關(guān)基因,構(gòu)建預(yù)后風險模型,初步探究風險評分與患者預(yù)后和免疫浸潤的關(guān)系。同時,將風險評分模型與臨床病理參數(shù)進行整合,構(gòu)建可視化的Nomogram模型,方便臨床應(yīng)用。

1 材料與方法

1.1 數(shù)據(jù)來源 從UCSC xena(https://xenabrowser.net/datapages/)中的TCGA-UCEC隊列下載子宮內(nèi)膜癌轉(zhuǎn)錄組RNA seq數(shù)據(jù)和生存信息,導入R軟件(版本4.0.2),將數(shù)據(jù)讀入并進行歸一化處理。使用分子特征數(shù)據(jù)庫(MSigDB)(https://www.gsea-msigdb.org/gsea/msigdb)和PubMed數(shù)據(jù)庫相關(guān)文獻報道(PubMed號:33828074和34992421)查找焦亡相關(guān)表型基因集,獲得22種焦亡基因納入后續(xù)分析。從cbioportal 數(shù)據(jù)庫(www.cbioportal.org)選擇癌癥類型為子宮內(nèi)膜癌并點擊TCGA-PanCancer Atlas模塊,獲得509例具有突變信息的UCEC樣本,進行后續(xù)分析。利用STRING數(shù)據(jù)庫(https://cn.string-db.org)獲得PPI互作網(wǎng)絡(luò)關(guān)系。

1.2 方法

1.2.1 焦亡相關(guān)預(yù)后模型的構(gòu)建 排除臨床數(shù)據(jù)和預(yù)后信息缺失的樣本,共納入533例有臨床病理信息和生存數(shù)據(jù)的UCEC樣本。其中,Ⅰ期331例,Ⅱ期50例,Ⅲ期123例,Ⅳ期29例。根據(jù)病灶部位及病情嚴重程度,排除沒有相關(guān)數(shù)據(jù)的樣本,最終分選出3種類型的UCEC樣本,即11例原發(fā)型UCEC,37例復(fù)發(fā)型UCEC和34例轉(zhuǎn)移型UCEC。根據(jù)病理分級分為G1級109例、G2級117例和G3級307例。使用R包glmne進行最小絕對收縮選擇算子(LASSO)回歸篩選焦亡相關(guān)預(yù)后基因,得出交叉驗證誤差最小的模型。預(yù)測模型的風險評分=Exp1*β1+ Exp2*β2+ Exp3*β3+…… +Expx*βx,其中Exp為相對應(yīng)預(yù)后的基因表達量,β為lasso回歸系數(shù)。使用R包“timeROC”構(gòu)建受試者工作曲線(ROC)評估風險評分模型預(yù)測能力的準確性。將風險評分和臨床病理參數(shù)(年齡、病理分期)納入以構(gòu)建預(yù)測子宮內(nèi)膜癌預(yù)后的Nomogram, 使用“cox.step”函數(shù)篩選出最佳變量組合用于建模。應(yīng)用C指數(shù)、受試者工作曲線(receiver operating character,ROC曲線)和校準曲線分析法(calibration)對該模型進行評估。

1.2.2 高風險組和低風險組預(yù)后和免疫浸潤比較 將風險評分的中位值作為截斷值,把子宮內(nèi)膜癌患者分為高風險組和低風險組,使用 R軟件的survival 包和 survminer 包繪制 Kaplan-Meier 生存曲線,評估兩組患者的預(yù)后差異?；赥CGA-UCEC隊列中的RNA seq數(shù)據(jù),通過ssGSEA算法計算樣本的28種免疫細胞豐度,比較高風險組和低風險組免疫細胞浸潤程度的差異。

1.3 統(tǒng)計學處理 使用R軟件(版本4.0.2)進行數(shù)據(jù)分析和圖形繪制。兩組間比較使用Wilcoxon秩和檢驗。生存曲線(Kaplan-Meier)使用log-rank檢驗。P<0.05為差異有統(tǒng)計學意義。

2 結(jié) 果

2.1 TCGA-UCEC隊列中22種焦亡相關(guān)基因轉(zhuǎn)錄組特征 根據(jù)MSigDB數(shù)據(jù)庫和PubMed數(shù)據(jù)庫相關(guān)文獻報道,共獲得22種焦亡相關(guān)基因, TCGA-UCEC隊列中包含了533例子宮內(nèi)膜癌樣本和35例癌旁樣本,導入R軟件對22種焦亡相關(guān)基因mRNA進行表達差異分析。結(jié)果顯示,13種基因在UCEC組織中的表達高于癌旁組織(P<0.05),3種基因在UCEC組織中的表達低于癌旁組織(P<0.05),6種基因在癌與癌旁組織中表達無差異(P>0.05)。見圖1A。

為了明確22種焦亡相關(guān)基因的表達水平是否受到轉(zhuǎn)錄組突變的影響,通過cbioportal數(shù)據(jù)庫對22種基因進行了突變分析。結(jié)果顯示,具有突變信息的509例子宮內(nèi)膜癌患者中有288例患者存在基因突變(57%),其中,TP53基因發(fā)生突變的樣本例數(shù)(194例)占總樣本(509例)的38%,在22種基因中突變頻率最高,突變類型主要為錯義突變(missense mutation)和截斷突變(truncating mutation),其次是CASP8基因,突變頻率為12%,提示焦亡相關(guān)基因的轉(zhuǎn)錄組變異可能影響基因的表達水平(圖1B)。進一步使用STRING數(shù)據(jù)庫構(gòu)建22種焦亡相關(guān)基因的蛋白質(zhì)-蛋白質(zhì)相互作用網(wǎng)絡(luò)(圖2),發(fā)現(xiàn)22種基因在蛋白質(zhì)水平上存在密切聯(lián)系。

圖1 TCGA-UCEC隊列中22種焦亡相關(guān)基因轉(zhuǎn)錄組特征

圖2 22種焦亡相關(guān)基因的蛋白質(zhì)-蛋白質(zhì)相互作用網(wǎng)絡(luò)

2.2 構(gòu)建預(yù)測UCEC預(yù)后的風險評估模型 對癌與癌旁組織有差異表達的16種焦亡基因進行進一步分析,篩選具有預(yù)后價值的相關(guān)基因。為了避免模型的過度擬合,通過LASSO回歸分析,基于lambda.min方法確定最佳lambda值為9(圖3A、B),共獲得9種與預(yù)后有關(guān)的焦亡基因(CASP3、CASP5、CYCS、GPX4、GZMB、IL1A、IRF2、TP53 和TNF)。計算每例UCEC患者的風險評分。風險評分=(-0.0169*CASP3表達量)+(0.7129* CASP5表達量)+(0.2239*CYCS表達量)+(-0.2027* GPX4表達量)+(-0.1743*GZMB表達量)+(-0.0185* IL1A表達量)+(-0.3798* IRF2表達量)+(-0.0752* TP53表達量)+(0.0631*TNF表達量)。

隨著風險評分的增加,患者生存時間縮短,死亡人數(shù)明顯增多(圖3C、D)。提示風險評分可能是預(yù)測患者預(yù)后的危險因素。通過時間依賴性ROC曲線評價風險評分模型預(yù)測效能,ROC曲線的1、3、5年AUC值分別為0.627、0.722、0.770,表明該風險評分對UCEC患者具有良好的預(yù)后預(yù)測能力(圖3E)。

圖3 LASSO回歸構(gòu)建UCEC預(yù)后風險的預(yù)測模型A:LASSO系數(shù)分布特征;B:選擇LASSO模型中最佳調(diào)整參數(shù) λ (lambda);C:風險評分和總體生存期的分布特征;D:不同生存狀態(tài)下的風險評分分布特征(生存狀態(tài):0代表存活;生存狀態(tài):1:代表死亡);E:時間依賴性ROC曲線評估風險評分模型的預(yù)后價值

2.3 高低風險組中臨床病理分型、預(yù)后特征和免疫細胞浸潤比較 將UCEC患者風險評分的中位值作為截斷值,將533例UCEC患者分為高風險組(261例)和低風險組(272例)。分析不同臨床病理分型在兩組間的分布特征,Ⅱ期、Ⅲ期和Ⅳ期UCEC樣本在高風險組中的比例高于低風險組,而Ⅰ期UCEC樣本在高風險組中的比例低于低風險組。原發(fā)型UCEC和轉(zhuǎn)移型UCEC在高風險組的比例低于低風險組,而復(fù)發(fā)型UCEC樣本在高風險組的占比高于低風險組。此外,對腫瘤細胞的分化程度進行分級,結(jié)果顯示低風險組中G1級和G2級樣本占比高于高風險組,G3級在低風險組占比低于高風險組。見表1。Kaplan-Meier生存曲線顯示,高風險組和低風險組生存率存在顯著差異,低風險組患者總生存期顯著高于高風險組(P<0.001),預(yù)后較好(圖4A)。

表1 不同臨床病理分型在高低風險組間的占比[n(%) ]

使用ssGSEA對兩組進行免疫浸潤分析,發(fā)現(xiàn)與高風險組相比,低風險組中激活型B細胞、激活型CD8+T細胞、激活型樹突狀細胞、CD56 bright NK細胞、CD56 dim NK細胞、中央記憶型CD4+T細胞、效應(yīng)記憶型CD8+T細胞、嗜酸性粒細胞、未成熟B細胞、巨噬細胞、肥大細胞、髓源性抑制細胞、單核細胞、自然殺傷細胞、自然殺傷T細胞、中性粒細胞、輔助性T細胞、濾泡輔助性T細胞、Ⅰ型輔助性T細胞和17型輔助性T細胞的免疫細胞亞群顯著上調(diào)(P<0.05),效應(yīng)記憶型CD4+T細胞在高風險組中的比例高于低風險組(P<0.05),激活型CD4+T細胞、中央記憶型CD8+T細胞、γδT細胞、未成熟樹突狀細胞、漿細胞樣樹突狀細胞、調(diào)節(jié)性T細胞和輔助型T細胞2在兩組中占比無差異(P>0.05)。見圖4B。

圖4 高低風險組間預(yù)后特征和免疫細胞浸潤A:高風險組與低風險組間總體生存期比較;B:單樣本基因集富集分析比較高低風險組之間的免疫細胞浸潤

2.4 構(gòu)建Nomogram模型 通過逐步回歸方法確定納入Nomogram模型的最佳變量,最終確定風險程度、年齡和臨床病理分期3個變量構(gòu)建預(yù)后Nomogram(圖5A)。Nomogram預(yù)后模型的C指數(shù)(C-index)為0.681(95%CI為0.618～0.742,P<0.001),1、3、5年總體生存率ROC曲線的AUC值分別為0.625、0.701、0.762(圖5B),表明構(gòu)建的Nomogram預(yù)測模型具有較高的準確性。預(yù)測1、3、5年總生存率的校準曲線(calibration)同樣顯示Nomogram預(yù)后模型與實際觀測結(jié)果一致性較高,提示Nomogram的預(yù)測效能良好。見圖5C～E。

圖5 構(gòu)建UCEC患者預(yù)后Nomogram模型A:預(yù)測UCEC患者預(yù)后的Nomogram;B:時間依賴性ROC曲線驗證Nomogram的預(yù)后價值;C～E:分別為Nomogram模型1年、3年、5年的校準曲線

3 討 論

焦亡不僅導致細胞死亡,而且會引發(fā)一系列炎癥反應(yīng)[12]。細胞壞死即可誘發(fā)大量炎性細胞因子的釋放,長期暴露于炎性環(huán)境能增加腫瘤發(fā)生的風險[13]。本研究中,從TCGA數(shù)據(jù)庫下載533例UCEC樣本,通過比較焦亡基因在腫瘤和正常組織中的表達水平及進行Lasso回歸分析,獲得了9個預(yù)后相關(guān)因子,繼而構(gòu)建了風險評分模型,并初步探究了風險評分與患者預(yù)后和免疫浸潤的關(guān)系,進一步結(jié)合其他臨床病例參數(shù)構(gòu)建UCEC焦亡相關(guān)基因預(yù)后模型。

本研究從22個焦亡基因中篩選出反映UCEC患者預(yù)后的9個焦亡關(guān)鍵基因,包括CASP3、CASP5、CYCS、GPX4、GZMB、IL-1A、IRF2、TP53和TNF,其中大部分在其他研究中已有相關(guān)報道。如CASP3編碼的半胱氨酸蛋白酶-3(caspase-3)在食管鱗狀細胞癌中表達上調(diào),且caspase-3陽性表達患者的預(yù)后良好,5年生存率可達25%,而caspase-3陰性患者的5年生存率僅為6%[14]。谷胱甘肽過氧化物酶4 (GPX4) 是一種保護脂質(zhì)過氧化的抗氧化酶,Kang等[15]發(fā)現(xiàn),其在免疫細胞中表達上調(diào),可對抗GSDMD介導的細胞焦亡。顆粒酶B (GzmB)是主要由自然殺傷細胞和細胞毒性T細胞分泌的蛋白質(zhì),GzmB通過直接裂解gasdermin 家族成員GSDME和激活caspase-3誘導腫瘤細胞焦亡,從而達到抗腫瘤作用[16]。Zhang 等[17]發(fā)現(xiàn),特異性敲除TP53會阻斷脂多糖 (LPS) 誘導的肺癌細胞的焦亡。以上研究提示,這些焦亡相關(guān)基因可能在腫瘤發(fā)生發(fā)展中扮演重要角色。目前CYCS、IL-1A、IRF2在UCEC中的作用尚無有關(guān)報道,本研究發(fā)現(xiàn),其與UCEC的預(yù)后顯著相關(guān),推測這些分子可能是預(yù)測子宮內(nèi)膜癌患者預(yù)后的潛在生物標記物。

以往研究發(fā)現(xiàn)焦亡相關(guān)炎性小體在腫瘤細胞中廣泛分布,如NLRP1、NLRP3和CARD等,炎性小體可通過調(diào)節(jié)半胱氨酸蛋白酶-1(caspase-1)觸發(fā)腫瘤細胞的焦亡[18]。細胞焦亡能導致炎癥小體活化及炎癥介質(zhì)釋放,后者可招募免疫細胞浸潤并殺傷腫瘤細胞。本研究結(jié)果顯示,具有抗腫瘤免疫反應(yīng)功能的CD8+T細胞、樹突狀細胞、CD56 bright NK細胞、CD56 dim NK細胞等免疫細胞亞群在低風險組上調(diào),提示低風險組免疫微環(huán)境表現(xiàn)為“熱腫瘤”,具有較高的免疫原性。同時,預(yù)后風險模型可在一定程度上反映免疫細胞浸潤情況,并可能為免疫治療提供新的切入點。

Cao等[19]發(fā)現(xiàn)了與胃癌臨床病理特征、預(yù)后和免疫微環(huán)境浸潤顯著相關(guān)的9個焦亡基因標記,通過無監(jiān)督聚類方法確定了3種不同的焦亡基因亞型,患者免疫微環(huán)境狀態(tài)、免疫治療反應(yīng)和預(yù)后特征在3種亞型之間有顯著差異。在肺腺癌中,Liu等[20]構(gòu)建了與焦亡相關(guān)的3種lncRNAs的預(yù)后模型,此模型與患者的總生存期、免疫細胞浸潤豐度及免疫治療反應(yīng)密切相關(guān)。本研究構(gòu)建風險評分模型,根據(jù)風險評分中位值將UCEC樣本分為高風險組和低風險組,結(jié)果顯示,低風險組患者的總生存期顯著高于高風險組,而高風險組表現(xiàn)出更差的臨床分期、病理分級和更高的復(fù)發(fā)風險。為了方便臨床應(yīng)用,將此預(yù)后模型進一步可視化為Nomogram,應(yīng)用ROC曲線和校準預(yù)測曲線證實Nomogram具有良好的預(yù)后預(yù)測效能。

本研究基于TCGA大樣本多中心的數(shù)據(jù)構(gòu)建了Nomogram預(yù)后模型,確保了該模型的穩(wěn)定性,但仍存在一定局限性,如子宮內(nèi)膜癌相關(guān)臨床指標(病理分期TNM及腫瘤標志物CA125和CA199)未納入TCGA數(shù)據(jù)庫,因此無法評估以上臨床指標對UCEC預(yù)后的影響。此外,焦亡相關(guān)基因?qū)CEC的預(yù)后價值及潛在調(diào)控作用尚需臨床和基礎(chǔ)實驗加以佐證。