实用肿瘤杂志   2026, Vol. 41 Issue (3): 233-246 本刊论文版权归本刊所有,未经授权,请勿做任何形式的转载

文章信息

周美琪, 邱吉利, 龚晓楠, 陈嘉妮, 胡跃
Zhou Meiqi, Qiu Jili, Gong Xiaonan, Chen Jiani, Hu Yue
基于内质网应激相关基因特征的乳腺癌预后模型构建与功能解析
Construction and functional analysis of a prognostic model for breast cancer based on endoplasmic reticulum stress-related gene signatures
实用肿瘤杂志, 2026, 41(3): 233-246
Journal of Practical Oncology, 2026, 41(3): 233-246

基金项目

浙江省自然科学基金(ZCLQN25H1607);北京市希思科临床肿瘤学研究基金(Y-QL2019-0393);中国医药卫生事业发展基金(chmdf2025-xrky07-12)

通信作者

胡跃,Email:huyuezju@zju.edu.cn

文章历史

收稿日期:2025-09-17
基于内质网应激相关基因特征的乳腺癌预后模型构建与功能解析
周美琪 , 邱吉利 , 龚晓楠 , 陈嘉妮 , 胡跃     
浙江大学医学院附属第二医院乳腺外科,浙江 杭州 310009
摘要目的 基于癌症基因组图谱(The Cancer Genome Atlas, TCGA)和基因表达综合数据库(Gene Expression Omnibus, GEO)数据库,构建一个基于内质网应激相关基因(endoplasmic reticulum stress-related genes, ERSRGs)的乳腺癌预后风险模型,并系统解析其生物学功能与临床意义。方法 整合TCGA数据库与GEO数据库数据集GSE20685和GSE42568中共1 469例乳腺浸润癌(breast invasive carcinoma, BRCA)组织和113例癌旁组织的转录组与临床数据(TCGA-GEO BRCA队列),筛选差异表达的ERSRGs,并采用LASSO Cox回归构建预后模型。通过Kaplan-Meier法作生存分析、受试者工作特征(receiver operating characteristic, ROC)曲线、功能富集、免疫浸润、基因组变异和药物敏感性分析以评估模型性能。结果 从分子特征数据库v7.0获得272个ERSRGs,基于TCGA-GEO BRCA队列在其中筛选出15个在乳腺癌组织中差异表达的基因[|log2 fold change (FC)| > 1.0,P < 0.05]。通过LASSO Cox回归在其中筛选出8个预后关键基因用于构建预后风险模型,分别为cAMP反应元件结合蛋白3样1(cAMP responsive element binding protein 3 like 1, CREB3L1)、基质细胞衍生因子2样1(stromal cell derived factor 2 like 1, SDF2L1)、磷酸肌醇-3-激酶调节亚基1(phosphoinositide-3-kinase regulatory subunit 1, PIK3R1)、失活同源物2相互作用蛋白(disabled homolog 2 interacting protein, DAB2IP)、蛋白磷酸酶1调节亚基15A(protein phosphatase 1 regulatory subunit 15A, PPP1R15A)、佛波醇-12-肉豆蔻酸-13-乙酸酯诱导蛋白1(phorbol-12-myristate-13-acetate-induced protein 1, PMAIP1)、丝裂原活化蛋白激酶激酶激酶5(mitogen-activated protein kinase kinase kinase 5, MAP3K5)和1, 4, 5-三磷酸肌醇受体1型(inositol 1, 4, 5-trisphosphate receptor type 1, ITPR1)。以LASSO Cox回归筛选所得的8个关键预后基因的回归系数与其归一化mRNA表达值相乘并求和,建立风险评分公式,并使用风险评分的中位数将TCGA-GEO BRCA队列的BRCA患者分为高风险组(n=734)和低风险组(n=735)。生存分析显示,高风险组患者总生存期(overall survival, OS)缩短(5年OS率:76.8% vs 88.8%,10年OS率:59.0% vs 77.6%;log-rank P < 0.01)。ROC曲线分析显示,风险模型预测1、3和5年OS的曲线下面积(area under the curve, AUC)分别为0.645、0.671和0.689。功能分析提示,高风险组差异表达基因富集于雌激素信号通路、E2F靶点、IL-17信号通路和炎性反应等。免疫浸润分析显示,高风险组免疫评分与基质评分均升高(均P < 0.01)。药物敏感性分析表明,高风险组对18种药物的半数抑制浓度(half maximal inhibitory concentration, IC50)值均高于低风险组(均P < 0.01),尤其是细胞周期蛋白依赖性激酶4/6(cyclin-dependent kinase 4/6, CDK4/6)抑制剂哌柏西利、哺乳动物雷帕霉素靶蛋白(mammalian target of rapamycin, mTOR)抑制剂替西罗莫司与AZD8055。结论 本研究成功构建一个基于ERSRGs的乳腺癌预后模型。该模型具有良好的预测性能,可能为乳腺癌患者的风险分层与个体化治疗提供参考。
关键词乳腺癌    内质网应激相关基因    预后风险模型    列线图    生物信息学    
Construction and functional analysis of a prognostic model for breast cancer based on endoplasmic reticulum stress-related gene signatures
Zhou Meiqi , Qiu Jili , Gong Xiaonan , Chen Jiani , Hu Yue     
Department of Breast Surgery, the Second Affiliated Hospital, Zhejiang University School of Medicine, Hangzhou 310009, China
Abstract: Objective To construct a prognostic risk model for breast cancer based on endoplasmic reticulum stress-related genes (ERSRGs) using The Cancer Genome Atlas (TCGA) and Gene Expression Omnibus (GEO) databases, and to systematically analyze its biological functions and clinical significance. Methods Transcriptomic and clinical data of 1 469 breast invasive carcinoma (BRCA) tissues and 113 adjacent normal tissues from TCGA database and the datasets, GSE20685 and GSE42568, from GEO database were integrated as the TCGA-GEO BRCA cohort. Differentially expressed ERSRGs were identified, and a prognostic model was constructed using LASSO-Cox regression. Model performance was evaluated by Kaplan-Meier survival analysis, receiver operating characteristic (ROC) curves, functional enrichment analysis, immune infiltration analysis, genomic variation analysis, and drug sensitivity analysis. Results A total of 272 ERSRGs were obtained from the Molecular Signatures Database (MSigDB) v7.0, of which 15 were differentially expressed in BRCA tissues from the TCGA-GEO BRCA cohort [|log2 fold change (FC)| > 1.0, P < 0.05]. Eight key prognostic genes were screened by LASSO-Cox regression to construct the prognostic risk model, including cAMP responsive element binding protein 3 like 1 (CREB3L1), stromal cell derived factor 2 like 1 (SDF2L1), phosphoinositide-3-kinase regulatory subunit 1 (PIK3R1), disabled homolog 2 interacting protein (DAB2IP), protein phosphatase 1 regulatory subunit 15A (PPP1R15A), phorbol-12-myristate-13-acetate-induced protein 1 (PMAIP1), mitogen-activated protein kinase kinase kinase 5 (MAP3K5), and inositol 1, 4, 5-trisphosphate receptor type 1 (ITPR1). A risk score formula was established using the normalized mRNA expression values of these genes weighted by the regression coefficients from the LASSO-Cox analysis, and the BRCA patients in the TCGA-GEO BRCA cohort were divided into a high-risk group (n=734) and a low-risk group (n=735) using the median risk score. Survival analysis showed that patients in the high-risk group had shortened overall survival (OS) (5-year OS rates: 76.8% vs 88.8%; 10-year OS rates: 59.0% vs 77.6%; log-rank P < 0.01). ROC curve analysis showed that the areas under the curves (AUCs) of the risk model for predicting 1-, 3-, and 5-year OS were 0.645, 0.671, and 0.689, respectively. Functional analysis suggested that differentially expressed genes in the high-risk group were enriched in estrogen response, E2F targets, interleukin-17 signaling pathway, and inflammatory response. Immune infiltration analysis showed that both immune scores and stromal scores were elevated in the high-risk group (both P < 0.01). Drug sensitivity analysis indicated that the high-risk group had higher half maximal inhibitory concentration (IC50) values for 18 drugs, especially the cyclin-dependent kinase 4/6 (CDK4/6) inhibitor palbociclib and the mammalian target of rapamycin (mTOR) inhibitors temsirolimus and AZD8055, compared to the low-risk group (all P < 0.01). Conclusions This study successfully constructed a prognostic model for breast cancer based on ERSRGs. The model has good predictive performance and may provide a reference for risk stratification and individualized treatment of breast cancer patients.
Key words: breast cancer    endoplasmic reticulum stress-related genes    prognostic risk model    nomogram    bioinformatics    

