基于失巢凋亡相关基因的前列腺癌预后风险模型的构建与分析

叶煜冰,杨泽锐,3,陈楚君,仇志坤,李 雄,2,3

(1.广东药科大学附属第一医院临床药学重点专科;2.广东药科大学基础医学院;3.国家药监局药物警戒技术研究与评价重点实验室,广东 广州 510000)

前列腺癌(prostate cancer,PCa)是男性中发病率较高的癌症之一[1]。在接受手术治疗或放化疗后,仍然有近三分之一的患者最终会出现生化复发(Biochemical recurrence,BCR)或临床复发,预示着潜在的临床转移和预后不良。能够早期发现BCR或临床进展是降低PCa死亡率的关键方法[2]。血清前列腺特异性抗原(prostate specific antigen,PSA)是目前公认的用来筛查和诊断PCa的生物标志物,但由于其特异性差,它也使得低风险PCa进行不必要的活检或严重过度治疗[3]。因此,迫切需要为PCa患者寻找新的早期预后标志物。

失巢凋亡是由于细胞对细胞外基质(ExtraCellular Matrix,ECM)的黏附丧失或细胞黏附不当而诱导出现的一种特殊的细胞死亡方式[4]。肿瘤细胞在丧失对ECM的粘附性后通过旁自分泌以及旁分泌机制抵抗凋亡继续存活,并重新获得附着能力促进肿瘤扩散、侵袭和转移[5]。研究证实,失巢凋亡相关基因(anoikis-related genes,ARGs)在包括前列腺癌中的许多肿瘤的进展中发挥作用,包括胃癌[6],头颈鳞癌[7],胶质瘤[8],乳腺癌[9]等。据报道,前列腺癌中跨膜丝氨酸蛋白酶4(transmembrane serine protease 4,TMPRSS4)能够通过上调SOX2(SRY-Box Transcription Factor 2,SOX2)促进失巢凋亡抵抗,提高了循环肿瘤细胞的存活率,促进肿瘤早期转移[10]。

尽管已证实失巢凋亡与肿瘤的预后相关,但基于ARGs的前列腺癌预后指标暂未得到分析。因此,本研究主要基于癌症基因组数据集(The Cancer Genome Atlas Program,TCGA),构建一个基于3个ARGs构成的预后风险模型并验证其预测能力,为预测前列腺癌患者的预后情况和寻找新的潜在治疗靶点提供新思路。

1 资料与方法

