原发性肝细胞癌焦亡相关LncRNA预后模型的构建和验证*

张海航,吴 祎,卢彦达

(海南医学院第一附属医院肿瘤内科,海南海口 570100)

原发性肝细胞癌 (hepatocellular carcinoma, HCC) 是世界发病率最高的恶性肿瘤之一,在全球,每年约有 714 850 例新病例和 664 700 例死亡[1]。乙型肝炎病毒或丙型肝炎病毒、大量饮酒、黄曲霉毒素污染的食品、吸烟、肥胖和 2 型糖尿病是 HCC 的重要危险因素[2]。 虽然HCC的诊疗有很多选择,如手术、介入、靶向治疗、放疗和化疗,但由于缺乏早期发现导致癌症患者病死率、复发和转移风险高。 因此,探索预测HCC的准确生物标志物和临床前期模型仍有必要。

肿瘤有多种策略来限制正常细胞保护机制的细胞死亡途径[3]。细胞死亡既是细胞增殖、应激反应和体内平衡的生理调节因子,也是一种肿瘤抑制机制[4]。细胞死亡在某些情况下(包括焦亡、自噬、铁死亡、坏死和凋亡),无疑对我们的健康有益,例如癌症治疗[5-8]。焦亡是一种炎症细胞死亡的新形式。它是由 caspase-1/4/5/11 触发的细胞程序性死亡的炎症形式。它被一些炎症小体激活,释放促炎因子IL-1β和IL -18以诱导细胞焦亡。 caspase-1 的激活将细胞焦亡与其他类型的细胞死亡区分开来。在此过程中,细胞焦亡导致细胞裂解、染色质碎裂并释放细胞内促炎成分[9]。因此,激活诱导细胞焦亡可能会促进HCC细胞死亡,克服现有化疗药物的耐药性,为癌症治疗开辟新的治疗前沿,从而发挥抗癌作用[10]。

长链非编码RNA(long non-coding RNA,lncRNA)是长度不少于200个核苷酸的非编码RNA[11]。 虽然lncRNA不编码蛋白质,但其参与各种生物学调控过程,并与肿瘤的发生、发展和转移有关[12, 13]。目前与焦亡相关lncRNA研究较少。最近的一项研究表明,lncRNA-D16366 在人类 HCC 组织中低表达,因此它可能是该疾病的独立预测因子。一项相关研究表明,SNHG7/miR-34a/SIRT1 轴有助于肝癌期间 NLRP3 依赖性细胞焦亡[14]。敲除RP1 - 85F18.6可以抑制结肠癌增殖和侵袭以及增加结肠癌细胞的焦亡[15]。同样,XLOC_000647 通过降低NRLP3 的表达来抑制胰腺癌细胞的增殖、侵袭[16]。然而,系统地评估焦亡相关的 lncRNA 特征与 HCC 患者总生存期相关性的研究仍然很少。本项研究基于癌症基因组图谱 (The Cancer Genome Atlas, TCGA) 数据构建细胞焦亡相关的 lncRNA 风险模型,并探讨细胞焦亡相关的lncRNA表达在 HCC 预后中的作用。

1 材料与方法

