基于细胞焦亡基因前列腺癌预后模型的构建及肿瘤微环境分析△

邵波，万水，田申，马园园，陈丹霞，杜沂宸

1 昭通市中医医院泌尿外科，云南 昭通 657000

芜湖市中医医院2 泌尿外科，3 皮肤科，4 肛肠科，安徽 芜湖 241000

前列腺癌（prostate cancer，PCa）是男性泌尿系统常见的恶性肿瘤之一，在中国发病率呈上升趋势，严重威胁中老年男性的生命健康[1]。PCa 发生于前列腺上皮[2]，是一种雄激素依赖性恶性肿瘤，病理类型包括导管腺癌、腺癌、鳞状细胞癌等，其中90%以上为腺癌[3]，其发病机制目前尚不清楚。细胞焦亡是近年来被发现并证实的一种新的细胞程序性死亡方式，可参与肿瘤微环境的形成，具有促进和抑制肿瘤发生、发展的双重效应，其特征为细胞膜穿孔、细胞肿胀及细胞破裂等[4-5]。目前，关于PCa 与细胞焦亡基因关系的研究较少。因此，本研究基于细胞焦亡基因构建PCa 预后模型并分析其肿瘤微环境，现报道如下。

1 资料与方法

1.1 数据获取

从癌症基因组图谱（The Cancer Genome Atlas，TCGA）数据库（https://portal.gdc.cancer.gov/）中获取502 例PCa 患者和51 例健康者的基因转录组数据、基因突变数据。从基因表达综合（Gene Expression Omnibus，GEO）数据库（https://www.ncbi.nlm.nih.gov/geo/）中获取96 例PCa 患者的临床特征及细胞焦亡基因表达数据。

1.2 获取差异表达的细胞焦亡基因

1.2.1 PCa 与细胞焦亡基因关系分析 将PCa 患者的基因突变数据导入R 语言“clusterProfiler”包中，分析肿瘤样品中基因突变的类型、频率及其碱基改变情况。通过perl 软件进一步整理基因突变数据，计算肿瘤突变负荷及拷贝数频率变化情况，分析拷贝数频率在染色体上的变化。将从TCGA数据库下载的502例PCa患者和51例健康者的转录组数据进行差异分析，提取差异表达的细胞焦亡基因。

1.2.2 预后相关基因的筛选 将从GEO 数据库下载的96 例PCa 患者的临床特征及细胞焦亡基因表达数据合并，设置参数P﹤0.001，筛选出与预后相关的细胞焦亡基因。根据风险比（hazard ratio，HR）将筛选出的预后相关基因分为高风险基因和低风险基因。HR≥1 为高风险基因，HR﹤1 为低风险基因。

1.2.3 PCa 分型及验证 根据细胞焦亡基因的表达量，设置聚类算法clusterAlg=“km”，对PCa 进行聚类分析。根据R 语言“ConsensusClusterPlus”包累积分布函数（cumulative distribution function，CDF）图的K 值确定最佳亚型分组，以CDF 图中各颜色图谱界限最清楚时的K 值为最佳分组。应用主成分分析（principal component analysis，PCA）验证上述聚类分析结果是否可靠。对细胞焦亡基因亚型的基因表达量、不同亚型患者生存情况进行差异分析，来验证分型是否可靠。

1.2.4 差异交集基因的获取及富集分析 对不同亚型PCa 患者细胞焦亡基因的表达量进行差异分析，以|logFC|﹥1，P﹤0.05 为过滤条件筛选差异交集基因。使用R 语言“clusterProfiler”包对差异交集基因的分子功能和功能途径进行基因本位（Gene Ontology，GO）功能注释和京都基因与基因组百科全书（Kyoto Encyclopedia of Genes and Genomes，KEGG）富集分析。

1.3 预后相关差异表达基因的获取、分析及验证

将上述差异交集基因、PCa 患者的临床特征进行单因素Cox 分析，以P﹤0.05 为过滤条件获得与预后相关的显著差异表达基因及其表达量。将显著差异表达基因的表达量进行聚类分析，算法及判断条件同1.2.3。将不同亚型的基因表达数据和临床数据进行差异分析，以P﹤0.001 为过滤条件分析不同亚型PCa 患者的生存情况与细胞焦亡基因的关系。

1.4 预后模型的筛选、构建及验证

