梁德宏，王艷，張全穎，陳燃，謝朝陽，羅海清，李祥勇

(1.廣東醫(yī)科大學(xué)生物化學(xué)與分子生物學(xué)教研室/廣東省教育廳醫(yī)學(xué)生物活性分子研究重點實驗室，廣東湛江524023；2.廣東醫(yī)科大學(xué)醫(yī)學(xué)技術(shù)學(xué)院，廣東東莞523808；3.廣東醫(yī)科大學(xué)附屬醫(yī)院腫瘤醫(yī)院，廣東湛江524001)

腎透明細胞癌(kidney clear cell carcinoma，KIRC)起源于腎小管上皮，是腎細胞癌的主要組織病理亞型，具有更高的腫瘤復(fù)發(fā)率和轉(zhuǎn)移率[1]。早期KIRC可通過手術(shù)或消融策略成功治愈，但晚期、轉(zhuǎn)移性KIRC患者的總體預(yù)后仍然較差[2]。隨著免疫治療的不斷發(fā)展，其應(yīng)用于KIRC的臨床收益獲得持續(xù)提升，但大多數(shù) KIRC患者仍沒有獲得持久的療效收益，甚至出現(xiàn)免疫相關(guān)不良事件，這主要歸因于腫瘤免疫微環(huán)境(tumor microenvironment, TME)的異質(zhì)性和復(fù)雜性[3-4]。研究表明，具有強免疫原性的KIRC，其預(yù)后模型主要基于病理分期等常規(guī)臨床指標(biāo)，為最大化免疫治療的臨床獲益，有必要進一步探索KIRC免疫相關(guān)生物標(biāo)志物和TME的特征[5]。因此，在本研究中，筆者擬通過生物信息學(xué)的方法篩選出免疫相關(guān)的差異表達基因，以此構(gòu)建和驗證KIRC免疫相關(guān)預(yù)后模型，分析TME特征，以期為KIRC預(yù)后管理和免疫治療決策的制定提供重要證據(jù)。

1 材料與方法

1.1數(shù)據(jù)獲取 基于UCSC Xena平臺( https://xenabrowser.net/datapages/)，下載源自癌癥基因組圖譜(The Cancer Genome Atlas, TCGA)的KIRC隊列(593例樣本)和國際癌癥基因組聯(lián)盟(International Cancer Genome Consortium，ICGC)的KIRC隊列(158例樣本)的mRNA表達數(shù)據(jù)以及對應(yīng)患者的臨床信息，包括人口統(tǒng)計學(xué)指標(biāo)(年齡、性別等)、預(yù)后數(shù)據(jù)(生存時間、死亡狀態(tài)等)和病理參數(shù)等。其中，TCGA的KIRC數(shù)據(jù)集作為目標(biāo)數(shù)據(jù)集，ICGC的KIRC數(shù)據(jù)集作為驗證數(shù)據(jù)集。

1.2篩選免疫相關(guān)的差異表達基因 基于R包limma，對目標(biāo)數(shù)據(jù)集進行差異表達基因(differentially expressed genes，DEGs)分析，篩選標(biāo)準(zhǔn)為矯正后P值小于0.05和差異倍數(shù)|Fold Change|大于1.5，即被視為DEGs。將得到的DEGs與從免疫學(xué)數(shù)據(jù)庫和分析平臺(Immunology Database and Analysis Portal，ImmPort)數(shù)據(jù)庫(https://immport.niaid.nih.gov/)獲取的參考免疫相關(guān)基因集取交集，得到免疫相關(guān)的DEGs，用于后續(xù)分析。

1.3構(gòu)建免疫相關(guān)預(yù)后風(fēng)險模型 基于R包Survival，對免疫相關(guān)DEGs與患者總生存時間(overall survival，OS)進行單因子回歸分析，篩選顯著關(guān)聯(lián)(P<0.05)的基因作為候選預(yù)后因子。接著，通過LASSO cox回歸分析，進一步篩選出重要的預(yù)后相關(guān)基因及其回歸系數(shù)，以此建立KIRC預(yù)后風(fēng)險評估模型(風(fēng)險分數(shù)=∑差異基因的回歸系數(shù)λi×歸一化處理后的基因表達量βi)，并利用R包survivalROC繪制ROC曲線來評估預(yù)后模型的性能。以真陽性率與假陽性率之間的差值最大處作為最佳風(fēng)險評分臨界值，將腫瘤樣本分為高風(fēng)險組和低風(fēng)險組，通過R包SUR-Miner 繪制兩組的生存曲線。

1.4腫瘤微環(huán)境的評估分析 采用CIBERSORT算法分析22種腫瘤浸潤免疫細胞在高風(fēng)險組和低風(fēng)險組的比例，以此評估和比較腫瘤免疫浸潤細胞類型的構(gòu)成和組間差異(通過非配對t檢驗進行)[6]。

1.5構(gòu)建和驗證列線圖 基于R包rms，結(jié)合臨床信息(年齡、性別和腫瘤分期)，對預(yù)后風(fēng)險評分模型進行列線圖模型構(gòu)建，以此分析其在不同患者臨床特征的預(yù)后價值。臨床信息、風(fēng)險評分與OS之間的關(guān)系通過森林圖形式展示，其中，一致性指數(shù)(C-index)表明列線圖的預(yù)測準(zhǔn)確性。

1.6統(tǒng)計學(xué)分析 采用R軟件對數(shù)據(jù)進行處理和統(tǒng)計分析，R包survivalROC繪制ROC曲線并評估預(yù)后模型的性能，Wilcoxon檢驗用于檢驗兩組間的差異，生存分析曲線的組間差異采用Log-Rank檢驗。以P<0.05為差異具有統(tǒng)計學(xué)意義。

2 結(jié)果

2.1篩選關(guān)聯(lián)預(yù)后的免疫相關(guān)DEGs 納入TCGA的KIRC隊列為目標(biāo)數(shù)據(jù)集，共計593個樣本(表1)，對腫瘤組織樣本和癌旁組織樣本的mRNA表達譜進行分析，共發(fā)現(xiàn)694個DEGs，其中326個基因表達上調(diào)(圖1A)。將所得DEGs與ImmPort數(shù)據(jù)庫的參考免疫基因集(1 509個)取交集后，篩選出57個免疫相關(guān)的DEGs(圖1B)，進一步通過單因子回歸和LASSO cox回歸分析，得到11個與OS顯著相關(guān)的免疫相關(guān)DEGs(P<0.05，表2)，以及各個基因?qū)?yīng)的回歸系數(shù)(圖1C)。NOS1、ANGPTL3和AVPR1B基因高表達與良好的預(yù)后顯著相關(guān)，而其他8個基因則相反。根據(jù)多因素回歸分析模型，Akaik信息標(biāo)準(zhǔn)(AIC)為705.42，C指數(shù)為0.71。值得注意的是，僅ANGPTL3與OS顯著相關(guān)(P=0.035)，其風(fēng)險率為0.89(95%CI：0.80～099)，見圖1D。