乳腺癌是全球女性最常见的恶性肿瘤,也是女性癌症相关死亡的主要原因之一[1]。在中国,乳腺癌发病率持续上升,2022年粗发病率达59.70/10万,位居女性恶性肿瘤发病首位[2]。随着早期诊断、精准分型和个体化治疗的进步,乳腺癌死亡率已显著降低[3]。乳腺癌的分子分型不仅能够区分患者预后,更是指导精准治疗的关键。Oncotype DX和MammaPrint等多基因检测工具已在预后预测与治疗决策中广泛应用[4]。然而,现有模型多集中于增殖和激素受体等经典通路,对内质网应激(endoplasmic reticulum stress, ERS)这一在肿瘤发生、代谢重编程、免疫逃逸和治疗抵抗中起关键作用的生物学过程尚未系统整合。

ERS作为细胞应激的核心,在维持蛋白质稳态、决定细胞命运和重塑肿瘤微环境中扮演关键角色[5]。鉴于肿瘤细胞与正常细胞对ERS敏感性的差异,ERS通路中的关键调控分子已成为潜在的抗肿瘤治疗新靶点。在乳腺癌进展过程中,ERS通过激活未折叠蛋白反应(unfolded protein response, UPR)调控细胞存活与死亡。研究表明,ERS关键通路(如PERK/eIF2α/ATF4和IRE1α/XBP1)的异常激活可促进耐药、转移与免疫逃逸;而靶向ERS可通过诱导凋亡或自噬抑制肿瘤生长[6]。然而,目前尚缺乏基于ERS相关基因(ERS-related genes, ERSRGs)的系统性预后模型,也未能将其与基因组变异、免疫微环境和药物敏感性进行多维度整合分析。

本研究通过R语言整合癌症基因组图谱(The Cancer Genome Atlas, TCGA)与基因表达综合数据库(Gene Expression Omnibus, GEO)数据库,构建基于ERSRGs的预后风险评分模型并深入探索其与肿瘤突变负荷(tumor mutation burden, TMB)、免疫浸润、信号通路活性和药物敏感性的关联,并通过蛋白质-蛋白质相互作用(Protein-Protein Interaction, PPI)网络挖掘核心调控机制,全面证实ERSRGs在乳腺癌预后、药物敏感性和免疫微环境中的关键作用,以期为乳腺癌的风险分层与精准治疗提供新的分子工具,深化对ERS分子机制的理解,并为个体化预后预测提供潜在生物标志物。

1 资料与方法 1.1 数据获取