1.1 数据处理通过将TCGA-UCEC(https://xena.ucsc.edu/)的数据整合到癌症基因图谱(TCGA)(https://portal.gdc.cancer.gov),获得了419名PCa患者的基因表达谱和临床数据,包括BCR状态、BCR生存时间、年龄、T分期、N分期、PSA值和Gleason评分。从GeneCards(https://www.genecards.org/)中获得778个ARGs,选择相关性评分>0.4的基因。通过在基于基因表达水平值的交互式分析平台(Gene Expression Profiling Interactive Analysis,GEPIA)中以“PRAD”为关键词,染色体分布选择“both”,其他为默认选项,获得了PCa中差异表达基因(Differentially expressed gene,DEGs),建立火山图。最后,通过与DEGs数据集取交集获得差异表达的ARGs。

1.2 GO和KEGG通路富集分析为进一步了解与前列腺癌相关的失巢凋亡基因所涉及的生物学功能和通路,使用R包“clusterProfler”进行基因本体论(Gene Ontology,GO)和京都基因组百科全书(Kyoto Encyclopedia of Genes and Genomes,KEGG)通路富集分析。设定分析的物种为“智人”,显著性水平P值<0.05。

1.4 建立基于ARGs的列线图通过将419名患者的临床信息与风险评分相结合,并使用R语言中的“survival”等包进行独立预后分析。此外,基于风险评分和其他临床病理学特征创建了基于ARGs的列线图以预测病例的临床结果,绘制校准曲线,通过折线图的预测能力与观察到的生存结果进行比较来估计ARGs的可靠性。

1.5 免疫浸润分析基于TCGA训练集数据,使用CIBERSORT(Cell-type Identification By Estimating Relative Subsets Of RNA Transcripts,CIBERSORT)工具量化了22种类型的免疫细胞比例。将P<0.05设置为阈值。此外,采用单样本基因集富集分析(single sample gene set enrichment analysis,ssGSEA)方法反映前列腺癌中免疫微环境(tumor microenvironment,TME)的状态。

1.6 统计分析所有数据分类和分析由R软件(4.2.1)完成。应用Kaplan-Meier(K-M)分析评估两组之间的生存差异。使用受试者工作特征曲线(receiver operating characteristic curve,ROC)分析测试了ARGs模型的可靠性。所有分析的统计学显著性均设置为P<0.05。

2 结果

2.1 差异表达失巢凋亡基因的确定和功能分析从GEPIA数据库中共收集3005个差异表达基因DEGs(图1a),接下来将这些DEGs与从GeneCards中提取的445个失巢凋亡基因ARGs取交集,获得了82个在前列腺癌中差异表达的失巢凋亡基因(differentially expressed anoikis-related genes,DE-ARGs)(图1b)。为了进一步探索PCa中发生失巢凋亡的可能途径,我们采用clusterProlifer包对82个DE-ARGs的GO和KEGG通路分析进行处理。结果显示,在生物过程方面(biological process,BP),DE-ARGs主要与细胞黏附调控、阿米巴型细胞迁移、细胞-底物黏附等有关(图1c)。从细胞成分(cellular component,CC)来看,DE-ARGs在含有胶原蛋白的细胞外基质、细胞前缘和局部黏附中具有极大的参与性。在分子功能(molecular function,MF)方面,82个DE-ARGs主要与整合素结合、生长因子结合和蛋白酶结合。KEGG富集分析显示,这些差异基因主要富集于PI3K-Akt信号通路、局部黏附,肿瘤中的MicroRNAs,Rap1信号通路和Ras信号通路等(图1d)。

图1 差异表达的失巢凋亡相关基因(DE-ARGs)

2.2 ARGs预后标志物的构建和验证我们将419名PCa患者随机分为训练集(n=210)和验证集(n=209)。采用单变量Cox方法在训练集中获得6个与BCR生存相关的ARGs:GNE(UDP-N-acetylglucosamine 2-epimerase/N-acetylmannosamine kinase,GNE)、凋亡抑制因子5(baculoviral IAP repeat containing5,BIRC5)、前列腺素内过氧化物合酶2(prostaglandin-endoperoxide synthase 2,PTGS2)、糖基化转移酶(O-GlcNAc transferase,OGT)、肿瘤转移相关基因1(metastasis-associated protein 1,MTA1)、肺腺癌转移相关转录物1(Metastasis-associated lung adenocarcinoma transcript 1,MALAT1),见图2a。然后使用LASSO线性回归和多元Cox分析来确定预后标志物最具预测性的基因,发现3个潜在的预测因子(GNE、PTGS2、OGT)(图2b～d)。采用以下风险评分公式,风险评分=[GNE×(-0.2806))+(PTGS2×(-0.3445))+(OGT×(1.0868)]。利用风险评分的中值将病例分为高风险和低风险组。箱式图结果显示6个基因在正常组织和PCa组织中基因表达存在差异,P<0.05(图2e)。具体数值分别为:GNE(5.1e-3),BIRC5(5.2e-19),PTGS2(4.7e-11),OGT(4.0e-15),MTA1(2.4e-4),MALAT1(1.9e-5)。

图2 差异失巢凋亡相关标志物的鉴定

2.3 验证3个ARGs在无BCR生存预测中的作用在训练集中,我们指定的ARGs可以区分PCa样本的风险评分和临床状态(图3a～b)。 K-M分析显示,与低风险组相比,高危组的无进展生存期(recurrence free survival,RFS)较差(图3c)。其预测患者1年,3年和5年复发率的风险评分的ROC曲线下的面积(area under curve,AUC)分别为0.76,0.74和0.68(图3d)。应用测试集和整个队列来验证ARGs的预测性能。PCa患者在测试集和整个队列中的风险评分和临床状态的布局如图3e～f和图3i～j,两个风险组之间的预后差异也在测试集和整队列中进行验证,如图3g和图3k所示,高风险组表现出更差的预后,差异具有显著性意义,p值分别为4.996e-2,3.796e-2。ROC曲线表明ARGs在测试和整组中的预测可靠性。结果显示测试集中预测患者1年,3年和5年生化复发率的风险评分的AUC值分别为0.58,0.68和0.61,整个队列中与预测患者1年,3年和5年复发率的风险评分的AUC值分别为0.62,0.73,0.71,见图3h和图3l。表明该风险模型可用于预测预后。

图3 PCa中ARGs的预测价值

2.4 列线图的确定为了确定ARGs模型的独立意义,我们应用单变量和多变量Cox分析显示,训练集中PSA值和风险评分均可作为PCa独立预后因素,PSA值(P<0.001),风险评分(P=0.001),见图4a～b。为了进一步挖掘ARGs的预后价值,结合临床因素如年龄、T分期、N分期、PSA数值和Gleason评分。我们构建了用于预测PCa患者1、3、5年生存率的诺模图(图4c～d),训练集和测试集列线图的C-index分别为0.784和0.702,表明上述列线图相当令人满意的预测性能。此外,校准曲线的结果表明上述诺模图有良好的预测性能。(图4e～h))