表1 TCGA KIRC數(shù)據(jù)集的樣本信息[n(%)]

表2 與KIRC預(yù)后相關(guān)的基因列表(單因素回歸分析)

注：A，TCGA KIRC數(shù)據(jù)集的差異基因分析，橫坐標(biāo)為差異倍數(shù)的對數(shù)值(以log2為底)，縱坐標(biāo)為校正后P值的對數(shù)值(以 log10為底)的相反數(shù)；B，差異表達基因(DEGs)和ImmPort參考免疫基因集的交集分析；C，基因的回歸系數(shù)柱狀圖；縱坐標(biāo)為基因，橫坐標(biāo)為回歸系數(shù)；D，多因素回歸分析的森林圖；縱坐標(biāo)為風(fēng)險因子(基因)，橫坐標(biāo)為風(fēng)險因子的風(fēng)險率。

2.2免疫相關(guān)預(yù)后模型的搭建 基于11個顯著與OS相關(guān)的免疫相關(guān)DEGs及其回歸系數(shù)構(gòu)建的免疫相關(guān)預(yù)后風(fēng)險模型，將風(fēng)險得分的閾值設(shè)為-0.02，以此將所有KIRC樣本分為高風(fēng)險組和低風(fēng)險組(圖2A)。結(jié)合預(yù)后數(shù)據(jù)分析，發(fā)現(xiàn)高風(fēng)險組的死亡比例顯著高于低風(fēng)險組，提示風(fēng)險評分對評估腎透明細胞癌患者生存狀況具有一定的應(yīng)用價值(圖2B)。值得注意的是，相比高風(fēng)險組，低風(fēng)險組具有更高的總生存時間(P<0.001)，表現(xiàn)出NOS1、ANGPTL3和AVPR1B表達量上調(diào)，而其他8個基因表達量下調(diào)(圖2C～D)。

注：A，高風(fēng)險組和低風(fēng)險組的風(fēng)險曲線，縱坐標(biāo)代表風(fēng)險評分，橫坐標(biāo)代表風(fēng)險評分等級；B，生存狀況分布，縱坐標(biāo)代表生存時間(d)，橫坐標(biāo)代表風(fēng)險評分等級；C，11個免疫相關(guān)DEGs在高風(fēng)險組和低風(fēng)險組的比較；D，總生存時間在高風(fēng)險組和低風(fēng)險組的比較。

