基于细胞毒性T淋巴细胞相关基因肝细胞癌预后模型的构建及功能分析

叶必成，缪夏晔，刘树青

1徐州医科大附属淮安医院消化内科，江苏淮安223001；2扬州大学医学院

肝癌全球发病率在恶性肿瘤中居第六位，是肿瘤相关死亡的第四大原因［1］。肝细胞癌（HCC）是最常见的原发性肝癌，早期HCC患者的主要治疗手段为外科手术切除，然而大部分患者在诊断时已为晚期，其对应的治疗手段非常局限，因此预后极差。对于HCC患者预后的预测及治疗方案的确定，临床通常用AJCC分期、肝功能指标及AFP进行评估［2-3］，然而相同分期的HCC患者预后也可能不同，且分子遗传学也证明不同亚群的HCC患者预后有差异［4］。因此，有必要探寻一种新型的生物学标志物以有效评估HCC患者预后。CTL作为T细胞的重要组成部分，其表面受体主要为CD8，是抵抗肿瘤进展的主要免疫细胞［5］。先前一项研究通过规律间隔成簇短回文重复序列在小鼠癌细胞株（肾癌、乳腺癌、黑色素细胞瘤、结直肠癌）中确定了182个与CTL相关的基因（CRG），CRG增强或减弱了CTL对癌细胞的杀伤［6］。然而，CRG在HCC中的作用尚未被系统研究。本研究旨在观察CRG在HCC的表达情况以及探索其与预后的关系，并基于这些基因构建预后预测模型，以期有助于评估患者预后及临床精准化治疗。

1 资料与方法

1.1 资料来源 从TCGA（https：//portal.gdc.cancer.gov/repository）下载371例HCC患者共424个样本（包含50个癌旁样本）的资料，其中3例HCC组织重复测序2次。下载时间：2017年6月—2021年1月。排除生存时间为0患者6例，共纳入HCC患者365例。患者男246例、女119例；年龄16～90岁，中位年龄61岁；肿瘤分级G1级55例、G2级175例、G3级118例、G4级12例、未知5例；肿瘤分期Ⅰ期170例、Ⅱ期84例、Ⅲ期83例、Ⅳ期4例、未知24例；存在血管浸润106例、无血管浸润205例、未知54例；AFP≤200 ng/mL 201例、AFP＞200 ng/mL 75例、AFP未知89例。

从ICGC（https：//dcc.icgc.org/projects/LIRI-JP）下载231例HCC患者的RNA-seq数据和临床资料，用作模型验证。下载时间：2019年3月—2021年1月。该队列患者未存在重复测序以及生存资料缺失，既往有乙型肝炎或丙型肝炎病史。患者男170例、女61例；年龄31～89岁，中位年龄69岁；肿瘤分期Ⅰ期36例、Ⅱ期105例、Ⅲ期71例、Ⅳ期19例。对于重复测序的HCC组织，取其基因表达谱数据的平均值。若有多个探针对应同一个基因，则该基因的表达量为这些探针的均值。以上数据均为公开数据。

1.2 差异表达CRG的筛选 应用“base”R软件包中的“wilcox.test”函数对371例HCC组织和50例癌旁组织的CRG进行表达差异分析，筛选标准：错误发 现 率（FDR）＜0.05且|log2Fold Change|≥1，FoldChange为差异倍数。

1.3 与总生存期（OS）相关CRG的筛选 将具有生存信息的365例HCC患者与各个CRG的表达数据进行合并，并用“base”R软件包中的“coxph”函数对CRG进行基于OS以及生存状态的单因素回归分析，计算各个CRG的P值以及风险比（HR），P＜0.05为与OS相关CRG。

1.4 预后预测模型的构建与验证 利用R软件“base”包中的“intersect”函数确定既为差异表达又与OS相关的CRG，并使用“glmnet”R软件包中的“glmnet”函数和“cv.glmnet”函数对这些基因进行基于OS以及生存状态的Lasso-Cox回归分析，其中maxit=1 000。计算经Lasso-Cox回归所确定的各个基因的系数，并乘以这些基因对应的表达水平，所得结果之和为风险评分。以风险评分的中位数为截断值，将TCGA队列以及ICGC队列的患者分为高危组和低危组。利用“survival”包的“survdiff”对高、低危组患者进行基于Kaplan-Meier法的生存分析，并进行Log-rank检验评估两组的预后差异。利用“Rtsne”的“prcomp”和“Rtsne”函数对Lasso-Cox回归所确定的基因进行PCA及t-SNE分析，以探索高、低危组患者分布情况。利用“timeROC”R软件包的“timeROC”函数进行时间相关的ROC曲线分析，计算曲线下面积（AUC）。