图4 基于ARGs列线图的构建

2.5 ARGs模型的免疫活性分析为了探究ARGs与免疫细胞浸润之间的相关性,我们使用CIBERSORT算法来分析训练集中ARGs在22种免疫细胞中的比例(图5a)。此外,我们还使用ssGSEA进一步比较分析TCGA训练集和测试集队列中低、高风险人群中间16种免疫细胞的富集分数和13种免疫相关通路的活性。结果显示,在训练集中,高风险人群免疫细胞的浸润程度普遍较低,活化的树突状细胞,CD8+T细胞,未成熟的树突状细胞,肥大细胞,自然杀伤细胞,Th1和Th2细胞以及T调节细胞在高低风险组的分布有差异,尤其是肥大细胞和自然杀伤细胞,见图5b。与低风险组相比,高风险组的一些免疫通路下调,如抗原呈递细胞共抑制、CC趋化因子受体、Ⅱ型干扰素反应等(图5c),差异有统计学意义,如表1所示。

图5 训练集中ARGs的免疫蓝图分析

表1 ARGs在16种免疫细胞和13种免疫相关通路分析结果

3 讨论

前列腺癌是危害男性健康的恶性肿瘤之一,容易发展为去势抵抗性前列腺癌。这种疾病在早期阶段很难被诊断出来,可靠的预后生物标志物可以帮助准确评价前列腺癌患者的预后情况。失巢凋亡是一种特殊的凋亡细胞死亡方式,发生失巢凋亡抵抗的细胞可以促进肿瘤转移的发生,是肿瘤侵袭和转移的重要机制之一[11]。这种独特的死亡方式吸引了大量研究,例如,色素框同源物2(chromobox homolog 2,CBX2)被确定为高级浆液卵巢癌中发生失巢凋亡逃逸和扩散的驱动因素。CBX2的缺失可以抑制细胞增殖,降低干性,并增加顺铂敏感性[12]。环状RNAs(CircCEMIP)通过增强前列腺癌细胞的保护性自噬来促进失巢凋亡抵抗[13]。同样的,转录因子BRN2表达能够诱导c-MET水平和STAT3的磷酸化增加,从而增加在黑色素瘤中的抗失巢凋亡[14]。在肺腺癌、肝细胞癌中,癌胚抗原相关细胞黏附分子6(carcinoembryonic antigen related cell adhesion molecule6,CEACAM6),GTP酶激活蛋白1(IQ motif containing GTPase-activating protein 1,IQGAP1)的表达均能通过激活Src-FAK信号通路来抑制失巢凋亡。表明它们可作为肿瘤转移和预后指标之一[15-16]。然而,失巢凋亡及其相关基因在PCa中的具体作用尚不清楚。

在本研究中,首先通过TCGA等数据库获得前列腺癌患者的基因表达和临床数据,利用GeneCards数据库确定PCa中相关的失巢凋亡基因。为了检测差异表达失巢凋亡基因影响的潜在生物学途径,我们通过GO和KEGG富集分析发现相关基因主要激活了各种肿瘤通路。GO分析结果显示主要与细胞迁移和细胞黏附有关,KEGG分析主要在PI3K-AKT信号通路,局部粘附,Rap1信号通路,肿瘤中的MicroRNAs中富集。这些通路都参与PCa的发生发展。例如,MicroRNA-151通过靶向PI3K/AKT通路调节人前列腺癌细胞的生长,化疗敏感性和转移[17]。Rap1信号传导调节整合素和钙粘蛋白,两者在细胞对ECM的黏附和细胞间黏附中起重要作用[18]。在前列腺癌细胞中,信号诱导增殖相关1(signal-induced proliferation-associated 1,SIPA1)通过抑制Rap1介导的细胞对ECM的粘附促进肿瘤细胞侵袭和转移[19]。Talin1作为一种粘附复合物蛋白通过局部粘附信号传导和失巢凋亡抵抗促进肿瘤侵袭和转移[20],因此,推测失巢凋亡基因可能主要通过激活以上通路促进PCa生化复发的发生。