2.3免疫相關(guān)預(yù)后模型的性能評估和驗證 繪制用于評估免疫相關(guān)預(yù)后風(fēng)險模型性能的ROC曲線，發(fā)現(xiàn)該模型在評估3年、4年、5年生存概率的AUCROC分別為0.669、0.707和0.750(圖3A)，且風(fēng)險評分在KIRC不同腫瘤分期中表現(xiàn)出顯著的差異，傾向于臨床分期越高，則風(fēng)險評分越高(圖3B)。同時，基于ICGC的KIRC隊列數(shù)據(jù)集，選擇其中具有mRNA表達數(shù)據(jù)和預(yù)后信息的腫瘤樣本共158例，以此作為預(yù)后模型的驗證集，其中男性100例，女性58例，年齡大于60歲共101例，小于60歲共57例。結(jié)果顯示，在驗證集中，該免疫相關(guān)預(yù)后風(fēng)險模型在評估3年、4年、5年生存概率的AUCROC分別為0.628、0.653和0.677(圖3C)，高風(fēng)險評分與較短生存期顯著相關(guān)(P<0.000 1)，見圖3D。

注：A，基于目標(biāo)數(shù)據(jù)集搭建預(yù)后模型用于預(yù)測3/4/5年生存率的ROC曲線；B，目標(biāo)數(shù)據(jù)集中不同臨床分期的腫瘤樣本的風(fēng)險評分比較；C，預(yù)后模型在驗證集中表現(xiàn)出來的預(yù)測3/4/5年生存率的ROC曲線；D，基于預(yù)后模型給驗證集樣本計算的風(fēng)險評分的生存分析曲線。

經(jīng)過多因素回歸模型分析，發(fā)現(xiàn)年齡(≥60歲)、腫瘤分期(Ⅲ和Ⅳ期)和基于免疫相關(guān)預(yù)后模型得到的風(fēng)險評分均與OS顯著相關(guān)，風(fēng)險評分與腫瘤Ⅲ期的風(fēng)險率相當(dāng)，且可通過列線圖估算3年、4年和5年的生存概率(圖4A～B)。通過繪制校準(zhǔn)曲線驗證列線圖的性能，可以看出預(yù)測曲線接近理想曲線(圖4C)，表明構(gòu)建的列線圖模型性能良好，其預(yù)測準(zhǔn)確性較先前的免疫相關(guān)預(yù)后模型得到進一步提高(C指數(shù)：0.78 vs 0.71)。值得注意的是，風(fēng)險評分對生存時間的風(fēng)險率(1.93，95%CI：1.59～2.3)比任一單個免疫相關(guān)DEGs都高。

注：A，基于目標(biāo)數(shù)據(jù)集，將患者年齡、性別、病理分期、風(fēng)險評分分別作為變量預(yù)測患者3年、4年和5年生存率的列線圖；B，多因素回歸分析的森林圖；C，分別針對3年、4年和5年的預(yù)后模型繪制的校準(zhǔn)曲線，灰色直線代表理想狀態(tài)，藍色折線代表實際狀態(tài)。

2.4探索KIRC腫瘤免疫微環(huán)境 基于mRNA表達譜，通過CIBERSORT算法來估算每例KIRC樣本中22種免疫細胞的比例，結(jié)果顯示，12種免疫細胞類型在高風(fēng)險組和低風(fēng)險組間的占比差異具有統(tǒng)計學(xué)意義(P<0.05)，其中，漿細胞(plasma cell)、調(diào)節(jié)T細胞(regulatory T cell，Tregs)、單核細胞(Monocyte)等6種細胞的組間差異最為顯著(P<0.000 1，圖5A)。值得注意的是，在高風(fēng)險組中表現(xiàn)為更高浸潤程度的調(diào)節(jié)T細胞、M0型巨噬細胞(macrophages M0)、活化的記憶CD4+T細胞(T cells CD4 memory activated)與較差的預(yù)后相關(guān)，同時，在高風(fēng)險組中浸潤程度更低的休眠樹突狀細胞(dendritic cells resting)則與較好的預(yù)后相關(guān)(圖5B～E)。

注：A，高風(fēng)險組(藍色)和低風(fēng)險組(紅色)免疫細胞浸潤程度的比較，縱坐標(biāo)為細胞比例，橫坐標(biāo)為免疫細胞類型；B，M0型巨噬細胞浸潤程度高低在生存時間上的比較分析；C，調(diào)節(jié)T細胞浸潤程度高低在生存時間上的比較分析；D，活化的記憶CD4+T細胞浸潤程度高低在生存時間上的比較分析；E，休眠樹突狀細胞浸潤程度高低在生存時間上的比較分析。

此外，進一步對高風(fēng)險組和低風(fēng)險組間免疫檢查點及其配體的表達水平進行分析，結(jié)果表明，相比于低風(fēng)險組，高風(fēng)險組中PD-1、CTLA-4和LAG3的表達水平顯著上調(diào)(P<0.000 1),而PD-L1顯著下降(圖6A～D)。