1.4.1 预后模型基因的筛选及风险评分公式的验证 将502 例PCa 患者随机分为训练集（n=376）和测试集（n=126），根据训练集显著差异表达基因的表达量进行Lasso 回归分析，并做10 倍交叉验证，输出风险评分。取训练集风险评分的中位值将训练集、测试集分别分为高风险组和低风险组。将高风险组和低风险组进行生存情况、基因表达差异分析，由测试集分析结果验证训练集所构建的风险评分计算公式是否可靠。

1.4.2 预后影响因素分析、预后模型的构建及验证 采用多因素Cox 回归模型分析PCa 患者预后的独立影响因素，并构建预后模型、绘制列线图，预测PCa 患者第6、9、12 年的生存情况，并采用一致性分析验证该模型与PCa 患者的实际生存情况是否相符。采用受试者工作特征（receiver operating characteristic，ROC）曲线及曲线下面积（area under the curve，AUC）验证预后模型的可靠性，AUC≥0.65 时，表明预后模型的准确度较高。

1.5 肿瘤微环境分析

从ImmPort 数 据库（https://www.immport.org/）中下载免疫细胞数据，应用R 语言“cibersort”包分析免疫细胞含量，数据下载截止日期为2023 年1月12 日。应用R 语言通过“quantile”包分析参与预后模型构建的基因与免疫细胞含量的相关性。应用R 语言“estimate”包得出肿瘤微环境基质细胞评分、免疫细胞评分并评估PCa 样品的肿瘤纯度。从PubMed 数据库（https://pubmed.ncbi.nlm.nih.gov/）中下载泛癌干细胞数据，分析干细胞含量与风险评分之间的相关性，数据下载截止日期为2023 年1 月12 日。对高风险组和低风险组样本进行肿瘤突变负荷分析，数据来源于1.2.1 部分。

2 结果

2.1 差异表达的细胞焦亡基因

2.1.1 PCa与细胞焦亡基因的关系 突变频率最高的细胞焦亡基因分别为斑点型BTB/POZ 蛋白（speckle type BTB/POZ protein，SPOP）（11%）、肿瘤蛋白p53（tumor protein p53，TP53）（11%）、肌联蛋白（titin，TTN）（10%）；TP53、带电多泡体蛋白7（charged multivesicular body protein，CHMP）7、颗粒酶A（granzyme A，GZMA）、CHMP4C、消皮素（gasdermin，GSDM）C、含核苷酸结合寡聚化结构域1（nucleotide binding oligomerization domain containing，NOD）1 等51 个焦亡基因存在拷贝数变异；各突变基因几乎分布于人体所有染色体。差异分析得到37 个显著差异表达的细胞焦亡基因，其中B 细胞淋巴瘤/白血病-2 相关X 蛋白（B cell lymphoma/leukemia- 2- associated X protein，BAX）、CHMP2A等10 个基因在PCa 样本中表达量高，胱天蛋白酶（caspase，CASP）1、CHMP2B等27 个基因在正常样本中表达量高。（图1）

图1 PCa患者与健康者焦亡基因表达差异箱线图

2.1.2 预后相关基因的筛选 根据细胞焦亡基因的表达量及临床数据分析，获得10 个与预后相关的细胞焦亡基因，其中9 个高风险基因为BAX、CASP3、CASP6、CASP8、黑色素瘤缺乏因子2（absent in melanoma 2，AIM2）、白细胞介素（interleukin，IL）1B、IL6、含PYD 和CARD 域（PYD and card domain containing，PYCARD）、肿瘤坏死因子（tumor necrosis factor，TNF），1个低风险基因为IL1A。

2.1.3 PCa 分型及验证 聚类分析显示，根据CDF 图将K 值确定为3，将PCa 分为A、B、C 三个亚型，采用PCA 验证分型结果可靠。差异分析结果显示，不同亚型之间存在显著差异，数据可用于后续分析。

2.1.4 差异交集基因的获取及富集分析 通过差异分析获得A、B、C 三个亚型的交集基因共109个。GO 及KEGG 富集分析结果显示，差异交集基因主要通过细胞外基质组织（extracellular matrix organization，EMO）、细胞外结构组织（extracellular structure organization，ESO）、肽基酪氨酸磷酸化（peptidyl-tyrosine phosphorylation，PTP）、肽基酪氨酸修饰（peptidyl-tyrosine modification，PTM）等通路来激活蛋白酪氨酸激酶（protein tyrosine kinase，PTK）、跨膜受体蛋白激酶（transmembrane receptor protein kinase，TRPK）等分子的功能，从而影响PCa 的发生、发展。（图2）

图2 GO功能注释及KEGG富集通路分析

2.2 预后相关的显著差异表达基因的获取、分析及验证