1.5 基因本体（GO）富集分析、京都基因与基因组百科全书（KEGG）富集分析及单样本基因集富集分析（ssGSEA） 利用“clusterProfiler”R软件包的“enrichGO”函数以及“enrichKEGG”函数对TCGA队列的高、低危组之间的差异表达基因进行GO以及KEGG富集分析，其中GO富集分析包括生物学过程（BP）、细胞组成（CC）以及分子功能（MF）。利用“gsva”R软件包的“gsva”函数对TCGA队列的高、低危组进行ssGSEA，以量化免疫细胞浸润及免疫相关功能激活情况。

1.6 统计学方法 采用R4.0.2软件。差异表达基因（DEG）的鉴定采用Wilcoxon检验，并用BH法对P值进行矫正，筛选标准为FDR＜0.05且|log2 Fold Change|≥1。如上文无特殊规定，则P＜0.05为差异有统计学意义。

2 结果

2.1 与OS相关的差异表达CRG 共筛选出37个与OS相关的差异表达DRG。见图1。

图1 与OS相关差异表达CRG的筛选

2.2 基于TCGA队列建立预后模型 Lasso-Cox回归将上述37个DRG进行再次筛选，并使用交叉验证建立模型，最终建立了基于6个CRG（ATG10、HDAC1、PIGU、AHSA1、CAD、CEP55）的预后预测模型，其风险评分的计算公式为0.384×ATG10表达值+0.076×HDAC1表 达值+0.270×PIGU 表 达 值+0.153×AHSA1表达值+0.208×CAD表达值+0.169×CEP55表达值。

2.3 预后预测模型的外部验证结果 上述公式计算ICGC队列的风险评分，并以风险评分的中位数作为截断值将ICGC队列的患者分为高危组和低危组。考虑到OS超过5年的患者仅2例，因此仅评估该模型对患者术后1～4年生存率的预测价值。ROC曲线显示，该模型对患者术后1～4年预后预测的曲线下面积分别为0.70、0.70、0.74、0.74。

2.4 GO富集分析、KEGG富集分析及ssGSEA结果

2.4.1 GO富集分析结果 差异表达的DRG在淋巴细胞介导免疫、补体激活、细胞吞噬、离子通道活性相关分子功能等生物学过程富集。见表1。

表1 GO富集分析结果（BP、CC及MF的前5位）

2.4.2 KEGG富集分析结果 差异表达的DRG基因在细胞因子-细胞因子受体的相互作用、细胞周期及部分癌症通路显著富集。见表2。

表2 KEGG富集分析（前10位）

2.4.3 ssGSEA结果 高危组活化的树突细胞（aDC）、未成熟的树突细胞（iDC）、浆细胞样树突细胞（pDC）、抗原提呈细胞共刺激（APC co-stimulation）、主要组织相容性复合体Ⅰ型分子（MHC class I）富集评分（0.62、0.61、0.72、0.98分）高于低危组（中位数分别为0.60、0.60、0.68、0.97分），比较差异有统计学意义（P均＜0.05）；高危组Ⅰ型干扰素反应（Type I IFN Reponse）、Ⅱ型干扰素反应（TypeⅡIFN Reponse）、NK细胞富集评分（0.82、0.77、0.60分）较低危组（0.83、0.80、0.62分）低，高危组调节性T细胞（Treg）富集评分（0.81分）较低危组（0.80分）高（P＜0.05），比较差异有统计学意义（P＜0.05）。

3 讨论

越来越多的证据表明，肿瘤微环境是影响肿瘤发生发展的重要因素［7］。肿瘤微环境主要由癌细胞、基质细胞、免疫细胞、趋化因子和细胞因子共同构成［8］，其中CTL是抵抗肿瘤进展的主要免疫效应细胞。然而，由于癌细胞不具有很强的免疫原性，因此CTL在肿瘤微环境中是被抑制的［9］。目前已开发了诸如免疫检查点抑制剂（PD-1、PDL-1和CTLA4）等免疫疗法用于促进CTL特异性免疫应答，然而对于HCC患者，PD-1及PDL-1的有效应答率却不到20%［10］。因此识别HCC中影响CTL功能的关键基因至关重要。本研究综合分析182个CRG，最终筛选了6个与OS显著相关的关键CRG，并基于CRG构建了HCC预后预测模型。该模型对患者术后1～4年预后预测的曲线下面积均≥0.70。提示该模型具有良好的预测效能。