注：A，PD-1的表達量與不同風(fēng)險評分等級的比較分析；B，PD-L1的表達量與不同風(fēng)險評分等級的比較分析：C，CTLA-4的表達量與不同風(fēng)險評分等級的比較分析；D，LAG3的表達量與不同風(fēng)險評分等級的比較分析。

3 討論

為探索免疫相關(guān)分子標(biāo)志物在KIRC預(yù)后中的預(yù)測價值，筆者從TCGA的KIRC數(shù)據(jù)集中篩選出11個顯著關(guān)聯(lián)KIRC患者總生存時間的免疫相關(guān)DEGs。其中包括鮮有被報道的ANGPTL3，在本研究中無論單因素還是多因素回歸分析結(jié)果，均表現(xiàn)出與KIRC預(yù)后顯著的關(guān)聯(lián)性，且其高表達被證實可抑制腫瘤的增殖和轉(zhuǎn)移，進而關(guān)聯(lián)更好的預(yù)后[7-8]。此外，單因素回歸分析結(jié)果顯示，NOS1和AVPR1B基因的高表達均與良好的預(yù)后相關(guān)，但在多因素回歸分析中不存在相關(guān)性，表明這2個基因?qū)IRC預(yù)后的預(yù)測價值仍需進一步探索。

鑒于單個免疫相關(guān)DEGs對預(yù)后結(jié)果的影響程度有限，筆者以這11個免疫相關(guān)DEGs構(gòu)建了預(yù)后模型，并通過ICGC的KIRC數(shù)據(jù)集進行驗證。結(jié)果顯示，基于該模型得到的風(fēng)險評分越高，腫瘤分期等級越高，生存時間越短，且風(fēng)險評分在評估3年、4年、5年生存概率上也表現(xiàn)出一定的預(yù)測性能。同時，風(fēng)險評分對KIRC預(yù)后的影響程度相比于單個免疫相關(guān)DEGs更強(風(fēng)險率：1.93 vs 1.37)，且與腫瘤分期Ⅲ期相當(dāng)(風(fēng)險率：1.93 vs 1.97)，這進一步支持該模型的臨床應(yīng)用潛力。

鑒于KIRC中TME異質(zhì)性和復(fù)雜性可影響免疫治療的臨床獲益，故而筆者進一步分析了風(fēng)險評分與TME的關(guān)聯(lián)。結(jié)果表明，組間顯著差異的免疫細胞類型高達12種，其中，漿細胞、調(diào)節(jié)T細胞和M0型巨噬細胞在高風(fēng)險組中表現(xiàn)出更高的浸潤程度，且均與較差的預(yù)后相關(guān)。其中，調(diào)節(jié)T細胞被認為是KIRC中有效抗腫瘤免疫的主要障礙，可通過負調(diào)節(jié)腫瘤內(nèi)部或外部的免疫和非免疫細胞，限制抗腫瘤免疫反應(yīng)，進而不利于KIRC患者的預(yù)后結(jié)果[9]。同時，相比低風(fēng)險組，在高風(fēng)險組中，浸潤程度更高的活化的記憶CD4+T細胞和M0型巨噬細胞與較差的預(yù)后相關(guān)，而浸潤程度較低的休眠樹突狀細胞，則與較好的預(yù)后相關(guān)，這與先前的研究結(jié)論相一致[10]。此外，筆者還分析了免疫檢查點及其配體的表達情況，發(fā)現(xiàn)高風(fēng)險組中PD-1、CTLA-4和LAG3顯著上調(diào)，而PD-L1顯著下調(diào)。以上結(jié)果表明與風(fēng)險評分相關(guān)聯(lián)的調(diào)節(jié)T細胞、休眠樹突狀細胞、記憶CD4+T細胞、M0型巨噬細胞和免疫檢查點可能是KIRC的獨立預(yù)后因素，有望為預(yù)后管理和免疫治療提供證據(jù)。

本研究存在以下局限性，一是研究數(shù)據(jù)均源自TCGA數(shù)據(jù)庫，樣本來源以白人為主，該結(jié)果是否適用于黃種人尚需驗證；二是研究數(shù)據(jù)的來源人群以60歲以上的男性為主，或許僅限這小部分人群可得到結(jié)果復(fù)現(xiàn)和應(yīng)用；三是該模型的應(yīng)用價值還需要進一步的臨床研究驗證。

綜上所述，本研究聚焦于KIRC患者群體，基于11個關(guān)聯(lián)OS的免疫相關(guān)DEGs構(gòu)建和驗證了1個性能良好的KIRC免疫相關(guān)預(yù)后風(fēng)險模型，并在此基礎(chǔ)上，探索不同風(fēng)險評分下的TME特征，有望為KIRC患者的臨床風(fēng)險分層和免疫治療決策提供支持。