单因素Cox 分析得到与预后相关的显著差异表达基因109 个。对不同亚型PCa 患者的细胞焦亡基因表达量、生存情况进行分析，得出不同亚型PCa 患者的生存情况与细胞焦亡基因表达量之间的关系与2.1.2 及2.1.3 相符，说明筛选出的显著差异表达基因具有可靠性，可用于后续预后模型构建。

2.3 预后模型的构建及验证

2.3.1 预后相关基因的筛选 将显著差异表达的基因纳入Lasso 回归模型进行筛选，采用10 倍交叉验证筛选出B 细胞淋巴瘤/白血病3（B-cell leukemia/lymphoma 3，BCL3）、骨形成蛋白2（bone morphogenetic protein 2，BMP2）、C-C 趋化因子受体5（C-C motif chemokine receptor 5，CCR5）、几丁质酶3 样蛋白2（chitinase 3 like 2，CHI3L2）、肝癌缺失蛋白1（deleted in liver cancer 1，DLC1）、主要组织相容性复合体Ⅱ类DQ alpha 2（major histocompatibility complex, class Ⅱ, DQ alpha 2，HLADQA2）、分泌球蛋白家族3A 成员1（secretoglobin family 3A member 1，SCGB3A1）、serpin 家族E 成员1（serpin family E member 1，SERPINE1）、溶质载体家族14 成员1（solute carrier family 14 member 1，SLC14A1）9 个与预后相关的基因，并构建预后模型。风险评分计算公式：风险值=DLC1表达量×（-3.01982305649243）+SERPINE1表达量×（ -1.79410577740395 ） +BCL3表达量 ×（ -1.33453048935853 ） +CCR5表达量 ×（ -4.46426540725949 ）+CHI3L2表达量 ×（ - 4.92826100054168）+BMP2表达量 ×（ 5.20551784803045 ）+SCGB3A1表达量 ×（ - 4.55414507965178）+SLC14A1表达量×（2.10817308576526）+HLA-DQA2表达量 ×（2.6657243859992），分别计算每例患者的风险评分。按照风险评分的中位值，将HR≥4、HR﹤4 的PCa 患者分别作为高风险组（n=376）与低风险组（n=126）。差异分析显示，训练集与测试集预后情况大致相同，且随着风险评分的升高，患者死亡数增加。由此可得，风险计算公式可靠，可用于预后模型构建。

2.3.2 预后影响因素分析、预后模型的构建及验证 多因素Cox 分析结果显示，年龄、肿瘤分期（T 分期、N 分期）、风险程度均是PCa 患者预后的独立影响因素。绘制预后模型列线图，并与患者实际生存情况进行一致性分析，结果显示，该模型预测PCa 患者第6、9、12 年的生存情况与患者的实际生存情况比较相符（图3）。ROC 曲线显示，该模型预测PCa 患者第6、9、12 年生存情况的AUC分别为0.803、0.722、0.689，表明该模型预测预后的准确度较高。

图3 PCa患者预后模型列线图

2.4 肿瘤微环境分析

2.4.1 免疫细胞相关性分析 免疫细胞含量与细胞焦亡基因存在相关性（正相关或负相关）（图4）。对高风险组与低风险组进行差异分析，结果显示，高风险组患者的基质细胞评分明显高于低风险组，差异有统计学意义（P﹤0.01）；两组免疫细胞评分、综合评分比较，差异均无统计学意义（P﹥0.05）（图5）。

图4 免疫细胞与参与预后模型构建焦亡基因关系图

图5 高风险组（n=376）和低风险组（n=126）肿瘤微环境差异分析

2.4.2 干细胞相关性分析 相关性分析结果显示，干细胞含量与风险评分呈正相关（P﹤0.01）。（图6）

图6 干细胞含量与风险评分的相关性

2.4.3 肿瘤突变负荷差异分析 差异分析结果显示，高风险组患者的肿瘤突变负荷明显高于低风险组（P﹤0.01）。（图7）

图7 高风险组（n=376）和低风险组（n=126）肿瘤突变负荷差异分析

3 讨论

PCa 属于泌尿系高发恶性肿瘤，与细胞焦亡密切相关。细胞焦亡基因主要通过雄激素受体（androgen receptor，AR）、CASP3、CASP1/GSDMD 等信号通路影响PCa 的发生和发展[6-8]。同时，也可以由紫杉醇、双氢青蒿素等药物来诱导PCa 相关的细胞焦亡[9-10]。细胞焦亡是一种不同于细胞凋亡的细胞死亡方式，通过活化CASP 蛋白来裂解GSDMD 蛋白，诱导炎症小体的表达以及炎症因子的分泌，最终导致细胞死亡[11]。炎症小体、炎症因子和GSDMD 蛋白等与肿瘤的发生、侵袭和转移密切相关[12]。

