胰腺癌是一种严重威胁患者生命的高度恶性肿瘤,是癌症相关死亡的第七大原因
[1]。由于胰腺属于腹膜后位器官,且位于腹腔深部,背靠脊柱,故胰腺癌患者的典型症状往往出现较晚。90%的胰腺癌在发生远处转移后才被患者察觉,已处于晚期,预后普遍较差,其5年总生存(OS)率不足10%
[1-2]。胰腺癌具有致密的结缔组织增生性基质,由癌症相关成纤维细胞产生,肿瘤组织伴有明显的纤维化,密度增加。肿瘤细胞的快速生长,其对氧气和养分需求较大,与之相应的血管供给却往往不足以满足需求,从而导致微环境缺氧状态
[3]。
氧化应激在肿瘤的发生发展中起着至关重要的作用,是细胞对活性氧(ROS)产生的一种自我保护机制。当ROS生成和抗氧化防御之间的平衡被破坏时,就会发生氧化应激
[4]。正常细胞在有氧条件下主要依靠有氧氧化产生能量,而肿瘤细胞即使在氧含量充足的情况下,仍然倾向于通过糖酵解途径产生能量,这种现象称为“Warburg效应”,这是肿瘤细胞生长、存活和能量供应的关键特征
[5]。此外,研究
[6]发现,ROS积累可进一步刺激有氧和无氧糖酵解,增加糖酵解产物乳酸,从而使细胞更加依赖能量供应,加速肿瘤进展。升高的ROS水平可激活信号通路,如缺氧诱导因子1(HIF-1)通路,进而影响乳酸代谢,使细胞进入催化厌氧糖酵解状态
[7]。Morusin通过ROS介导抑制FOXM1/c-Myc信号通路,在前列腺癌细胞中诱导凋亡和Warburg效应
[8]。乳酸代谢也可以影响氧化应激。在神经母细胞瘤细胞中,L-乳酸通过促进轻微的ROS爆发,促进细胞防御机制,在体外和体内对抗氧化应激具有保护作用
[9]。乳酸通过激活自噬机制,保护视网膜免受氧化应激引起的退化
[10]。乳酸运输稳态的破坏扰乱了心肌细胞内乳酸-丙酮酸的平衡,导致氧化应激和炎症反应,加剧了心肌损伤
[11]。因此,氧化应激和乳酸相互作用和相互调控,二者合并起来分析将产生更为完整的分析视角,从而筛选出关键的生物标志物。
本研究旨在用生物信息学的分析方法,采用公共数据库的测序数据及临床信息,对氧化应激联合乳酸代谢在胰腺癌中潜在作用进行了综合分析。筛选并得到胰腺癌预后相关的生物标志物,以揭示潜在的致病机制,为胰腺癌的个体化诊疗提供新的理论依据和选择。
1 资料与方法
1.1 数据收集和预处理
来自肿瘤基因组图谱(The Cancer Genome Atlas,TCGA)数据库(
https://portal.gdc.cancer.gov)的胰腺癌(PAAD)研究项目提供了预处理的RNA测序(RNA-seq)信息,包括180个肿瘤和3个正常样本。剔除了缺少样本表达及生存信息的6个样本后,获得了177例患者的临床数据。
为了测试预后特征并避免数据偏差,在基因表达全景图谱(Gene Expression Omnibus,GEO)数据库(
https://www.ncbi.nlm.nih.gov/geo)选择了数据集GSE57495,包含GSE57495的180个胰腺癌肿瘤样本。数据合并及去除批次效应的操作依赖于“inSilicoMerging”R包。来自TCGA的RNA-seq数据以及GSE57495的序列矩阵使用“log2(标准化基因表达+1)”方法处理。表达量含有“NA”的样本被去除,缺失值使用“impute”R包进行补充。
所有使用的数据均可公开获取,严格遵守TCGA和GEO数据库发布的政策和指南。在MesigDB上检索以“乳酸代谢,氧化应激”为检索词,得到其相关的基因。
1.2 预后模型构建
本研究首先使用“limma”R包进行差异表达基因的分析。特别利用了包中的“Wilcox”功能来执行非参数的Wilcoxon秩和测试。通过这种方法,能够识别出在胰腺癌样本与对照样本之间表达差异显著的基因。笔者设定了严格的筛选标准,包括调整后的P<0.05,以及至少2倍的表达变化,从而确保筛选出的差异表达基因与乳酸代谢和氧化应激高度相关。通过“Wilcox”测试,筛选出显著差异表达的基因,这些基因在乳酸代谢和氧化应激途径中表现出显著的上调或下调。
基于训练队列的OS,对筛选得到的差异表达基因进行了单变量Cox回归分析,来评估胰腺癌中乳酸代谢与氧化应激相关基因的预后预测价值。每个基因的表达水平被视为独立变量,计算其与患者OS的相关性和风险比(HR),识别出具有统计学意义的预后标志物,即那些HR显著不等于1的基因。
为了进一步细化候选基因列表,使用Venn图将单变量Cox分析中得到的预后相关基因和差异表达基因进行交叉比较,重叠基因被认为是与胰腺癌预后密切相关的差异表达基因。预后特征的训练数据集来自TCGA训练集中的肿瘤组织。随后采用了最小绝对收缩和选择算子(LASSO)回归方法鉴定最具预测力的基因。该过程确定了每个基因的最优惩罚参数λ,以建立风险评分计算系数,通过交叉验证确定了最优的λ值,这一过程提高了模型的预测精度。最后,根据LASSO回归的结果,利用变量的线性组合构建了一个线性风险评分模型。模型是采用如下公式构建:。其中“N”为基因数,“exp”为基因表达量,“coef”是通过LASSO回归得到的该基因的对应系数。将患者以风险评分得分情况分为高风险组和低风险组,并使用Kaplan-Meier生存曲线对比组间生存差异情况。这种方法不仅可以揭示不同风险组患者的生存预期差异,还能通过Log-rank检验来评估这种差异的统计学意义。此外,使用“timeROC”R包(4.21版)评估模型的预测性能,特别是通过测量受试者工作特征(ROC)曲线下面积(AUC),并以此评估模型对患者5年OS率和其他指标的预测准确率。
1.3 预后模型的验证
从GEO数据库中获取GSE57495数据集,验证上述建立的胰腺癌预后风险模型。通过建立的模型的基因系数算出胰腺癌患者的风险评分,基于此划分高风险患者组和低风险患者组,同时利用Kaplan-Meier曲线图来比较组间的预后生存异同。
1.4 列线图构建
使用风险评分预测胰腺癌患者的OS率,使用列线图结合临床信息预测OS率。校准曲线可以评估预测的OS率和观察到的OS率之间的一致性。为了明显基因的准确度,采用了一致性指数(C指数)曲线。
1.5 风险组间不同的肿瘤免疫微环境模式
CIBERSORT利用一系列参考转录组信息来估计大量肿瘤样本表达数据中发现的22种不同免疫细胞的比例。这种估计是利用线性支持向量回归原理进行的。为了评估与免疫相关的功能,借助“GSVA”R包进行了ssGSEA分析。ESEIMATE算法用于分析免疫细胞和基质细胞在肿瘤微环境中的浸润情况。
2 结 果
2.1 胰腺癌中氧化应激与乳酸代谢差异基因识别
为鉴定并得到胰腺癌样本对比正常组织中存在差异表达基因,本研究使用limma差异分析算法得到共有49个基因在180个肿瘤组织和3个正常组织中差异表达,且具有统计学意义(
P<0.001,且log|FC|>2)(
图1A-B)。对
n个乳酸代谢及氧化应激相关基因进行单因素Cox回归,并根据TCGA队列的OS选择99个预后相关基因(
图1C)。将差异基因差异表达基因与预后基因使用Venn图相交,最终得到7个差异-生存相关基因(
图1D)。计算7个差异-生存相关基因的Pearson相关系数,以0.2为相关系数来评估其潜在相互作用(
图1E),结果显示,
MUCI与
MUC5AC、
KRT18、
SDC1呈正相关,与
ATPAF2呈负相关;
KRT18与
AREG、
MUC1呈正相关,与
DDC呈负相关;
ATPAF2与
DDC呈正相关。
2.2 预后模型构建
采用LASSO算法,结合RNA测序数据及患者预后生存信息,筛选得到6个基因:
MUC1、
KRT18、
SDC1、
ATPAF2、
DDC、
AREG,并组成联合乳酸代谢与氧化应激为视角的胰腺癌预后预测模型。6个基因的风险评分计算方法:风险评分=(0.054 713 97)×MUC1表达值+(0.208 527 45)×KRT18表达值+(0.035 224 30)×SDC1表达值+(0.109 608 38)×AREG表达值+(-0.029 191 22)×DDC表达值+(-0.543 015 84)×ATPAF2表达值(
图2)。通过建立的模型计算TCGA中胰腺癌患者的风险评分,并分为高风险(
n=251)和低风险(
n=252)两组(
图3A),结果显示高风险组患者处于死亡状态的分布更高(
图3B)。主成分分析(principal component analysis,PCA)表明预后特征能有效地将患者分为高风险和低风险两组(
图3C)。Kaplan-Meier曲线,高风险组患者的OS率明显降低(
P<0.01)(
图3D)。ROC分析结果显示,此模型对胰腺癌患者1、2、3年的OS率的预后预测具有较高的准确率,AUC分别为0.710、0.674、0.649(
图3E)。
图3F展示了高低风险组之间6个模型相关基因的表达情况及差异。
2.3 GEO样本的模型效能验证
GSE57495胰腺癌患者样本经过中位评分被分为高风险组和低风险组(
图4A)。结果显示高风险组患者处于死亡状态的分布更高(
图4B)。PCA表明,预后特征能有效地将患者分为高风险和低风险两组(
图4C)。AUC在1、2、3年分别达到0.600、0.546和0.618(
图4D)。验证队列的Kaplan-Meier生存曲线显示,低风险组样本的OS较好,与训练队列相同(
图4E)。风险热图分析显示,高风险组中的患者6个建模基因的表达更为显著(
图4F),这与之前得到的结果基本一致。
2.4 列线图构建及评估
为形成可靠的临床决策,促进个体化诊疗,纳入风险评分、性别、分期、分级、年龄等多种临床参数构建列线图模型。校准曲线显示了1、3、5年的模型预测的生存结果与实际生存结果之间的一致性,临床分期与年龄以及风险评分被视为有意义的临床预测指标。结合三者可以对胰腺癌患者预后情况进行更为精准的预测(
图5A)。校准曲线证明本列线图的预测效果与实际情况接近(
图5B)。
2.5 临床相关性分析
分析结果显示,风险评分年龄(以65岁为分界点)、性别与临床分期之间差异无统计学意义(均
P>0.05)。由此可见,不论年龄、性别或者临床分期结合本模型相关的风险评分对于胰腺癌患者的预测无意义,但病理分级的1~3级之间,两两差异存在统计学意义(均
P<0.05)(
图6)。这或许提示本风险评分对于早期胰腺癌患者的预后预测更为精准,而晚期胰腺癌患者的预后则需要纳入更多临床及试验数据进行分析,以提升本模型的预测效能。
2.6 CIBERSORT评估免疫细胞组成情况
为进一步阐明风险评分与免疫反应之间的内在联系,使用CIBERSORT工具来估计胰腺癌样本中22种免疫细胞类型的丰度(
图7A)。评估胰腺癌中高风险样本和低风险样本之间免疫细胞浸润的差异,结果显示,效应T细胞、CD8
+ T细胞、M0巨噬细胞和中性粒细胞在不同风险组之间存在明显差异(
图7B)。在这些免疫细胞差异分析中,以M0巨噬细胞在高风险组样本中的浸润高于低风险组,这可能提示高风险样本中较差的生存表现与M0巨噬细胞的浸润情况存在相关性。六种免疫状态中,包括C1(伤口愈合),C2(IFN-γ显性),C3(炎症),C4(淋巴细胞耗尽),C5(免疫安静)和C6(TGF-β显性)。C1和C3之间,C2和C3之间存在统计学意义(
图7C),这提示胰腺癌患者在伤口愈合,IFN-γ显性,炎症方面存在相关性。
3 讨 论
氧化应激和乳酸代谢分别都在肿瘤微环境中发挥作用。作为预后较差的肿瘤类型,胰腺癌的研究价值和重要性毋庸置疑。本研究利用生物信息学中的LASSO/Cox算法分析了其在胰腺癌中的潜在作用,并筛选出了与预后相关的生物标志物。认识到这些基因在胰腺癌生物学中的作用,并评估该模型对胰腺癌患者预后的预测性能,关键基因在胰腺癌个体化治疗策略等其他潜在应用中的价值。
胰腺癌是一种高度恶性肿瘤,5年OS率较低,需要精确的预后模型来评估患者状况、指导治疗并改善治疗效果,有助于选择最佳治疗策略,也可以预测患者对治疗的反应和生存期,这将帮助患者获得更有效的治疗并提高生存质量。借助对预后模型的研究,可以更深入地了解到胰腺癌的发病机制、高危因素等,从而开发出更精确、靶向性更强的治疗策略。根据构建的预后模型,TCGA胰腺癌样本和GEO验证队列被分为高风险和低风险两组,这证明本模型具有良好的判别能力。该模型对患者预后的预测能力得到多种模型功效测试算法、生物信息学分析和多角度验证的肯定。结果证实,该模型可以独立、稳定地预测胰腺癌患者的预后,并具有良好的预测效果,还可以被视为一个独立的预后因素。在此模型的基础上,建立了一个包含风险评分、性别、分期、分级和年龄等因素的列线图。风险评分是预测胰腺癌患者预后的独立危险因素。通过测量胰腺癌患者中与该预后特征相关的6个核心基因的表达水平,可以对其预后作出一定的预测。
氧化应激与乳酸代谢是相互作用和相互调控的
[12-14],所涉及的关键基因很多。本研究中涉及的关键基因对胰腺癌预后的影响已被大量报道。MUC1,即黏附蛋白1,是一种高度糖基化的跨膜黏蛋白,位于上皮细胞腔面,为黏蛋白家族中最易被识别的跨膜蛋白成员之一
[15]。MUC1蛋白具有重要的生物学功能,主要是通过参与不同信号通路调节肿瘤的各个方面,如细胞生长、增殖、转移、凋亡和发育过程等
[16],并与其他几种致癌受体和途径相互作用,如HER2、β-catenin、NF-κB和雌激素受体(ERα)
[17]。MUC1蛋白可能对细胞间的黏附、生长因子受体及其最终的信号通路发挥多方面作用,增强肿瘤细胞对环境信号的响应能力。针对MUC1亚基细胞外区域的单克隆抗体MUC1-C(即抗hMUC1抗体)在实验中被证明能够对胰腺癌的发展产生显著影响
[16],具体表现为触发胰腺癌细胞中EGFR的内化,导致EGF诱导的ERK磷酸化受到抑制,从而抑制胰腺癌细胞的增殖和迁移。角蛋白18(KRT18)是一种特定类型的角蛋白,在上皮细胞中表达的中间纤维蛋白家族成员,KRT18与角蛋白8(KRT8)紧密结合形成角蛋白细胞骨架,起着支持和保护细胞结构的作用
[18]。已被证实KRT18可以作为胰腺癌肿瘤进展的潜在标志物
[19]。在胃癌组织中,KRT18表达上调,研究发现其过表达与不良预后密切相关,可能通过调控细胞周期、增殖和凋亡途径促进胃癌的恶性生物学行为
[20]。SDC1(也称为CD138、免疫球蛋白超家族的1因子)是一种跨膜肿瘤相关糖蛋白,广泛表达在多种肿瘤类型中
[21]。SDC1在肿瘤细胞中的主要作用是调节细胞核心蛋白聚糖的表达和分布以及细胞外基质的构成,这些糖蛋白质与细胞外基质相关的生物分子(如生长因子、趋化因子和蛋白酶等)相互作用,从而影响肿瘤细胞间通讯和信号传导,进一步调节肿瘤生长、迁移和免疫逃逸的能力
[22]。过度表达也能促进胰腺癌和乳腺癌中癌细胞的生长和增殖
[23]。根据Chen等
[21]的研究,SDC1可以通过与CCL5相互作用,进而促进肿瘤细胞的迁移能力。两性调节蛋白(amphiregulin,简称AREG)属于表皮生长因子(EGF)家族成员,可与EGF受体(EGFR)结合,激活EGFR介导的信号传导通路
[24]。在先前的一项研究中发现,在胰腺癌中,AREG表达水平明显上升,且与不良的临床预后紧密相关;AREG通过上调和激活多种信号通路(如EGFR、RAS、RAF、MEK、MAPK、PI3K/Akt等)促进胰腺癌细胞的生长、侵犯和转移。AREG的抑制或下调可显著减慢胰腺癌细胞的增殖、迁移和侵袭能力
[25-26]。在卵巢癌中,AREG可以通过调控AREG-EGFR-ERK信号通路,促进癌细胞的生长和存活,增强细胞对化疗药物的抵抗力;研究发现,抑制或干扰AREG的表达可以降低卵巢癌细胞的耐药性,从而提高化疗的疗效
[27]。多巴胺脱羧酶(dopa decarboxylase,简称DDC)是一种关键的生物合成酶,在神经递质合成过程中具有重要作用,DDC能将多巴胺前体氨基酸(3,4-二羟基苯基固氨酸,L-DOPA)转化为多巴胺,同时也负责将5-羟色胺前体氨基酸(色氨酸)转换为5-羟色胺(血清素)
[28-29]。在胰腺癌细胞中,DDC的过表达可显著增强细胞生长、侵袭和转移能力;当使用小干扰RNA(siRNA)或其他方法抑制DDC表达后,胰腺癌细胞的肿瘤生成能力显著降低
[30]。这证实了DDC在恶性胰腺癌的生长、侵犯和转移中发挥了关键作用。这些发现提示DDC在胰腺癌的生长、侵犯和转移过程中发挥着重要作用。
ATPAF2基因被认为在细胞内负责催化ATP的合成和降解
[31]。尽管目前没有直接证据将ATPAF2与胰腺癌关联,但若考虑到ATPAF2在线粒体呼吸链复合体I中的职能,它可能在肿瘤生长、侵犯和转移过程中发挥关键作用,在肿瘤环境中,氧化应激水平升高,可能影响ATPAF2的表达水平。因此,研究ATPAF2在不同肿瘤类型中的表达差异可能有助于评估肿瘤发育、恶变和预后。
尽管关于胰腺癌预后模型的研究有价值,但本研究仍然有一定的局限性,在研究中发现与胰腺癌预后有关的基因,需通过细胞、动物实验和临床试验进一步验证研究发现的生物学功能和意义;另外,由于当前的技术和成本效益问题,基于测序的转录组数据构建的预后预测模型无法在临床普遍应用。针对这些局限性,后续研究将通过加强深度学习、细胞与动物实验等方式,持续优化胰腺癌预后模型,提高其预测准确率和实用性。有理由相信,随着技术的不断进步和成本的逐渐降低,未来这些基于高通量测序的预后预测模型将逐渐走向临床应用。
本研究分析了来自在线开放数据库的数据,首次尝试用系统的方法评估和总结与胰腺癌相关基因的预后意义。本研究建立的预后模型可以预测胰腺癌患者的生存和预后,可能是胰腺癌患者潜在的免疫治疗靶点。