基于铁死亡相关基因的生物信息学分析构建结直肠癌预后模型

陈梁玥 段怡平 柳家翠 马甜甜 朱翠雯 李东旭 刘永豪 李慧 喻明霞★

2019年癌症统计数据显示消化系统肿瘤新增病例数和死亡人数均位居世界第一，而结直肠癌（Colorectal Cancer，CRC）两者均位居消化系统榜首［1］。根据我国恶性肿瘤死亡率排名，在男性中居第5 位，女性中居第4 位［2］。CRC 起病隐匿，早期症状不明显且不典型，缺乏特异性，有时与其他消化系统疾病共存，因此其发现时多为中晚期，且部分患者出现转移灶［3］。尽管目前对CRC 的治疗方式有手术切除、放化疗、基因治疗等多种手段，但其五年生存率仍不足50%［4］。因此，通过CRC 预后相关的生物标志物对CRC 患者的预后进行早期评估，成为了必然。

铁死亡（Ferroptosis）是一种非凋亡形式的细胞死亡，其特征是铁依赖性脂质氢过氧化物的积累达到细胞致死水平［5］。自铁死亡机制被提出以后，人们在多种疾病发生与发展的关键过程中证实了铁死亡的作用，如阿尔兹海默病、急性肾损伤、肿瘤的增殖与凋亡等［6］。Zhang［7］等发现，苯并吡喃衍生物能通过调节AMPK/mTOR/p70S6k 信号通路的活性下调SLC7A11的表达，降低半胱氨酸和谷胱甘肽的含量来诱导结肠癌细胞的铁死亡。而Malfa［8］等证实桦木树提取物通过促进细胞微环境氧化，导致血红素加氧酶⁃1 表达水平增高提高细胞内游离铁浓度，从而诱导铁死亡。此外，XIAN［9］等研究表明circABCB10可以通过调控miR⁃326/CCL5 轴促进直肠癌细胞铁死亡，从而成为潜在的治疗靶点。

本研究通过分析TCGA 数据库中CRC 的高通量测序数据，通过生物信息学方法筛选与CRC 预后相关的铁死亡基因，并加以分析构建临床预后模型，以期为CRC 的诊疗及预后提供评估参考。

1 材料与方法

1.1 数据采集

从TCGA（http：//gdc⁃portal.nci.nih.gov/）数据库下载CRC 患者肿瘤组织癌旁组织的mRNA 高通量测序数据和临床资料，其中包括568 个肿瘤组织和44 个正常组织。对基因表达数据进行归一化处理，并Log2 转化后用于后续分析。排除了临床信息缺失或患者生存时间为0 的样本，最终纳入了538 名CRC 患者。

使用Strawberry Perl 5.32.1 软件将原始基因表达谱数据整理成基因表达矩阵文件。

1.2 方法

1.2.1 铁死亡相关基因的提取和处理

基于美国国家生物信息中心数据库挖掘调控肿瘤发生发展的60 个铁死亡相关基因（Iron death⁃related genes，IRGs），用于后续对照分析。通过R语言“limma”包对肿瘤和癌旁组织进行差异分析，筛选出铁死亡相关的差异基因（Differentially ex⁃pressed gene，DEGs）。筛选阈值为伪发现率（false discovery rate，FDR）<0.05。运用R 语言“Venn 包”对具有差异和价值的IRGs 取交集，并绘制热图。

1.2.2 预后模型的构建及评价

为了达到最大化减小过度拟合风险的目的，采用LASSO Cox 回归分析并结合R 语言“glmnet”包再次筛选基因。回归分析中的独立变量为交集基因的表达量，而响应变量是TCGA 队列中患者的总体生存率和生存状态。根据每个基因的标准化表达水平及其相应的回归系数计算患者的风险评分。风险评分表达式如下：风险评分（Risk Score），其中N 是基因数，Expi 是基因的表达水平，Ci 是LASSO 回归分析得到的回归相关系数。根据预后模型得到的风险评分求其中位数，将纳入的CRC 患者划分为高、低风险组。对两组采用Kaplan⁃Meier 方法绘制生存曲线，同时使用R 语言“time ROC”包来评估预后模型的预测能力。根据下载整理的临床病理信息，采用单因素和多因素COX 回归对预后模型进行分析，以验证构建的预后模型是否能作为独立的预测因子。

1.2.3 功能富集分析和免疫功能状态分析

通过R 语言“cluster profiler”和“limma”包对高、低风险组中的DEGs 进行基因本体论（GO）和京都基因和基因组百科全书（KEGG）分析。此外，基于R 语言“gsva”包，通过单样本基因集富集分析（ssGSEA）对16 个免疫细胞浸润进行打分，并分析13 个免疫相关的通路。

1.3 统计学分析

应用R（4.0.2）和SPSS Statistic25.0 软件进行统计学分析；使用Wilcoxon 非参数秩和检验分析肿瘤样本和癌旁样本中的基因表达差异；采用χ2检验对计数资料之间的差异进行分析，以P<0.05 为差异有统计学意义。

2 结果

2.1 确定与铁死亡相关的差异预后基因

大部分IRGs 在肿瘤组织和临近组织间存在差异表达（83.3%，50/60），见图1A。单因素COX 分析显示，有6 个基因与预后相关，并绘制了基因表达热图（图1B⁃C），同时进行了LASSO 回归分析，见表1。