1.1数据收集 从 TCGA 数据库 (https://portal.gdc.cance r.gov/) 获取 mRNA 表达文件和相应的临床病理学数据。 收集了HCC患者临床病理资料,包括年龄(Age)、性别(Gender)、临床分期(Stage)、病理分级(Grade)、T分期、N分期、M分期、生存状态和生存时间。 结果获得了涉及374个HCC组织和50个正常组织的mRNA表达。 去除3个复发样本后,最终获得421个样本,包括肝细胞癌组织样本(n=371)和正常肝细胞组织样本(n=50)的数据。 此外,基于先前的文献,从 TCGA -HCC中提取了 33 个焦亡相关基因表达矩阵[17-20]。

1.2GO和KEGG功能注释 采用R(3.6.3)程序(根据FDP<0.05,|log2FC|≥1)筛选HCC组织和正常肝组织中差异表达的焦亡基因相关的LncRNA。 使用“clusterProfiler” R软件包[21]进行功能注释,GO分析揭示了差异表达基因(differentially expressed genes,DEG)在生物过程、细胞组成和分子功能中的作用,KEGG分析显示了DEG的富集途径。

1.3焦亡相关lncRNA 的鉴定 使用R中的“limma”包[22],根据焦亡基因表达水平进行Pearson相关分析(∣correlation coefficient | >0.4 andP<0.001)鉴定焦亡相关lncRNA。

1.4焦亡相关lncRNA 风险模型的构建和验证 使用R(v3.6.3)中的“glmnet”包中构建焦亡相关lncRNA预后模型[23-24]。风险评分计算公式为:

其中Coefi表示该系数,xi是每个焦亡相关lncRNA的FPKM值。计算所有患者的风险评分。然后通过生存分析、风险图和受试者工作特征(receiver operating characteristic ,ROC)曲线检验风险评分的有效性。使用单因素和多因素Cox回归分析以验证风险评分的独立作用。最后,采用R软件包评估风险评分和临床病理特征、免疫细胞浸润的相关性。

1.5开发列线图来预测HCC患者的预后风险 为了建立一个预测列线图,我们使用多因素Cox回归分析来识别与复发显著相关的危险因素并将其结合起来构建了一个整合预后特征的列线图,来预测HCC患者的1 年、2 年和 3 年生存率。

1.6肿瘤免疫微环境(tumor immune microenvironment,TIME)特征的风险评分 为揭示风险评分与肿瘤浸润免疫细胞之间的相关性,采用TIMER、CIBERSORT、CIBERSORT-ABS、QUANTISEQ、MCPCUNTER、XCELL和EPIC七种方法评估免疫浸润情况,ssGSEA 用于量化两组之间的肿瘤浸润免疫细胞亚群以及评估其免疫功能。 使用Spearman相关性分析探索风险评分与免疫浸润状态之间的相关性。

1.7风险特征在免疫检测点治疗中的作用 基于先前的研究,免疫检测点阻断的靶点的表达有助于免疫治疗的有效进行[25]。 为了进一步探讨不同风险组在免疫治疗中的潜在作用,分析了高低风险组与30个枢纽靶标的表达值之间的相关性。

1.8统计学分析 总生存期(OS)是指从诊断日期到死亡日期的时间间隔。 通过 Kaplan-Meier 对数秩和检验绘制生存曲线。 使用“GSVA”包将ssGSEA量化两组之间的肿瘤浸润免疫细胞亚群进行比较。通过 Cox 回归模型进行单因素和多因素 Cox 回归分析,以验证风险特征的独立预后性能。卡方检验分析了与焦亡相关lncRNA亚组与免疫细胞浸润评分之间的相关性,Spearman分析计算相关系数。 使用ROC曲线评估 1 年、2 年和 3 年OS的预后价值。P<0.05 被认为具有统计学意义。 R 软件(3.6.3 版)用于所有统计分析。

2 结果

2.1GO和KEGG功能注释 从 TCGA-HCC 数据集中提取了33个焦亡相关基因的表达矩阵。基于GO的分析显示过表达的基因主要参与焦亡,BP参与了白细胞介素-1β的产生、焦亡、白细胞介素-1等的产生。 MF主要调节炎症小体复合物、胞质部分、膜筏合成途径的产生。 KEGG富集分析显示,差异基因主要参与NOD样受体信号通路、Toll样受体信号通路、TNF信号通路、IL-17信号通路。见图1。

图1 焦亡相关差异基因的GO和KEGG分析

2.2筛选焦亡相关lncRNA 焦亡基因相关的一个或多个相关(∣correlation coefficient |> 0.4和P<0.001)的lncRNA被定义为焦亡相关lncRNA, 最终筛选出518个焦亡相关lncRNA。见图2。

图2 焦亡基因和焦亡相关lncRNA的共表达网络

2.3焦亡相关lncRNA的风险模型的构建 通过对518 个焦亡相关 lncRNAs进行单因素 Cox 回归分析,确定了 54 个预后相关的焦亡相关 lncRNA,多因素 Cox 构建基于 14 个焦亡相关 lncRNAs风险模型。根据风险评分的公式计算每一位患者的风险值。风险值= MKLN1-AS*0.698+ NRAV*0.190+ LINC00205*0.220+PTOV1-AS1*(-0.347)+SNHG20*(-0.430)+FAM111A-DT*(-0.454)+SNHG4*0.389+ AC074117.1*0.407+AC091057.1*(-1.474)+PXN-AS1*0.190+AL031985.3*0.358+ AC025178.2*0.481+ AC016747.1*(-0.108)+ AC099850.4*0.152。

2.4风险模型独立预后效能评估 根据风险值的中位值,将HCC样本分为高风险组和低风险组。单变量和多变量Cox分析显示:焦亡相关IncRNA 风险模型 (HR: 1.287,95 %CI: 1.220-1.357)和肿瘤分期(HR: 1.680,95 %CI: 1.369-2.062)是 HCC 患者的独立预后因素(表1)。Kaplan-Meier 分析显示低危组的患者预后优于高危组(P<0.001,图 3A)。从风险生存状态图中发现,患者的风险模型与 HCC 患者的生存率成反比(图3B)。 1 年、2 年和 3 年生存率的 AUC值、风险因素的 AUC 值和风险因素DCA在预测HCC预后方面表现出优于传统临床病理学特征的性能(图3C,D,E)。风险模型中14个焦亡相关lncRNAs和焦亡基因的表达关系密切相关(图3F)。结合临床病理学特征和风险评分的列线图是稳定且准确(图 4)。因此,可应用于 HCC 患者的临床管理。

表1 HCC患者OS的独立预后因素分析

图3 焦亡相关 lncRNAs预后效能评估

图4 临床病理因素和焦亡相关lncRNA风险评分的列线图

2.5免疫和基因表达 基于TIMER、CIBERSORT、CIBERSORT-ABS、QUANTISEQ、MCPCUNTER、XCELL和EPIC算法的低危组和高危组细胞成分和细胞免疫反应如图5A所示。免疫细胞亚群之间的相关性分析和基于 TCGA-HCC 数据的 ssGSEA 的相关函数显示,细胞溶解、MHC I 类和Ⅰ型及Ⅱ 型 IFN 反应在低风险和高风险组之间显著不同(图 5B)。 鉴于检查点抑制剂的免疫疗法的重要性,我们进一步探讨了两组之间免疫检查点表达的差异。 我们发现程序性死亡 1(Programmed death 1, PD-1,也称为 PDCD1)、程序性死亡配体 1(Programmed death ligand 1, PD-L1,也称为 CD274)、含黏蛋白结构域的分子 3(Mucin domain‐containing molecule‐3, TIM-3、 也称为HAVCR2),吲哚胺2,3-双加氧酶2(Indoleamine 2,3‐dioxygenase 2, IDO2)等,在两组患者中的表达存在显著差异(图5C)。

图5 免疫细胞浸润和免疫检测点表达

3 讨论

焦亡可影响肿瘤的增殖、侵袭和转移,并受到一些非编码RNA等分子的调控。因此,了解焦亡相关lncRNA调控肿瘤细胞的机制以及焦亡在癌症中的潜在作用有助于预防和治疗癌症。我们首先从TCGA数据集中提取了33个与焦亡相关的基因,基于GO的分析显示,主要参与炎症小体复合物的焦亡和调控诱导。KEGG分析进一步揭示,这些基因主要参与NOD样受体信号通路、Toll样受体信号通路、TNF信号通路、IL-17信号通路。最近的一项研究发现,来自人诱导多能间充质干细胞的外泌体可能靶向抑制NOD样受体蛋白3(NOD-like receptor protein 3,NLRP3)炎症小体通路,以抑制LPS/ATP诱导的肺泡巨噬细胞(AM)[26]。据报道,Toll样受体(Toll-like receptors, TLRs)识别危险相关信号并激活细胞中的核因子-κB(nuclear factor-κB, NF-κB)信号通路,诱导炎症小体的激活并发挥焦亡的作用[27]。重组人白细胞介素-1受体拮抗剂(RhIL-1Ra)可能通过减少肿瘤坏死因子-α(Tumor necrosis factor-α, TNF-α)和白细胞介素-17(Interleukin-17, IL-17) 的分泌以抑制伴刀豆球蛋白A参与肝细胞焦亡[28]。

最近的一项研究发现,HCC组织和细胞系中MKLN1-AS的上调与HCC 患者的总生存期较短有关。在功能上,敲除MKLN1-AS会破坏HCC细胞的增殖、迁移和侵袭,并诱导细胞凋亡[29]。lncRNA LINC00205的过表达增加了HCC细胞的致瘤性,这是由miR-184和EPHX1诱导的细胞凋亡引起的[30]。 NHG20的表达与肿瘤的大小和临床分期相关。在SK-Hep-1细胞中敲除SNHG20可显著抑制细胞增殖、迁移和侵袭[31]。转录组学分析显示,MBNL3诱导lncRNA PXN-AS1第4外显子的纳入, PXN-AS1可作为HCC患者诊疗的重要靶点[32,33]。

将焦亡相关的 lncRNA 的不同表达分为高风险和低风险两类,以探索它们在 HCC 中的潜在作用。免疫检查点抑制剂 (Immune checkpoint inhibitors, ICIs) 可与细胞焦亡结合增强抗肿瘤活性[34]。尽管高 ICI 评分亚组具有较高的免疫评分,但在低 ICI 评分亚组中细胞溶解、I 型 IFN 反应和Ⅱ型 IFN 反应的分布明显更高。然而,具有高 ICI 评分的 HCC 患者显示出更高的 MHC I 类表达水平。先前的研究表明,肿瘤浸润性 T 细胞和 B 细胞的紧密接近表明它们之间的功能相互作用与增强的局部免疫激活有关并有助于改善 HCC 患者的预后[35]。此外,较高的风险评分与预后受损相关,包括Ⅱ型 IFN 反应和 I 型 IFN 反应的活性以及溶细胞的分数。因此,高危患者的焦亡受损和抗肿瘤免疫减弱可能是其预后不良的原因。

分析高低风险组中焦亡相关 lncRNA 的功能。 KEGG通路分析显示,高危组差异化焦亡相关lncRNA主要富集于焦亡相关信号通路,如转录因子、NOD样受体信号通路、mTOR信号通路和Toll样受体信号通路。 因此,本研究认为与焦亡相关的lncRNA可能在焦亡过程中起关键作用。总的来说,焦亡相关lncRNA在肝癌预后中的作用相对缺乏,本研究开发的新的预测预后的模型可能为未来的癌症控制提供有价值的建议。