为获得PCa的生化复发预后指标,我们将所有患者分为训练集和测试集。进行逐步回归后,在训练集中开发了3个基于ARGs的预后模型,并将患者分为高风险组和低风险组。数据显示两组的预后存在差异。并在测试集和TCGA整个队列确认了3个标志物的准确性。为了进一步研究ARGs模型的功能,我们生成了一个,结合了年龄、分级、阶段和风险评分等因素的列线图。同时校准图显示列线图与预后指标的预测具有很好的拟合。

我们提出的ARGs模型与PCa患者的生存结果密切相关。包括OGT、PTGS2、GNE。OGT广泛存在于细胞核和细胞质中,参与基因表达、代谢反应、细胞应激反应等与细胞营养状况相关的反应的调控[21]。有研究发现OGT在胰腺癌肿瘤组织中表达水平大幅升高,高水平的OGTmRNA表达与PDAC患者的总体生存率低下有关[22]。此外,OGT还能促进胶质母细胞瘤中脂质代谢和肿瘤的发展,可能作为潜在的治疗靶点[23]。相反,OGT下调促进卵巢癌中顺铂耐药的发生,可作为化疗耐药新的分子标志物[24]。而我们的结果发现OGT在高风险组中表达升高。PTGS2是前列腺素的一种限速酶,PTGS2/NF-κB信号通路参与了胶质瘤细胞的放射耐受性和PCa细胞的增殖,迁移和侵袭进展[25-26]。在胃癌中表达上调,与患者生存率低有关。可用作辅助治疗靶点来逆转胃癌细胞的化学耐药性[27]。GNE作为一种激酶,是唾液酸合成过程中的主要调节因子之一[28]。GNE表达的改变可导致唾液酸化模式的改变,引起包括癌症,神经退行性疾病和感染性疾病等多种疾病表型[29]。研究发现,GNEmRNA和蛋白表达在人淋巴瘤细胞和人髓系白血病细胞中表达显著降低,富含CpG的GNE启动子异常高甲基化,导致GNE表达沉默[28,30],关于GNE和OGT基因在前列腺癌中的作用尚未明确,值得开展进一步的研究。

免疫抑制细胞可以诱导肿瘤微环境中肿瘤免疫逃逸的发生,进而促进肿瘤转移和进展。但前列腺癌是典型的“冷肿瘤”之一,肿瘤内杀伤性T细胞浸润少,抑制性免疫细胞和基质细胞比较多[31]。有研究报道Th1细胞与PCa进展最相关,T调节细胞与PCa患者预后之间存在正相关关系[32]。Th2细胞和中枢记忆性T细胞与根治性前列腺切除术后前列腺癌复发有关,是复发的独立保护因子[33-34]。在这项研究中,我们发现几种免疫细胞的组成在高风险组中表达低,包括Th2细胞、CD8+T细胞,肥大细胞,NK细胞和Treg细胞等。而PCa的典型肿瘤微环境表现为T细胞密度降低,尤其是在CD8细胞亚群中。这种从TME中排除T细胞有助于形成免疫特权位点[35]。总体而言,高风险组患者的免疫细胞在肿瘤组织中浸润的比例较低,表明高风险组患者可能抑制前列腺癌中的免疫反应。

综上所述,我们的研究发现了可能参与PCa生化复发过程的失巢凋亡相关基因的预后价值,这些基因对预测患者的BCR有重要价值。尽管我们提出的ARGs在预测PCa患者的预后方面显示出强大的性能,但本研究具有一定的局限性。首先,当前研究中分析的所有PCa病例临床队列仅从TCGA网站获得。ARGs模型的效力也需要由外部队列进行验证。其次,这是一项基于TCGA数据库的回顾性研究,缺乏吸烟史等其他因素可能影响前列腺癌预后的信息。此外,应基于实验分析来阐释ARGs基因的潜在分子机制。