表1 6 个铁死亡相关基因的详细信息Table 1 Information of the 6 genes associated with iron⁃death

2.2 构建预后模型和评价

PCA 和t⁃SNE 分析表明，高风险和低风险两组患者分布在两个方向，见图2A⁃B。Kaplan⁃Meier 生存曲线提示低风险组总体生存率显著高于高风险组，差异有统计学意义（P<0.05），见图2C。并绘制受试者工作曲线（AUC1=0.594，AUC2=0.708，AUC3=0.599，图2D）。同时，对患者进行了风险评分（图3A⁃B）。单因素和多因素COX 分析提示，该预后模型可以作为一个独立的预后因子（P<0.05）（图3C⁃D）。

图2 模型的预测效能检验Figure 2 Predictive power test of the model

图3 预后模型的生存状态及COX 回归分析Figure 3 Survival status and COX regression analysis of the prognostic model

2.3 功能分析

GO 分析中，DEGs 富集在金属离子的隔离以及离子门控通道活性等几个铁相关的通路，此外，DEGs 也明显富集在与免疫相关的通路中，如补体的激活等。KEGG 分析中，在细胞因子与细胞因子受体相互作用等与免疫反应相关的功能通路富集。见表2。

表2 GO 富集分析和KEGG 通路分析Table 2 GO enrichment analysis and KEGG pathway analysis

ssGSEA 分析表明，巨噬细胞、调节性T 细胞等在两组之间比较，差异有统计学意义（P<0.05）。两组之间的免疫功能比较，除Ⅰ型主要组织相容性复合物（MHC）外，差异有统计学意义（P<0.05），见图4。

图4 不同风险组之间ssGSEA 评分Figure 4 ssGSEA scores between different risk groups

3 讨论

在本研究中，通过文献查阅和数据库挖掘，系统的分析了60 个IRGs 在CRC 组织中的表达与其总体生存率之间的关系，首次构建了一个关于6个IRGs 的预后模型。该模型将538 例CRC 患者分为高、低风险组，两组患者DEGs 在铁相关和免疫相关的通路中被富集。一些基因可能通过调节某种激活剂或抑制剂进而调控肿瘤的铁死亡［10］，但这些基因的表达与CRC 患者的生存预后之间的相关性仍不清楚。但值得注意的是，大约83.3%的IRGs 在CRC 组织和癌旁组织之间存在有差异表达，单因素COX 分析显示有6 个基因与总体生存预后相关。这表明铁死亡在CRC 中具有潜在的作用，因此，利用IRGs 建立CRC 预后模型是有意义的。

本研究构建的预后模型由6 个IRGs 组成，包括AKR1C1、ATP5MC3、HMGCR、CRYAB、FDFT1、ALOX12。人Aldo⁃keto 还原酶1 家庭成员C1（Human Aldo⁃keto reductase 1 family member C1，AKR1C1）的过度表达可能导致肿瘤细胞的耐药［11］，并通过催化醛酮转化为相应的醇，解毒活性脂质过氧化物分子，进而影响铁死亡进程［12］。ATP5MC3又称ATP5G3，编码线粒体ATP 合成酶亚单位，参与肿瘤的多种生物学过程，并抑制了erastin 诱导的铁死亡［13］。3⁃羟基⁃甲基戊二醇⁃CoA还原酶（3⁃hydroxy⁃methylpentanediol⁃CoA reduc⁃tase，HMGCR）在甲戊酸的合成中发挥重要作用，Shimada 等人发现抑制HMGCR可增强FIN⁃56 诱导的铁死亡［14］。CRYAB是小型热休克蛋白家族的一员，参与调节肿瘤中的几种信号通路，有报道称其可以调节铁的摄入量和GPX4 的丰度［15］。角鲨烯合成酶（Squalene synthase，FDFT1）是合成甾醇并最终合成胆固醇的关键酶，该基因的敲低与CRC 的恶性进展和不良预后相关，并通过负调控AKT/mTOR/HIF1α 通路抑制肿瘤发生［16］。ALOX12是脂氧合酶超家族成员之一，有氧条时催化多不饱和脂肪酸底物形成氢基二十碳四烯酸，有研究证明，ALOX12的失活减少了p53介导的铁死亡，并消除了p53依赖的肿瘤生长抑制［17］。综上，这些基因在一定程度上影响CRC 的进展，但其中的大部分基因是否通过铁死亡机制来干扰CRC 患者的预后还尚不明确。

近年来铁死亡在肿瘤领域的研究中较为热门，但肿瘤免疫与铁死亡之间的潜在调控机制尚不清楚。对高、低风险间的DEGs 进行GO 和KEGG 分析，发现这些DEGs 富集于一些与免疫相关的生物学通路和功能中。由此可推断，铁死亡与肿瘤的免疫功能之间可能存在密切关联。此外，在ssGESA 分析中，发现抗原提呈过程在高、低风险组之间存在显著差异，一种可能性是铁死亡细胞释放某些信号，吸引抗原提呈细胞，从而激活抗原提呈反应［18］。同时，本研究通过对比分析，高风险组中有较高的肿瘤相关巨噬细胞、调节性T 细胞表达，既往研究已经证实该类细胞水平的升高与CRC 肿瘤的不良预后相关［19］。

综上所述，本文经过多种生物学方法，构建了6 个与铁死亡相关基因的预后模型，可能为CRC患者的个体化治疗和评估提供参考。