本文筛选出9 个与预后相关的显著差异表达的细胞焦亡基因，包括DLC1、BCL3、SERPINE1、CCR5、CHI3L2、SCGB3A1、HLA-DQA2、SLC14A1、BMP2。DLC1为抑癌基因，通过与其伴侣磷脂酶Cδ1、真核延伸因子、鸟苷三磷酸（guanosine triphosphate，GTP）酶活化蛋白结合，影响细胞动力，或与张力蛋白、黏着斑激酶相互作用重塑黏着斑，加强黏附蛋白的作用，从而降低肿瘤细胞转移风险[13-14]。BCL3可阻断核因子κB（nuclear factor-κB，NF-κB）信号通路，影响细胞增殖、分化、凋亡以及免疫反应等，对多数肿瘤具有促进作用，但对如PCa 等炎症性肿瘤具有明显的抑制作用[15-16]。SERPINE1、CCR5、CHI3L2、SCGB3A1基因在肺癌、直肠癌、乳腺癌等恶性肿瘤中具有促进肿瘤发生、发展的作用[17-20]。根据本文分析结果，上述基因可以降低PCa 患者的死亡风险并延长患者的生存期。分析原因可能为，基因表达物分泌至肿瘤微环境中，与炎症因子相互作用，激活了细胞焦亡基因对肿瘤的抑制作用。因此，上调上述细胞焦亡基因的表达或可成为治疗PCa 的新方法。BMP2属于转化生长因子β超家族成员的低分子糖蛋白基因，可诱导髓源性抑制细胞增殖，促进肿瘤细胞浸润[21-22]。SLC14A1编码的尿素转运蛋白B（urea transporter B，UT-B）在肿瘤组织中大多处于碱基缺失状态，因此功能低下，增加微环境肿瘤负荷可促进肿瘤的发生、发展[23]。HLA-DQA2属于白细胞抗原的一种亚型，可限制CD4+辅助性T 细胞的相关反应，降低机体抗原呈递、免疫激活等功能的效应，对恶性肿瘤的预后产生不利影响[24]。

GO 功能注释和KEGG 富集分析结果显示，细胞焦亡基因通过EMO、ESO、PTP、PTM 等通路来激活PTK、TRPK 等分子的功能，从而调控PCa 的发生、发展。EMO、ESO 信号通路包括肿瘤坏死因子样弱凋亡诱导物（TNF-like weak inducer of apoptosis，TWEAK）/成纤维细胞生长诱导因子14（fibroblast growth factor-inducible 14，Fn14）、黏着斑激酶（focal adhesion kinase，FAK）/Rho 关联含卷曲螺旋结合蛋白激酶1（Rho associated coiled-coil containing protein kinase 1，ROCK1）/Yes 相关蛋白（Yes-associated protein，YAP）、WNT/β-联蛋白（βcatenin）等，对细胞成分、功能、增殖、分化以及构建起调节作用，失衡时可引起细胞增殖、分化紊乱及密集重叠生长，促进肿瘤的发生、转移和浸润[25-27]。PTP 通路属于肽基酪氨酸修饰的方式之一。酪氨酸在调控细胞活动过程起到分子开关的作用，可导致肿瘤扩散。几乎所有的多肽因子均通过PTP 通路来激活和刺激细胞的生长、增殖，促进肿瘤的发生[28]。PTK 可催化多种含Scr 同源结构域2（Src homology domain 2，SH2）的底物蛋白磷酸化，诱导细胞增殖和分化，促进PCa、膀胱癌等多种恶性肿瘤的发生、发展[29]。TRPK 属于跨膜蛋白，膜外信号物质与跨膜受体蛋白结合后，可激活其活性，促进跨膜受体蛋白功能，如NOTCH 受体蛋白易造成肿瘤细胞增殖、分化及逃逸[30]。

本文成功构建了基于细胞焦亡基因的PCa 预后模型，并进行了肿瘤微环境分析，发现了BCL3、BMP2、CCR5、CHI3L2、DLC1、HLA-DQA2、SCGB3A1、SERPINE1、SLC14A19 个细胞焦亡基因，可预测PCa 患者的预后和免疫治疗反应，有利于进一步指导临床实践。