从TCGA官网(https://portal.gdc.cancer.gov/)获取截至2025年3月21日的乳腺浸润癌(breast invasive carcinoma, BRCA)患者每百万映射读取中每千碱基转录本片段数(fragments per kilobase of transcript per million mapped reads, FPKM)格式的转录组数据,并转换为每百万转录本数(transcripts per million, TPM)值和临床信息(包括性别、年龄、TNM分期与生存数据),共涵盖1 038例肿瘤组织[年龄26~90岁,(58.5±13.2)岁;男性12例,女性1 026例;Ⅰ期173例,Ⅱ期591例,Ⅲ期233例,Ⅳ期20例,临床分期不明21例]与113例癌旁组织。同时,下载TCGA-BRCA患者的体细胞拷贝数变异(copy number variation, CNV)和体细胞突变数据用于计算TMB,并使用maftools包[7]进行可视化。从GEO数据库下载GSE20685[8]与GSE42568[9]数据集(平台:GPL570[HG-U133_Plus_2]),包含基因表达谱和临床病理信息。其中GSE20685数据集包含327例女性BRCA组织样本[年龄24~84岁,(47.9±10.7)岁],GSE42568数据集包含104例女性BRCA组织样本[年龄31~90岁,(58.8±11.8)岁]。这2个数据集与TCGA数据合并成TCGA-GEO BRCA队列(1 469例BRCA组织和113例癌旁组织),并应用limma[10]与sva[11]包进行归一化与批次效应校正。

1.2 风险模型构建

从分子特征数据库v7.0(Molecular Signature Database, MSigDB)[12]下载2个ERSRGs集合(GO RESPONSE TO ENDOPLASMIC RETICULUM STRESS和GO REGULATION OF RESPONSE TO ENDOPLASMIC RETICULUM STRESS),去重后得到272个ERSRGs。分析ERSRGs在癌与癌旁组织中的差异表达、共表达网络和基因组变异,从而解析ERSRGs在BRCA中的特征。基于TCGA-GEO BRCA队列,以|log2 fold change(FC)| > 1.0且错误发现率(false discovery rate, FDR) < 0.05为阈值,筛选出15个在BRCA组织中差异表达的ERSRGs。采用LASSO Cox回归对上述15个差异表达的ERSRGs进行特征筛选,通过10折交叉验证确定最优惩罚参数λ,最终筛选出8个关键预后基因。以各基因的回归系数与其归一化表达值乘积之和构建风险评分公式,并依据中位风险评分将患者划分为高风险组和低风险组。

1.3 差异表达分析与功能富集分析

为鉴定与ERS风险模型相关的基因,采用R的limma包[10]分析TCGA-GEO BRCA队列的BRCA组织的高低风险组间的差异表达基因(differentially expressed genes, DEGs),筛选标准为|log2 FC| > 0.5且FDR < 0.05,结果通过热图与火山图展示。基于筛选出的DEGs,利用clusterProfiler包[13]进行基因本体论(Gene Ontology, GO)与京都基因与基因组百科全书(Kyoto Encyclopedia of Genes and Genomes, KEGG)富集分析,设定FDR < 0.05为显著性阈值。同时,基于MSigDB的“h.all.v7.2.symbols.gmt”基因集,对全部基因表达数据开展基因集富集分析(Gene-Set Enrichment Analysis, GSEA)[12],校正P < 0.05为差异具有统计学意义。

1.4 免疫浸润分析

应用单样本GSEA(single-sample GSEA, ssGSEA)算法评估BRCA患者的免疫细胞浸润水平。分析采用Bindea等[14]定义的28种免疫细胞标志基因集,通过R的RSVA包[12]计算富集分数以量化各免疫细胞亚型的相对丰度。同时,使用ESTIMATE算法[15]计算肿瘤样本的免疫与基质评分,以评估肿瘤微环境。

1.5 CNV分析

为解析TCGA-BRCA患者不同风险组间的CNV特征,基于TCGA biolinks包获取的Masked Copy Number Segment数据,采用GISTIC 2.0算法(GenePattern平台)进行系统性分析。关键参数设置为置信度0.99并保留X染色体。最终结果通过R的Maftools包实现可视化。

1.6 药物敏感性分析

基于抗癌药物敏感性基因组学数据库(Genomics of Drug Sensitivity in Cancer, GDSC;https://www.cancerrxgene.org/[16]的细胞系基因突变与不同抗癌药物半数抑制浓度(half maximal inhibitory concentration, IC50)数据,利用R的oncoPredict包[17]预测并比较TCGA-GEO BRCA队列高低风险组对138种抗癌药物的敏感性差异。

1.7 临床预测模型构建

整合ERS相关的风险评分与临床病理特征,构建预测TCGA-GEO BRCA队列中BRCA患者总生存期(overall survival, OS)的临床预测列线图。采用Harrell's C-index评估模型的区分能力,并通过校准曲线检验其预测准确性。

1.8 核心基因筛选

基于STRING数据库[18](交互评分 > 0.4)构建PPI网络,并使用Cytoscape(v3.7.2)软件进行可视化。通过CytoHubba插件[19-20]中的最大集团中心性(Maximal Clique Centrality, MCC)算法筛选枢纽基因(hub genes),并取MCC值最高的前8个基因定义为最终的核心基因集。

1.9 统计学分析

采用R(v3.6.2)软件进行统计学分析。组间比较采用以下方法:连续变量符合正态分布时使用Student t检验,否则使用Mann-Whitney U检验;分类变量使用χ2检验或Fisher精确检验。生存分析基于Kaplan-Meier法(R的survival包),组间差异通过对数秩检验评估。利用单因素和多因素Cox回归模型识别独立预后因素。采用Pearson相关系数(R的corrplot包)评估差异表达ERSRGs在BRCA组织中两两之间的表达相关性。通过时间依赖性ROC曲线分析(pROC包)[21]评估模型的预测准确性。所有统计学分析均为双侧检验,以P < 0.05为差异具有统计学意义。

2 结果 2.1 ERSRGs在BRCA患者中的表达和突变概况

对GSE20685、GSE42568和TCGA数据集进行标准化、合并和批次效应校正后形成TCGA-GEO BRCA队列,系统分析ERSRGs在BRCA样本中的mRNA表达、单核苷酸变异和CNV。基因表达热图显示,ERSRGs在TCGA-GEO BRCA队列肿瘤组织与癌旁组织中存在差异表达(图 1A)。以|log2 FC| > 1.0且FDR < 0.05为阈值,共鉴定出15个差异表达的ERSRGs(图 1B)。

注  A:TCGA-GEO BRCA队列肿瘤组织(n=1 469)与正常乳腺组织(n=113)中ERSRGs的mRNA表达热图;B:TCGA-GEO BRCA队列样本中15个差异表达ERSRGs在肿瘤组织(n=1 469)与正常乳腺组织(n=113)中的mRNA表达箱线图;C:TCGA-BRCA患者(n=536)中ERSRGs的体细胞突变图谱;D:TCGA-BRCA患者(n=536)中ERSRGs在23对染色体上的拷贝数变异概览;**P < 0.01 图 1 BRCA患者ERSRGs表达谱与基因组变异 Fig.1 Overall expression profile and genomic alterations of ERSRGs in BRCA patients

对有突变数据的536例TCGA-BRCA样本单核苷酸多态性(single nucleotide polymorphism, SNP)分析显示,432例(80.6%)存在ERSRGs突变,突变频率以肿瘤蛋白p53(tumor protein p53, TP53)基因最高,其次为磷酸肌醇-3-激酶调节亚基1(phosphoinositide-3-kinase regulatory subunit 1, PIK3R1)基因(图 1C)。TCGA-BRCA样本CNV分析表明,ERSRGs变异以拷贝数扩增为主,少数基因存在拷贝数缺失,提示ERSRGs基因组不稳定性在BRCA发生和发展中可能发挥重要作用(图 1D)。

基于TCGA-GEO BRCA队列的相关性分析显示,15个差异表达的ERSRGs中有19对基因在BRCA组织中表达呈正相关(均r > 0.5,均P < 0.01;图 2A)。基于这些差异基因,LASSO Cox回归进一步筛选出8个预后关键基因(图 2B~2C),包括cAMP反应元件结合蛋白3样1(cAMP responsive element binding protein 3 like 1, CREB3L1)、基质细胞衍生因子2样1(stromal cell derived factor 2 like 1, SDF2L1)、PIK3R1、失活同源物2相互作用蛋白(disabled homolog 2 interacting protein, DAB2IP)、蛋白磷酸酶1调节亚基15A(protein phosphatase 1 regulatory subunit 15A, PPP1R15A)、佛波醇-12-肉豆蔻酸-13-乙酸酯诱导蛋白1(phorbol-12-myristate-13-acetate-induced protein 1, PMAIP1)、丝裂原活化蛋白激酶激酶激酶5(mitogen-activated protein kinase kinase kinase 5, MAP3K5)和1, 4, 5-三磷酸肌醇受体1型(inositol 1, 4, 5-trisphosphate receptor type 1, ITPR1)。将上述8个预后关键基因的回归系数与其归一化mRNA表达值(X)相乘并求和,建立风险评分公式:风险评分= 0.102 5×XCREB3L1+(-0.166 3)×XSDF2L1+(-0.155 8)× XPIK3R1+0.291 0×XDAB2IP+0.000 1×XPPP1R15A+(-0.102 8)×XPMAIP1+0.165 5×XMAP3K5+(-0.125 5)× XITPR1,计算每例患者的风险评分,并按其中位数(0.274 3)将TCGA-GEO BRCA队列的BRCA患者分为高风险组(n=734)和低风险组(n=735)。Kaplan-Meier生存分析显示,高风险组患者OS短于低风险组(log-rank P < 0.01,图 2D)。高风险组和低风险组5年OS率分别为76.8%和88.8%,10年OS率为59.0%和77.6%,高风险组OS率更低(HR=2.057,95% CI:1.603~2.640)。时间依赖性ROC曲线分析显示,模型预测1、3和5年OS的曲线下面积(area under the curve, AUC)分别为0.645、0.671和0.689(图 2E)。

注  A:基于TCGA-GEO BRCA队列BRCA组织(n=1 469)构建的差异表达ERSRGs的相关性分析表达矩阵图;B:LASSO Cox回归分析中15个BRCA组织差异表达ERSRGs的回归系数轨迹图;C:15个BRCA组织差异表达ERSRGs的LASSO Cox回归10折交叉验证曲线图,左侧虚线对应λ_min,右侧虚线对应λ_1se,以λ_min为最优参数筛选出8个预后关键基因;D:整合TCGA-GEO BRCA队列中高风险组(n=734)与低风险组(n=735)BRCA患者的Kaplan-Meier总生存曲线;E:基于ERS相关的风险评分预测BRCA患者总生存率的时间依赖性ROC曲线;F:ERS相关风险评分模型的患者风险评分分布图、生存状态散点图和高低风险组中模型基因mRNA表达热图;AUC:曲线下面积(area under the curve) 图 2 BRCA患者ERS相关预后风险评分模型的构建和评估 Fig.2 Construction and evaluation of an ERS-related prognostic risk scoring model for BRCA patients
2.2 高低风险组的DEGs分析

ERS相关的高低风险组间的差异表达分析共鉴定出241个DEGs,包括62个上调基因和179个下调基因(图 3A~3B)。GO功能富集分析表明,这些DEGs主要富集于多细胞生物稳态和胰岛素样生长因子受体信号通路等过程(图 3C)。KEGG通路分析显示,DEGs在IL-17信号通路、松弛素信号通路、雌激素信号通路和AMPK信号通路等多条肿瘤相关通路中显著富集(图 3D)。GSEA进一步提示,高风险BRCA组织中,内质网腔、滑面内质网、炎性反应、E2F靶点和G2M期检查点等通路显著激活,而ERS正向调控在低风险患者中显著富集(图 4A~4F)。

注  A:TCGA-GEO BRCA队列中高风险组(n=734)与低风险组(n=735)间DEGs的火山图;B:TCGA-GEO BRCA队列中高风险组(n=734)和低风险组(n=735)间DEGs的热图;C:DEGs的GO功能富集分析;D:DEGs的KEGG通路富集分析 图 3 基于TCGA-GEO BRCA队列患者ERS相关的预后风险模型DEGs的GO和KEGG富集分析 Fig.3 GO and KEGG enrichment analysis of DEGs based on the ERS-related prognostic risk model in the TCGA-GEO BRCA patients
注A  :内质网腔的GSEA富集曲线;B:滑面内质网的GSEA富集曲线;C:内质网应激反应的正向调控的GSEA富集曲线;D:炎性反应特征通路的GSEA富集曲线;E:E2F靶基因特征通路的GSEA富集曲线;F:G2M期检查点特征通路的GSEA富集曲线;GO:基因本体论(Gene Ontology);NES:标准化富集得分(normalized enrichment score) 图 4 基于TCGA-GEO BRCA队列患者ERS相关的预后风险模型DEGs的GSEA分析 Fig.4 GSEA analysis of DEGs based on the ERS-related prognostic risk model in the TCGA-GEO BRCA patients
2.3 PPI网络图和相关调控网络的构建

基于STRING数据库构建ERS相关的高低风险组DEGs的PPI网络,并导入Cytoscape进行可视化(图 5A)。进一步利用CytoHubba插件在网络中筛选局部高密度区域,在枢纽基因中识别出MCC值最高的8个PPI调控网络核心基因(图 5B)。Friends分析结果显示,在这些PPI调控网络核心基因中,胰岛素受体底物1(insulin receptor substrate 1, IRS1)处于关键中心位置(图 5C)。为验证其表达特征,基于TCGA泛癌表达谱数据分析结果显示,IRS1基因在多种肿瘤与对应癌旁组织间的表达差异均具有统计学意义(均P < 0.05,图 5D)。

注  A:基于STRING数据库构建的ERS相关预后风险模型DEGs的PPI网络;B:由CytoHubba插件识别ERS相关的枢纽基因子网络中前8个PPI调控网络核心基因;C:Friends分析法对ERS相关的PPI调控网络核心基因重要性的评估结果,突出显示IRS1处于网络的关键中心位置;D:基于TCGA泛癌数据集分析IRS1基因在多种肿瘤与癌旁正常组织中的差异表达;*P < 0.05;**P < 0.01 图 5 基于TCGA-GEO BRCA队列BRCA患者ERS相关预后风险模型DEGs的PPI网络分析 Fig.5 PPI network analysis of DEGs based on the ERS-related prognostic risk model in the TCGA-GEO BRCA patients
2.4 免疫细胞浸润分析

基于ERS相关风险评分对TCGA-GEO BRCA队列进行免疫浸润分析显示,高风险组患者的免疫评分与基质评分均高于低风险组(均P < 0.01,图 6A~6B)。吲哚胺2, 3-双加氧酶1(indoleamine 2, 3-dioxygenase 1, IDO1)、趋化因子配体10(chemokine C-X-C motif ligand 10, CXCL10)、穿孔素1(perforin 1, PRF1)和T盒转录因子2(T-box transcription factor 2, TBX2)等多种免疫相关基因在高风险组均高表达(均P < 0.05,图 6C)。ssGSEA分析进一步揭示28种免疫细胞亚型在两组间的浸润模式(图 6D);其富集分数间多呈正相关(图 6E)。其中,巨噬细胞、髓源性抑制细胞(myeloid-derived suppressor cell, MDSC)和调节性T细胞等多种抑制性免疫细胞在高风险组高表达(均P < 0.05,图 6F),提示高风险组的免疫微环境以免疫抑制为核心特征。

注  A:高低风险组BRCA患者免疫评分的比较;B:高低风险组BRCA患者基质评分的比较;C:免疫治疗关键分子在高低风险组间的表达差异;D:28种免疫细胞浸润水平的ssGSEA富集分数热图;E:ssGSEA分析不同免疫细胞浸润水平的相关性矩阵图;F:ssGSEA分析28种免疫细胞在高低风险组间浸润水平的差异;高风险组734例,低风险组735例;*P < 0.05;**P < 0.01 图 6 基于TCGA-GEO BRCA队列分析ERS相关预后风险评分与不同免疫细胞浸润的关联 Fig.6 Association between the ERS-related prognostic risk score and different immune cell infiltration in the TCGA-GEO BRCA patients
2.5 ERS相关的风险评分与基因组变异特征

对具有突变数据的536例TCGA-BRCA样本进行驱动基因突变谱分析显示,高低风险组间的常见突变基因存在明显差异(图 7A)。整体层面,高风险组具有更高的TMB(P=0.009 1,图 7B)和肿瘤新生抗原负荷(P=0.006 2,图 7C)。进一步利用GISTIC 2.0对TCGA-BRCA拷贝数数据进行分析显示,与低风险组比较,高风险组存在更为广泛的拷贝数扩增事件,而低风险组则相对富集于部分拷贝数缺失区域(图 7D~7E)。

注  A:高低风险组驱动基因突变谱;B:高低风险组肿瘤突变负荷比较;C:高低风险组肿瘤新生抗原负荷比较;D:基于TCGA-BRCA高风险组样本的GISTIC 2.0分析展示显著拷贝数扩增(红色)与缺失(蓝色)区域;E:基于TCGA-BRCA低风险组样本的GISTIC 2.0分析展示显著拷贝数扩增(红色)与缺失(蓝色)区域;高低风险组各268例 图 7 536例有突变数据的TCGA-BRCA患者不同ERS相关风险分组的遗传变异分析 Fig.7 Genetic-alteration landscape of 536 TCGA-BRCA patients with available somatic mutation data stratified by ERS-related risk
2.6 基于ERS相关风险评分的BRCA患者药物敏感性分析

基于GDSC数据库,利用R的oncoPredict包通过细胞系药物反应数据预测138种抗癌药物在TCGA-GEO BRCA队列BRCA患者中的IC50值。结果显示,高风险组对18种药物的IC50值更高(均P < 0.01,图 8),提示药物敏感性降低,包括细胞周期蛋白依赖性激酶4/6(cyclin-dependent kinase 4/6, CDK4/6)抑制剂哌柏西利、哺乳动物雷帕霉素靶蛋白(mammalian target of rapamycin, mTOR)抑制剂替西罗莫司与AZD8055。高风险组对长春瑞滨和鼠双微体2(mouse double minute 2, MDM2)抑制剂Nutlin.3a等部分化疗与靶向药物亦呈耐药趋势。

注  IC50:半数抑制浓度(half maximal inhibitory concentration);高风险组734例,低风险组735例;**P < 0.01 图 8 TCGA-GEO BRCA队列BRCA患者ERS相关预后风险分层的药物敏感性分析 Fig.8 Drug sensitivity analysis based on ERS-related prognostic risk stratification in the TCGA-GEO BRCA patients
2.7 临床预测模型的构建

分析TCGA-BRCA样本显示,ERS相关风险评分在年龄和性别方面比较,差异均无统计学意义(均P > 0.05,图 9A~9B),临床分期Ⅳ期患者风险评分高于Ⅰ~Ⅲ期者(均P < 0.05,图 9C)。多因素Cox回归分析显示,ERS相关风险评分是TCGA-BRCA患者OS的独立危险因素(P=0.008,图 9D)。整合风险评分与临床病理特征构建预测TCGA-BRCA患者1、3和5年OS概率的列线图(图 9E)。模型区分度C-index为0.691(95% CI:0.644~0.738),经Bootstrap法重复抽样1 000次进行内验证,校准曲线显示1、3和5年OS率的预测值与实际观测值吻合良好(图 9F)。

注  A:TCGA-BRCA队列中不同年龄组患者ERS相关的风险评分的分布;B:TCGA-BRCA队列中ERS相关的风险评分在不同性别中的分布;C:TCGA-BRCA队列中ERS相关的风险评分在不同临床分期中的分布;D:TCGA-BRCA队列ERS相关的风险评分和临床特征对BRCA患者预后的影响的多因素Cox回归分析;E:整合ERS相关的风险评分与临床特征的BRCA患者预后列线图;F:列线图预测BRCA患者1、3和5年总生存率的校准曲线(Bootstrap重复1 000次);AIC:赤池信息准则(Akaike information criterion) 图 9 ERS相关的风险评分对BRCA患者预后的预测能力分析 Fig.9 Prognostic power of the ERS-derived risk signature for BRCA patients
3 讨论

本研究构建了包含8个关键基因的ERS预后模型。尽管预测精度仍有提升空间,但其核心价值在于首次从ERS这一维度整合代谢重编程、治疗抵抗与免疫微环境调控,为BRCA预后预测提供全新机制视角。

模型中的8个基因通过协同调控ERS-UPR信号轴,驱动BRCA的恶性进展。CREB3L1作为UPR关键传感器,在不同亚型中的功能异质性提示ERS调控的复杂性[22-23]。而PIK3R1[24]和DAB2IP[25-26]侧重于通路调节,前者低表达与内分泌耐药相关,后者通过抑制RAC1/β-catenin轴减少化疗耐药。此外,PPP1R15A表达水平可能决定ERS持续时间和最终的细胞命运[27]。值得注意的是,ITPR1通过调控钙稳态连接ERS与自噬过程,而PMAIP1则整合ERS与线粒体凋亡通路[28-29],揭示ERS与细胞死亡调控之间的复杂网络[30]

ERS高风险评分与TMB、肿瘤新生抗原和CNV扩增密切相关。GSEA分析显示高风险组富集E2F靶点[31]、G2M检查点和IL-17信号通路[32],提示细胞周期失调与免疫逃逸是高危肿瘤的重要特征。值得注意的是,AMPK通路在不同应激条件下的功能转换值得进一步探讨[33]。PPI网络的构建和Friends分析揭示IRS1处于信号调控中心,作为胰岛素样生长因子1受体(insulin-like growth factor 1 receptor, IGF1R)和叉头框蛋白A1(forkhead box A1, FOXA1)的关键衔接蛋白[31, 34-35],在三阴性乳腺癌中低表达与不良预后独立相关[36]。而在激素受体阳性乳腺癌中则与癌干细胞样特征和内分泌耐药协同,展现出明显的亚型异质性[37],是潜在的预后标志物和治疗靶点。

药物敏感性分析高风险组对CDK4/6抑制剂(哌柏西利)和mTOR抑制剂(替西罗莫司和AZD8055)相对不敏感,推测可能与E2F通路异常激活和PI3K-AKT通路(如IRS1和PIK3R1)上调有关。这种分类分析不仅揭示ERS模型在预测传统靶向药物反应方面的价值,更为克服耐药性和开发新的治疗策略提供方向,尤其是近年来发病率显著增加的早发性乳腺癌,迫切需要新型药物改善预后[38]

在免疫微环境方面,高风险组呈现出显著的“免疫悖论”:即高水平的免疫浸润未能转化为生存优势,反而表现出以MDSC和巨噬细胞富集为核心的免疫抑制特征。这可能源于ERS诱导的IL-17通路激活及其驱动的促炎性抑制环境。同时,IDO1等检查点分子的同步高表达提示T细胞处于功能耗竭状态。因此,针对高风险患者,单纯免疫激活疗效可能受限,联合干预MDSC募集通路或IDO1等靶点或是提升个体化疗效的新策略[39-41]

本研究主要基于公共数据库进行回顾性分析,受数据集异质性影响,模型的稳健性仍需在多中心、大规模前瞻性临床队列中进一步验证。此外,目前研究尚处于描述性关联阶段,缺乏对核心基因的体内外功能实验支撑。未来研究将聚焦以下方向:通过细胞和动物模型探讨ERS相关基因在乳腺癌恶性进展与耐药中的具体分子通路;深入解析ERS信号轴与糖脂代谢重构及免疫抑制微环境间的复杂交互网络;探索风险评分与免疫及化疗响应的关联,为乳腺癌的精准分层与个体化治疗提供坚实的实验证据。

综上所述,本模型不仅为BRCA预后分层提供新工具,更通过系统解析ERS相关的基因调控网络,揭示其在重塑免疫微环境与介导治疗耐药中的关键作用,为精准治疗方案的制定提供重要线索,使其有望成为现行临床病理评价体系的有益补充。临床应用中,该ERS风险评分可作为一种辅助决策工具,提示医生关注高风险患者潜在的内分泌和化疗耐药性,并为开发针对免疫微环境调节的个体化治疗策略提供重要的理论依据。

利益冲突  所有作者声明无利益冲突

参考文献
[1]
Siegel RL, Kratzer TB, Giaquinto AN, et al. Cancer statistics, 2025[J]. CA Cancer J Clin, 2025, 75(1): 10-45.
[2]
Han BF, Zheng RS, Zeng HM, et al. Cancer incidence and mortality in China, 2022[J]. J Natl Cancer Cent, 2024, 4(1): 47-53.
[3]
Wilcock P, Webster RM. The breast cancer drug market[J]. Nat Rev Drug Discov, 2021, 20(5): 339-340. DOI:10.1038/d41573-021-00018-6
[4]
Delaloge S, Piccart M, Rutgers E, et al. Standard anthracycline based versus docetaxel-capecitabine in early high clinical and/or genomic risk breast cancer in the EORTC 10041/BIG 3-04 MINDACT phase Ⅲ trial[J]. J Clin Oncol, 2020, 38(11): 1186-1197. DOI:10.1200/JCO.19.01371
[5]
Zhang W, Shi Y, Oyang L, et al. Endoplasmic reticulum stress-a key guardian in cancer[J]. Cell Death Discov, 2024, 10: 343.
[6]
Xu D, Liu Z, Liang MX, et al. Endoplasmic reticulum stress targeted therapy for breast cancer[J]. Cell Commun Signal, 2022, 20(1): 174. DOI:10.1186/s12964-022-00964-7
[7]
Mayakonda A, Lin DC, Assenov Y, et al. Maftools: efficient and comprehensive analysis of somatic variants in cancer[J]. Genome Res, 2018, 28(11): 1747-1756. DOI:10.1101/gr.239244.118
[8]
Kao KJ, Chang KM, Hsu HC, et al. Correlation of microarray-based breast cancer molecular subtypes and clinical outcomes: implications for treatment optimization[J]. BMC Cancer, 2011, 11(1): 143. DOI:10.1186/1471-2407-11-143
[9]
Clarke C, Madden SF, Doolan P, et al. Correlating transcriptional networks to breast cancer survival: a large-scale coexpression analysis[J]. Carcinogenesis, 2013, 34(10): 2300-2308. DOI:10.1093/carcin/bgt208
[10]
Ritchie ME, Phipson B, Wu D, et al. Limma powers differential expression analyses for RNA-sequencing and microarray studies[J]. Nucleic Acids Res, 2015, 43(7): e47. DOI:10.1093/nar/gkv007
[11]
Chen SL, Yang D, Lei CX, et al. Identification of crucial genes in abdominal aortic aneurysm by WGCNA[J]. PeerJ, 2019, 7: e7873. DOI:10.7717/peerj.7873
[12]
Subramanian A, Tamayo P, Mootha VK, et al. Gene set enrichment analysis: a knowledge-based approach for interpreting genome-wide expression profiles[J]. Proc Natl Acad Sci U S A, 2005, 102(43): 15545-15550. DOI:10.1073/pnas.0506580102
[13]
Yu GC, Wang LG, Han YY, et al. clusterProfiler: an R package for comparing biological themes among gene clusters[J]. OMICS, 2012, 16(5): 284-287. DOI:10.1089/omi.2011.0118
[14]
Bindea G, Mlecnik B, Tosolini M, et al. Spatiotemporal dynamics of intratumoral immune cells reveal the immune landscape in human cancer[J]. Immunity, 2013, 39(4): 782-795. DOI:10.1016/j.immuni.2013.10.003
[15]
Yoshihara K, Shahmoradgoli M, Martínez E, et al. Inferring tumour purity and stromal and immune cell admixture from expression data[J]. Nat Commun, 2013, 4: 2612. DOI:10.1038/ncomms3612
[16]
Yang W, Soares J, Greninger P, et al. Genomics of Drug Sensitivity in Cancer (GDSC): a resource for therapeutic biomarker discovery in cancer cells[J]. Nucleic Acids Res, 2013, 41(Database issue): D955-D961.
[17]
Maeser D, Gruener RF, Huang RS. oncoPredict: an R package for predicting in vivo or cancer patient drug response and biomarkers from cell line screening data[J]. Brief Bioinform, 2021, 22(6): bbab260. DOI:10.1093/bib/bbab260
[18]
Szklarczyk D, Gable AL, Lyon D, et al. STRING v11: protein-protein association networks with increased coverage, supporting functional discovery in genome-wide experimental datasets[J]. Nucleic Acids Res, 2019, 47(D1): D607-D613. DOI:10.1093/nar/gky1131
[19]
Shannon P, Markiel A, Ozier O, et al. Cytoscape: a software environment for integrated models of biomolecular interaction networks[J]. Genome Res, 2003, 13(11): 2498-2504. DOI:10.1101/gr.1239303
[20]
Chin CH, Chen SH, Wu HH, et al. cytoHubba: identifying hub objects and sub-networks from complex interactome[J]. BMC Syst Biol, 2014, 8(4): S11.
[21]
Robin X, Turck N, Hainard A, et al. pROC: an open-source package for R and S+ to analyze and compare ROC curves[J]. BMC Bioinform, 2011, 12(1): 77. DOI:10.1186/1471-2105-12-77
[22]
Reed VL, Lalu E, Yoon L, et al. Uncovering a novel role of ROR1 in the epigenetic regulation of tumor suppressor gene CREB3L1 in triple-negative breast cancer cells[J]. Biomolecules, 2025, 15(5): 734. DOI:10.3390/biom15050734
[23]
Pan ZF, Xu T, Bao LS, et al. CREB3L1 promotes tumor growth and metastasis of anaplastic thyroid carcinoma by remodeling the tumor microenvironment[J]. Mol Cancer, 2022, 21(1): 190. DOI:10.1186/s12943-022-01658-x
[24]
Tsai YF, Lai JI, Liu CY, et al. Correlation between PIK3R1 expression and cell growth in human breast cancer cell line BT-474 and clinical outcomes[J]. World J Oncol, 2025, 16(1): 131-141. DOI:10.14740/wjon1986
[25]
Mukherjee A, Kakati RT, Van Alsten S, et al. DAB2IP loss in luminal a breast cancer leads to NF-κB-associated aggressive oncogenic phenotypes[J]. JCI Insight, 2024, 9(23): e171705. DOI:10.1172/jci.insight.171705
[26]
Xiong ZC, Yang L, Li N, et al. DAB2IP attenuates chemoresistance of triple-negative breast cancer through sequestration of RAC1 to prevent β-catenin nuclear accumulation[J]. Clin Transl Med, 2022, 12(12): e1133. DOI:10.1002/ctm2.1133
[27]
Marciniak SJ, Chambers JE, Ron D. Pharmacological targeting of endoplasmic reticulum stress in disease[J]. Nat Rev Drug Discov, 2022, 21(2): 115-140.
[28]
Zhang LL, Li SN, Shi JJ, et al. The LncRNA RMST-miR-4295-ITPR1 axis: a key mechanism in regulating autophagy in triple-negative breast cancer cells[J]. BMC Cancer, 2025, 25(1): 782.
[29]
Shang FJ, Nie HF, Du LY, et al. Inhibition of ATG5-mediated autophagy maintains PMAIP1 stability to promote cell apoptosis and suppress triple-negative breast cancer progression[J]. Discov Oncol, 2025, 16(1): 687.
[30]
Perez Kerkvliet C, Dwyer AR, Diep CH, et al. Glucocorticoid receptors are required effectors of TGFβ1-induced p38 MAPK signaling to advanced cancer phenotypes in triple-negative breast cancer[J]. Breast Cancer Res, 2020, 22(1): 39.
[31]
Gou X, Anurag M, Lei JT, et al. Transcriptional reprogramming differentiates active from inactive ESR1 fusions in endocrine therapy-refractory metastatic breast cancer[J]. Cancer Res, 2021, 81(24): 6259-6272.
[32]
Deng H, Muthupalani S, Erdman S, et al. Translocation of Helicobacter hepaticus synergizes with myeloid-derived suppressor cells and contributes to breast carcinogenesis[J]. OncoImmunology, 2022, 11(1): 2057399.
[33]
Pang BR, Wu H. Metabolic reprogramming in colorectal cancer: a review of aerobic glycolysis and its therapeutic implications for targeted treatment strategies[J]. Cell Death Discov, 2025, 11: 321.
[34]
Huggins RJ, Greene GL. ERα/PR crosstalk is altered in the context of the ERα Y537S mutation and contributes to endocrine therapy-resistant tumor proliferation[J]. NPJ Breast Cancer, 2023, 9(1): 96.
[35]
Tan CC, Lyu H, Ruan SB, et al. Subtype-specific HER3 enrichment in basal-like breast cancer is regulated via the GATA2/GATA3-FOXA1 axis[J]. Cancer Lett, 2025, 633: 218001.
[36]
Kong H, Pan M, Sun L. Identification and analysis of metabolic reprogramming-related genes in triple-negative breast cancer[J]. Clin Exp Med, 2025, 25(1): 332.
[37]
Dwyer AR, Truong TH, Kerkvliet CP, et al. Insulin receptor substrate-1 (IRS-1) mediates progesterone receptor-driven stemness and endocrine resistance in oestrogen receptor+ breast cancer[J]. Br J Cancer, 2021, 124(1): 217-227.
[38]
曲龙嘉, 郑亚迪, 罗姿麟, 等. 2000—2019年中国早发型癌症疾病负担及危险因素变化趋势[J]. 实用肿瘤杂志, 2024, 39(2): 97-103. DOI:10.13267/j.cnki.syzlzz.2024.017
[39]
高天, 李照廷, 肖敏. 免疫检查点抑制剂在激素受体阳性乳腺癌新辅助治疗中的研究进展[J]. 实用肿瘤杂志, 2025, 40(1): 1-6. DOI:10.13267/j.cnki.syzlzz.2025.001
[40]
Duan S, Li X, Song C, et al. Isoliquiritigenin inhibits triple-negative breast cancer progression via targeting the IRF5/SLC7A5/IDO1-mediated tryptophan metabolism pathway[J]. Oncol Res, 2025, 33(11): 3543-3556.
[41]
Shi A, An Y, Yun F, et al. IDO1, IL4I1: Novel immune checkpoints in breast cancer tumor-associated macrophages[J]. Breast Cancer, 2025, 17: 1279-1292.