本研究发现，近30%的CRG在肿瘤组织与癌旁组织之间差异表达，单因素Cox回归分析表明超过50%的CRG与OS相关，这进一步证明了CTL在HCC发生发展中起重要作用。最终筛选出6个与预后相关的关键CRG，即ATG10、HDAC1、PIGU、AHSA1、CAD、CEP55。除AHSA1和CAD以外，均存在相关研究报道，这些基因在HCC中高表达，并且与OS相关，与本研究的发现一致。ATG10，位于5q14.1，是细胞自噬启动的必要条件［11］。JO等［12］报道ATG10表达增加与血管侵袭和淋巴结转移相关。HDAC1是组蛋白去乙酰化酶家族成员之一，直接参与调控细胞自噬［13］。HDAC1的高表达促进了HCC细胞增殖并抑制HCC细胞凋亡［14］。PIGU是GPI-T复合物的一个重要亚基［15］，GPI-T复合物在肿瘤的发生发展中起至关重要的作用。WEI等［16］的研究发现，PIGU高表达不但促进了HCC细胞的增殖，迁移和侵袭，而且抑制了细胞凋亡。此外，他们还发现，PIGU的下调显著增强了NK-92细胞对HCC细胞的敏感性。AHSA1为HSP90的一个关键分子伴侣，参与致癌蛋白的成熟、稳定、激活［17］，其与HSP90构成的复合体也被证实是细胞自噬通路的关键蛋白［18］。本研究首次发现，AHSA1在HCC中高表达，并且与OS显著相关，其在HCC的具体生物学作用尚未明确，我们将进一步研究证实。CAD是一个具有三个酶结构域的多肽，包括氨基甲酰基磷酸合成酶、天冬氨酸转氨甲酰酶及二氢乳清酸酶，为嘧啶从头合成的关键酶［19］。CAD与肿瘤相关的报道较少，本研究首次发现，CAD在HCC中高表达，并且与OS相关。CEP55，又名C10orf3，为中心体蛋白相关蛋白家族成员，其主要生物学功能为中心体复制、细胞周期及胞质分裂［20］。自噬蛋白 NBR1通过识别 CEP55，进而招募自噬相关蛋白形成自噬体。LI等［21］发现，CEP55通过调节JAK2-STAT3-MMPs信号通路来促进HCC细胞迁移和侵袭。除PIGU和CAD以外，以上基因均有报道证明直接或间接参与细胞自噬，细胞自噬有助于癌细胞逃避CTL的攻击［6］。本研究进一步证实，在HCC中细胞自噬通路介导CTL相关免疫逃避的关键性，也揭示了细胞自噬在HCC在发生发展中起重要作用。然而，对以上CRG在HCC中作用的了解仍然有限，其影响HCC预后的分子机制以及对HCC患者临床治疗的意义有待进一步研究。

为进一步探索高危组患者预后差的潜在机制，本研究基于TCGA队列的高、低危组进行了GO富集分析、KEGG富集分析及ssGSEA。本研究发现，免疫相关通路、离子通道以及一些癌症通路显著富集。研究证明离子通道在T细胞的稳定与活化发挥重要作用［22］，这与本研究的发现一致。低危组和高危组的抗原提呈相关细胞及功能的富集得分存在显著差异，表明肿瘤相关抗原在高危组更容易识别，但高危组患者预后更差，可能是与高危组的CTL功能下降有关。高危组的Tregs富集得分相对于低危组更高。Tregs显著抑制CTL的功能，进一步证明了高危组的CTL功能被抑制。高危组的Ⅰ型干扰素反应、Ⅱ型干扰素反应及NK细胞富集得分较低危组低。因此，高危组患者抗肿瘤免疫力的减弱是其预后不良的主要因素。

综上所述，本研究构建的模型可有效预测HCC患者预后，有助于HCC患者预后预测体系的完善及临床个体化治疗方案的开发。