基于单细胞测序与遗传因果推断鉴定结直肠癌致病基因和治疗靶点

Identification of causal pathogenic genes and therapeutic targets in colorectal cancer by single-cell sequencing and genetic causal inference

摘要

目的 整合单细胞RNA测序（scRNA-seq）与遗传因果推断方法系统阐明结直肠癌（CRC）的转录组学特征和遗传致病机制，挖掘其潜在治疗靶点。 方法 利用GEO数据库GSE221575数据集对CRC肿瘤组织及其配对癌旁正常组织的scRNA-seq数据进行分析，解析肿瘤微环境的细胞图谱及转录状态；筛选肿瘤上皮细胞中的差异表达基因（DEGs），进而整合CRC全基因组关联研究（GWAS）与结直肠组织表达数量性状位点（eQTL）数据系统推断基因表达与CRC发病风险之间的因果关联。通过孟德尔随机化（MR）与共定位分析筛选具有因果效应的候选致病基因，结合拟时序轨迹推断、细胞间通讯分析及蛋白质-蛋白质互作（PPI）网络构建系统解析候选基因的生物学功能与药物靶点关联。 结果 scRNA-seq分析共鉴定出CRC肿瘤微环境中7类主要细胞亚群，证实其存在显著的转录异质性与细胞亚群组成差异；经两样本MR分析，初步筛选出27个与CRC发病风险存在因果关联的基因（校正后P0.05，逆方差加权法OR值95%置信区间不跨越无效线1）；进一步通过共定位分析，最终锁定7个高置信度致病基因（共享因果变异后验概率PPH40.8）：CNBP、C4orf3、CDH1、COX14、DARS、LIMA1和TMEM258。拟时序分析提示，上述核心基因主要参与上皮细胞的恶性演化与状态跃迁过程；细胞间通讯分析发现，内皮细胞是CRC肿瘤微环境中的核心信号枢纽，其与其他细胞亚群的信号互作数量最多、作用强度最高；PPI网络分析显示，EGFR为整个互作网络的核心节点，可直接连接候选致病基因与CRC临床常用药物的作用靶点。 结论 本研究通过整合单细胞转录组与遗传因果推断鉴定的7个CRC核心致病基因揭示了上皮细胞恶性演化与肿瘤微环境调控的关键机制；构建的EGFR核心靶点网络，为CRC遗传致病机制阐释与治疗靶点开发提供了坚实依据支撑。

Abstract

Objective To systematically elucidate the transcriptomic characteristics and genetic pathogenic mechanisms of colorectal cancer（CRC） and identify potential therapeutic targets by integrating single-cell RNA sequencing（scRNA-seq） with genetic causal inference methods. Methods The scRNA-seq data of CRC tumor tissues and paired adjacent normal tissues from the GEO dataset GSE221575 were analyzed to characterize the cellular landscape and transcriptional states of the tumor microenvironment. Differentially expressed genes（DEGs） in tumor epithelial cells were screened， and CRC genome-wide association study（GWAS） data were integrated with colorectal tissue expression quantitative trait loci（eQTL） data to systematically infer causal associations between gene expression and CRC risk. Mendelian randomization（MR） and colocalization analysis were employed to identify candidate causal genes， and pseudotime trajectory inference， cell-cell communication analysis， and protein-protein interaction（PPI） network construction were integrated to systematically elucidate the biological functions of candidate genes and their associations with drug targets. Results scRNA-seq analysis identified seven major cell subpopulations in the tumor microenvironment of CRC， suggesting significant trans⁃criptional heterogeneity and compositional differences. Two-sample MR analysis preliminarily identified 27 genes with causal associations to CRC risk（P_adjusted0.05， inverse variance weighted OR 95% CI not crossing the null value of 1）. Subsequent colocalization analysis further pinpointed seven high-confidence pathogenic genes（posterior probability of shared causal variants PPH40.8）： CNBP， C4orf3， CDH1， COX14， DARS， LIMA1， and TMEM258. Pseudotime analysis results showed that these core genes were mainly involved in malignant transformation and state transitions of epithelial cells. Cell-cell communication analysis results revealed that endothelial cells， the central signaling hub in the tumor microenvironment of CRC， exhibited the highest number and strength of interactions with other cell subpopulations. PPI network analysis result indicated that EGFR served as the core node， directly connecting candidate pathogenic genes with the targets of clinically used CRC drugs. Conclusion Seven core pathogenic genes for CRC identified by single-cell transcriptomics integrated with genetic causal inference can uncover the key mechanisms of malignant evolution of epithelial cells and tumor microenvironment regulation. The established EGFR-centered target network provides a solid foundation for elucidating the genetic pathogenesis of CRC and developing its therapeutic targets.

Graphical abstract

关键词

结直肠癌 / 单细胞RNA测序 / 孟德尔随机化 / 肿瘤微环境 / 治疗靶点

Key words

colorectal cancer / single-cell RNA sequencing / Mendelian randomization / tumor microenvironment / therapeutic target

结直肠癌（colorectal cancer， CRC）是全球发病率和死亡率双高的消化道恶性肿瘤，年新发病例约190万例、死亡超90万例，疾病负担持续攀升^［1，2］。CRC受人口老龄化、生活方式改变及肥胖相关代谢性疾病流行等因素影响，未来其全球负担将进一步加重^［3］。CRC的恶性进展与肿瘤细胞内在转录异质性、肿瘤微环境（tumor microenvironment， TME）免疫抑制重塑、上皮-间质转化（epithelial-mesenchymal transition， EMT）等核心机制密切相关，但其关键驱动机制尚未被系统阐明。当前CRC筛查、手术及系统治疗已取得诸多进展，但临床靶向治疗仍集中于EGFR、VEGF等少数靶点，以西妥昔单抗、阿柏西普及氟尿嘧啶类化疗药物为核心的晚期治疗方案，仅能使部分患者临床获益，且继发性耐药频发，患者总体预后不佳，复发转移与耐药问题突出^［4］。核心机制解析不足、精准诊疗靶点匮乏，仍是当前CRC研究领域的核心瓶颈。

在基础研究层面，传统批量转录组分析已广泛用于肿瘤机制研究，但其仅能获得大量细胞的平均表达信号，无法解析肿瘤组织内广泛存在的细胞间异质性，而该异质性正是调控肿瘤侵袭能力与治疗耐药性的核心因素，存在显著技术局限^［5，6］。单细胞RNA测序（single-cell RNA sequencing， scRNA-seq）技术则突破了这一瓶颈，实现了单细胞分辨率下肿瘤微环境细胞组成与功能状态的精准解析，完成了从“组织平均水平”到“单细胞水平”的技术跨越，为揭示CRC肿瘤异质性、识别恶性细胞特征性转录改变提供了技术手段^［7］。但其本质上仍属于观察性研究，仅能揭示基因表达与肿瘤表型的相关性，无法明确基因与CRC发病风险的因果关联，筛选出的候选基因往往缺乏坚实的遗传学因果证据支撑。

在群体遗传学研究领域，全基因组关联研究（genome-wide association studies， GWAS）已识别数百个CRC发病相关遗传位点，为解析CRC遗传致病机制提供了群体遗传学基础^［8］。表达数量性状位点（expression quantitative trait locus， eQTL）可揭示特定遗传变异对基因表达的调控效应，搭建起遗传位点与基因功能间的关联桥梁^［9］。孟德尔随机化（Mendelian randomization， MR）分析借助GWAS数据，以遗传变异为工具变量，可以从遗传学层面推断暴露因素与疾病结局间的因果关联，有效规避混杂偏倚与反向因果问题^［10］。MR联合共定位分析（colocalization analysis）更是为肿瘤致病基因的筛选提供了有力的方法学支撑^［11］。然而，现有遗传分析体系仍存在关键瓶颈：GWAS识别的绝大多数CRC发病相关遗传位点位于基因组非编码区，靶基因与细胞类型特异性的功能调控机制尚不明确^［8］。同时，传统遗传因果推断多基于组织批量水平数据，无法将遗传风险信号精准锚定到CRC发生发展的核心效应细胞群中，难以实现精细细胞类型层面的因果推断，这一核心研究空白，严重制约了CRC驱动机制的解析与诊疗靶点的挖掘。

因此，本研究拟通过整合scRNA-seq与群体遗传学多组学数据，突破传统研究无法兼顾单细胞异质性解析与遗传因果推断的技术瓶颈，筛选出与CRC发病风险具有显著因果效应的关键驱动基因，系统探究候选基因在CRC恶性进展中的生物学功能与核心调控通路，最终为CRC的精准诊疗提供新的理论依据与潜在干预靶点。

1 材料与方法

1.1 分析流程与数据来源

本研究整体分析流程如下：①对scRNA-seq原始数据进行质量控制、标准化与细胞分群分析，识别具有不同转录特征的细胞亚群；②基于细胞经典标志基因完成细胞类型注释，通过通路富集分析解析各细胞群的功能特征，重点关注在肿瘤发生发展中发挥核心驱动作用的上皮细胞亚群；③筛选肿瘤上皮细胞与正常上皮细胞间的差异表达基因（diffe⁃rentially expressed genes， DEGs），将其与CRC GWAS汇总数据、结直肠组织eQTL数据进行整合分析；④采用MR联合共定位分析，筛选同时获得遗传学与转录组学证据支持的CRC潜在致病基因；⑤通过拟时序轨迹分析、细胞间通讯网络分析、蛋白质-蛋白质相互作用（protein-protein interaction， PPI）网络构建等下游分析，解析候选基因的功能背景与作用机制。

本研究单细胞转录组数据来源于基因表达综合数据库（Gene Expression Omnibus， GEO，https：//www.ncbi.nlm.nih.gov/geo/），数据集编号为GSE221575，该数据集包含2例CRC组织及2例配对癌旁正常组织的scRNA-seq数据。CRC相关GWAS汇总数据来源于FinnGen R12队列，该队列共纳入11 790例CRC患者与378 749例健康对照人群。CRC相关eQTL数据则来源于表达数量性状位点基因组学联盟（expression Quantitative Trait Loci Genomics Consortium， eQTLGen Consortium，https：//www.eqtlgen.org/），该数据集为解析遗传变异对基因表达的调控关系提供了重要依据。

1.2 数据预处理与质量控制

原始测序数据（FASTQ文件）采用Cell Ranger软件（版本6.0）进行处理，测序序列比对至GRCh38版本人类参考基因组序列，通过唯一分子标识符（unique molecular identifier， UMI）计算每个细胞的基因表达矩阵。所有下游分析均在R软件（版本4.2.2）中完成，主要采用Seurat包（版本4.0.5），该软件包可实现单细胞数据分析中的标准化、缩放与聚类分析。本研究选取了4组代表性样本，分别为：直肠癌癌旁正常组织（Rectum-03N）、直肠癌肿瘤组织（Rectum-03T）、乙状结肠癌癌旁正常组织（Sigmoid-03N）、乙状结肠癌肿瘤组织（Sigmoid-03T）。从数据库下载数据后，采用严格过滤阈值对样本细胞进行质量控制（quality control， QC），细胞保留标准需同时满足以下3项：单个细胞检测基因数（nFeature_RNA）6 000；唯一分子标识符计数（nCount_RNA）1 000；线粒体基因转录本比例（percent.mt）20%。该标准可有效去除受环境RNA污染与完整性受损的细胞，确保后续分析基于高质量的单细胞数据^［12］。为评估上述质量控制流程的有效性，采用小提琴图可视化4组代表性样本的nFeature_RNA、nCount_RNA与percent.mt分布，过滤低质量细胞后进行下游分析。

1.3 降维与聚类分析

对经过质量过滤与标准化的单细胞数据进行对数转换与缩放处理，采用Seurat包中的FindVariableFeatures函数识别高可变基因（highly variable genes， HVGs）。基于筛选得到的HVGs进行主成分分析（principal component analysis， PCA）降维，通过方差贡献率评估与碎石图肘部法则评估，选取前10个主成分（principal components， PCs）用于下游分析。采用均匀流形逼近与投影（uniform manifold approximation and projection， UMAP）将高维数据投射至二维空间进行可视化。随后通过FindClusters函数，基于Louvain算法对降维后的细胞进行图聚类分析，分辨率设置为0.5；该参数设置可有效将细胞划分为不同细胞亚群，真实反映细胞内在的转录程序特征。

1.4 细胞类型注释与通路富集分析

首先基于细胞经典标志基因对细胞类型进行初步注释，采用SingleR软件包，以人类原代细胞图谱数据集（Human Primary Cell Atlas， HPCA）为参考，对细胞注释结果进行进一步验证。随后采用基因集变异分析（gene set variation analysis， GSVA），探究各细胞类型中显著富集的生物学通路^［13］。该分析使用分子特征数据库（Molecular Signatures Database， MSigDB，版本7.5.1）中的Hallmark基因集，基于GSVA R包完成。

1.5 差异表达基因筛选与孟德尔随机化分析

基于Wilcoxon秩和检验，采用Seurat软件中的FindMarkers函数进行差异基因表达分析。以校正后P0.05为标准筛选DEGs，从而识别上皮细胞恶性转化过程中转录水平发生显著改变的基因。将具有显著eQTL的基因与CRC GWAS汇总统计数据进行整合，探究基因表达的遗传调控与CRC发病的因果关联。研究采用R软件TwoSampleMR包进行两样本MR分析，主要采用逆方差加权法（inverse variance weighted， IVW）进行MR因果推断，结果经错误发现率（false discovery rate， FDR）校正，以校正后P0.05且IVW法比值比（odds ratio， OR）值的95%置信区间不跨越无效线1为筛选标准，最终获得与CRC发病风险存在显著因果关联的基因^［14］。为进一步强化因果证据的可靠性，采用R软件coloc包进行共定位分析，评估同一基因组位点内，调控基因表达的遗传变异与CRC发病风险相关的遗传变异是否共享同一因果变异；以基因表达与CRC发病共享同一因果变异的后验概率PPH40.8作为强共定位证据的“金标准”，最终锁定高置信度CRC因果致病基因^［15］。

1.6 拟时序轨迹推断与细胞间通讯分析

为探究CRC肿瘤微环境的细胞发育轨迹与通讯模式，采用Monocle 3软件对肿瘤上皮细胞进行拟时序分析，以经前期筛选的高可变基因构建细胞发育轨迹主图^［16］；同时采用CellChat软件分析细胞间配体-受体介导的信号相互作用，通过定量分析不同细胞群间相互作用的数量与强度，识别TME中的核心信号枢纽与关键调控通路，分析结果采用热图与网络图进行可视化呈现^［17］。

1.7 蛋白质-蛋白质相互作用网络的构建

为明确结直肠癌关键调控蛋白与治疗靶点间的关联，采用STRING数据库（版本11.5），构建候选致病基因与CRC临床常用药物（西妥昔单抗、阿柏西普、氟尿嘧啶、亚叶酸钙）已知靶点间的PPI网络^［18］。设置PPI互作的综合置信度（范围0~1）最低阈值为0.7，剔除网络中的孤立节点；采用Cytoscape 3.9.1软件对网络拓扑结构进行可视化，计算每个节点的度值并识别网络核心枢纽蛋白^［19］。该分析可将本研究新发现的候选致病基因映射至已知药物靶点的互作网络中。

2 结果

2.1 单细胞转录组数据的质量控制结果

对研究纳入的4组样本（Rectum-03N、Rectum-03T、Sigmoid-03N与Sigmoid-03T）进行质量控制后，采用小提琴图对各组细胞的质量控制指标分布特征进行可视化分析（见图1）。结果显示，直肠来源样本（含肿瘤及癌旁正常组织）nFeature_RNA与nCount_RNA分布更集中，提示其转录本捕获稳定性更佳，技术异质性更低；乙状结肠来源样本nFeature_RNA与nCount_RNA分布范围更广，部分细胞UMI计数超50 000，提示该组样本生物学异质性更高、细胞整体转录活性更强。4组样本的mt_percent整体分布趋势一致，乙状结肠来源样本mt_percent异质性更显著，提示不同解剖部位肿瘤微环境的细胞应激状态存在差异。本研究依据上述质控指标特征，设定mt_percent20%、nCount_RNA1 000、nFeature_RNA6 000的三重过滤标准，有效剔除低质量细胞与离群细胞，保障后续分析数据的可靠性。

2.2 单细胞转录组数据的降维与聚类分析结果

基于前期筛选得到的高可变基因，对质控后的数据进行主成分分析降维处理。PCA肘部图结果显示，前10个主成分后的方差解释率出现显著下降，且这10个主成分的累计方差贡献率70%（见图2A），其可作为有效主成分用于后续的降维分析。UMAP可视化结果显示，细胞在二维空间中形成界限清晰的聚类群体，共得到18个细胞亚群（聚类编号0~17，见图2B）。部分细胞亚群在UMAP空间中形成致密且分界明确的聚类，而肿瘤上皮细胞（聚类编号0、4、17）与部分免疫细胞亚群（6、13）间存在局部区域的重叠。按样本来源对UMAP投影结果进行着色分析显示，直肠与乙状结肠来源的细胞呈部分区分但相互重叠的分布模式，其中直肠来源样本的细胞群分布相对紧密，而乙状结肠来源样本的细胞群分布更为分散，该特征与质量控制分析中观察到的乙状结肠来源样本异质性更高的结果相一致；肿瘤样本与癌旁正常样本的细胞在UMAP空间中呈现明显的分离趋势（见图2C）。

2.3 单细胞转录组数据的细胞类型注释与功能富集分析结果

基于经典标记基因并结合参考数据库，对质控、降维聚类后的结直肠癌单细胞转录组数据进行系统的细胞类型注释，共鉴定出7类主要细胞亚群：上皮细胞、CD8⁺ T细胞、B细胞、成纤维细胞、单核细胞、造血干细胞及内皮细胞。其中上皮细胞在UMAP空间中呈现出最为集中的空间分布特征（见图3A）。对比肿瘤组织与癌旁正常组织的细胞分布特征发现，上皮细胞与免疫细胞可清晰区分为恶性与非恶性细胞群体，但部分细胞聚类的边界区域存在群体重叠现象（见图3B）。

采用GSVA对结直肠癌肿瘤微环境各细胞亚群的生物学通路富集特征进行分析，并构建功能富集热图。结果表明，不同细胞类型呈现出高度特异的通路活化模式，肿瘤微环境内部存在显著的功能异质性（见图3C）。肿瘤来源的上皮细胞显著富集于经典肿瘤相关信号通路，包括KRAS信号通路、上皮-间质转化（EMT）通路及TNFα/NF-κB信号通路；CD8⁺T细胞、B细胞等免疫细胞富集于免疫应答相关通路；成纤维细胞则富集于细胞外基质重塑相关通路。

2.4 CRC上皮细胞中潜在因果基因的鉴定与表达分析结果

研究首先对肿瘤组织与癌旁正常组织来源的上皮细胞进行差异表达分析，以校正后P0.05为阈值，筛选获得多个携带显著eQTL调控信号的基因。将上述eQTL关联基因与CRC GWAS汇总统计数据整合，经两样本MR分析，最终得到27个与CRC发病风险存在统计学显著关联的候选因果基因。27个候选基因的MR分析结果显示，其95%CI均未跨越OR=1的无效线，结合FDR校正后P0.05（见图4A）。进一步对上述基因进行共定位分析，最终锁定7个高置信度核心致病基因，包括CNBP、C4orf3、CDH1、COX14、DARS、LIMA1及TMEM258，其共享因果变异的后验概率均PPH40.8。点阵图可视化结果显示，7个核心基因在CRC肿瘤微环境各细胞亚群中的表达分布存在显著的细胞特异性差异，其中CNBP、CDH1、TMEM258在肿瘤上皮细胞中呈高表达，在免疫细胞与基质细胞中亦可检测到一定水平的表达（见图4B）。小提琴图验证结果显示，CNBP在肿瘤上皮细胞中表达上调幅度最显著；CDH1、LIMA1、DARS在肿瘤上皮细胞中均呈显著表达上调，变化趋势一致；COX14等其余核心基因表达上调幅度相对温和，整体呈一致上调趋势（见图4C）。

2.5 CRC细胞发育轨迹与细胞间通讯特征分析结果

使用Monocle 3软件对CRC肿瘤上皮细胞进行拟时序轨迹分析，结果显示肿瘤上皮细胞沿拟时序轨迹呈现连续的状态转换模式，而非离散的群体划分（见图5A）。进一步在同一UMAP二维嵌入空间及拟时序轨迹框架下对基因表达动态进行分析，筛选得到CDH1、DARS、TMEM258等分支特异性调控基因，该类基因在轨迹的特定分支中呈显著上调表达（见图5B）。

基于CellChat软件对CRC肿瘤微环境中所有细胞亚群进行细胞间通讯分析，构建全细胞亚群配体-受体介导的信号互作网络，结果显示上皮细胞、免疫细胞与基质细胞间存在广泛的信号通讯关联（见图5C）；其中内皮细胞与其他细胞亚群间的信号互作数量最多、作用强度最高，成纤维细胞与单核细胞同样为细胞间通讯网络的重要组成部分。

2.6 CRC候选致病基因与常用药物靶点PPI网络分析结果

使用STRING 11.5数据库，以0.7为最低相互作用置信度阈值，构建CRC临床常用药物靶点与MR联合共定位筛选的候选致病基因的PPI网络（见图6）。网络纳入节点包括：西妥昔单抗靶基因EGFR，阿柏西普靶基因VEGFB、PGF，氟尿嘧啶/亚叶酸钙靶基因TYMS、DHFR、SLC46A1、DPYD，以及候选致病基因CNBP、C4orf3、CDH1、COX14、DARS、LIMA1、TMEM258。PPI网络分析结果显示，EGFR在相互作用网络中处于核心枢纽地位，其与Fcγ受体家族成员（FCGR1A、FCGR2A、FCGR3A、FCGR3B）、补体成分（C1QA、C1QB、C1QC）均存在直接相互作用，同时与CDH1、LIMA1、DARS、TMEM258等候选致病基因存在直接分子相互作用。

阿柏西普的作用靶点VEGFB和PGF在PPI网络中形成局部功能作用簇，该簇与EGFR存在直接分子连接，且与血管生成相关功能模块紧密关联。氟尿嘧啶/亚叶酸钙的作用靶点在PPI网络中呈相对独立分布，其中TYMS、DHFR与EGFR存在直接分子相互作用。CNBP、C4orf3等候选致病基因位于PPI网络的边缘区域；CDH1、LIMA1在CRC肿瘤上皮细胞中呈显著富集表达，且与以EGFR为核心的药物靶点功能模块紧密连接。

3 讨论

近年来，单细胞测序技术的普及极大推动了CRC研究范式的革新，该技术可从高分辨率视角解析肿瘤异质性、微环境调控特征及潜在分子机制，为CRC的发病机制研究提供了全新视角^［20］。但目前多数研究仍局限于对肿瘤微环境的单细胞表型解析，鲜有研究将单细胞水平的基因表达特征与群体遗传学的因果推断相结合，难以从遗传根源上明确CRC发病的关键驱动基因。针对这一研究空白，本研究创新性整合scRNA-seq、GWAS、组织特异性eQTL、两样本MR及共定位分析技术，实现了单细胞异质性解析与群体遗传学因果推断的有机结合，系统解析了CRC TME的细胞组成、转录特征及细胞间通讯网络，筛选获得7个高置信度CRC致病基因，并构建了致病基因与临床药物靶点的PPI网络，为CRC发病机制阐释与精准治疗提供了全新的理论依据。

CRC并非均一的疾病实体，既往研究已明确，直肠与乙状结肠在胚胎起源、肠道菌群组成、脂质代谢特征等方面存在固有差异，这种差异可进一步转化为肿瘤细胞转录活性与微环境细胞组成的异质性，是影响CRC疾病进展与治疗响应的重要因素^［21］。本研究首先对纳入的scRNA-seq数据完成了严格的质量控制，结果显示直肠与乙状结肠来源样本的质控结果存在系统性差异：乙状结肠来源样本的检测基因数与UMI计数分布范围更广，提示该类样本具有更高的生物学异质性，也可能在检测过程中引入了更多技术噪音。这一结果与既往解剖部位相关的CRC异质性研究结论相呼应，同时也提示，在CRC相关单细胞研究中，必须将解剖部位作为关键混杂因素纳入分析体系，方可有效区分技术噪音与真实生物学异质性，精准揭示部位特异性的肿瘤发生发展机制。

既往研究已充分证实，CRC肿瘤微环境中，恶性上皮细胞是驱动肿瘤自主恶性增殖、侵袭转移的核心细胞，而免疫细胞与基质细胞则通过激活免疫应答、重塑细胞外基质等方式，共同塑造免疫抑制、促肿瘤生长的微环境^{［22，23］}。与之相符的是，本研究的UMAP可视化结果清晰呈现了肿瘤细胞与正常细胞的显著分离趋势，同时也印证了直肠与乙状结肠来源细胞的分布模式差异；基于经典标志基因，本研究完成了CRC肿瘤微环境7类主要细胞类型的注释，GSVA功能富集分析进一步显示，恶性上皮细胞显著富集KRAS、EMT等经典肿瘤相关通路，免疫细胞富集免疫应答相关通路，成纤维细胞则富集细胞外基质重塑相关通路，明确了不同细胞类型在肿瘤进展中的功能分工。在此基础上，本研究进一步发现，CNBP、CDH1、TMEM258在肿瘤上皮细胞中呈高表达特征，且在免疫细胞与基质细胞中均可检测到稳定的表达信号，提示这类基因可能广泛参与微环境的细胞间互作调控^［24］。同时本研究观察到，上皮细胞与部分免疫细胞群落在UMAP空间存在局部区域重叠，这一现象反映了肿瘤微环境中普遍存在的细胞可塑性，而既往研究也指出，这种连续的细胞状态特征，正是单一细胞群靶向治疗疗效受限的核心分子机制^［25］。

既往GWAS研究已识别出大量与CRC发病风险相关的遗传位点，但绝大多数位点位于非编码区，无法直接明确其靶基因与致病机制，难以实现从风险位点到功能基因的转化。为从遗传学层面精准筛选CRC发生的关键驱动基因，本研究整合单细胞转录组分析、GWAS、组织特异性eQTL与MR分析技术，从表型相关性中精准推断因果关联，识别出在肿瘤上皮细胞中异常表达、且与CRC遗传发病风险存在因果关联的27个候选基因；经共定位分析进一步验证，最终明确了CNBP、C4orf3、CDH1、COX14、DARS、LIMA1、TMEM258共7个高置信度致病基因。其中，CNBP在恶性上皮细胞中表达显著上调，且在微环境多类细胞中呈广泛表达模式，既往研究显示，CNBP可通过调控RNA稳定性与翻译过程，同时参与肿瘤细胞固有生物学行为调控与微环境细胞间通讯^{［26，27］}，这也为我们观察到的多细胞表达特征提供了功能解释。CDH1编码的E-钙黏蛋白是细胞间黏附的核心分子，经典EMT理论中其在肿瘤进展中呈表达下调趋势，但本研究观察到CDH1在肿瘤上皮细胞中表达显著上调，这一现象与部分CRC相关研究的结果一致，既往研究提示，这种异常上调可能与肿瘤早期细胞团黏附维持、亚群特异性表型或蛋白翻译后修饰相关，并非完全遵循经典EMT的表达模式^［28］。此外，本研究发现TMEM258、DARS等功能尚未被深入研究的基因，在微环境多细胞群中广泛表达，这为CRC分子机制研究提供了高优先级的候选靶点。

肿瘤细胞的表型可塑性是肿瘤适应微环境应激、产生药物耐药的核心机制之一，既往研究已证实，CRC上皮细胞在分化与去分化状态间的动态转换，是肿瘤内部异质性的重要来源^［29］。本研究拟时序轨迹分析结果显示，CRC上皮细胞沿连续轨迹发生动态演化，存在明显的状态转换与命运分支，直观揭示了肿瘤上皮细胞的可塑性特征。在此基础上，本研究进一步筛选得到CDH1、DARS、TMEM258等分支特异性调控基因，这类基因在轨迹的不同分支中呈显著上调表达，提示其可作为调控细胞命运决定的“分子开关”，介导CRC上皮细胞的谱系分叉，进而加剧肿瘤内部异质性，最终影响疾病进展与治疗响应^［30］。其中TMEM258作为内质网跨膜蛋白，既往研究已报道其主要参与内质网应激反应，而内质网应激与肿瘤细胞存活、增殖、化疗耐药及免疫逃逸密切相关，这也提示，靶向TMEM258及相关通路，有望为克服CRC药物耐药提供全新的研究方向^{［31，32］}。

CRC肿瘤微环境中，细胞间的异常信号通讯是驱动肿瘤进展、免疫逃逸与治疗耐药的关键环节。既往研究已证实，血管内皮细胞、成纤维细胞与单核/巨噬细胞是微环境细胞间通讯的核心参与者，其中血管微环境是协调肿瘤内部各类信号交流的关键平台^［33］。而成纤维细胞与单核/巨噬细胞可通过分泌大量信号分子，参与免疫抑制与基质纤维化微环境的构建，为肿瘤进展提供支撑^［34］。本研究基于CellChat的细胞通讯分析，完整描绘了CRC微环境中复杂的细胞间通讯网络，结果显示内皮细胞是介导信号传递的核心枢纽，其与其他细胞类型的互作数量最多、强度最高，成纤维细胞与单核细胞同样是通讯网络的重要组成部分，这一结果与既往研究的结论高度一致。在此基础上，本研究进一步整合CRC临床常用药物靶点与筛选获得的核心致病基因构建PPI网络，结果显示EGFR处于互作网络的核心位置，其不仅与免疫相关蛋白存在显著互作，还与CDH1、LIMA1等多个核心致病基因存在直接分子连接，是西妥昔单抗、阿柏西普、氟尿嘧啶等不同作用机制临床药物的信号整合枢纽，为临床联合用药方案的优化提供了分子机制解释。这一结果也提示，CRC治疗策略需从“靶向单个分子”向“干预整个信号网络”转变，本研究筛选的核心基因兼具作为疗效预测生物标志物与新型治疗靶点的潜力。

然而，本研究仍存在一定局限性：单细胞测序仅纳入2例患者的4例肿瘤组织及癌旁正常组织样本，样本量有限，虽通过自身配对设计和大样本遗传因果推断弥补了部分不足，但结论的外推性仍需更大样本独立验证；本研究基于公共数据库测序数据开展，尽管完成了严格的质量控制，样本来源与临床信息的异质性仍可能对结果产生影响^［35］；筛选得到的候选基因功能与分子调控机制，仍需体内外功能实验进一步验证；对免疫细胞、基质细胞的亚型细分与功能异质性解析仍有待深入。后续可通过类器官模型与基因编辑技术开展功能验证，借助空间转录组技术实现基因表达的原位验证，在前瞻性临床队列中完成标志物与联合用药方案的转化探索。

4 结论

本研究通过整合scRNA-seq、GWAS、组织特异性eQTL、MR及共定位分析等多组学技术，构建了从CRC遗传风险、单细胞功能解析到药物靶点挖掘的多维度研究框架，系统揭示了CRC TME的解剖部位特异性异质性特征，成功鉴定出7个与CRC发病存在因果关联的高置信度致病基因（CNBP、C4orf3、CDH1、COX14、DARS、LIMA1、TMEM258），明确了其在上皮细胞恶性演化中的作用，并构建了以EGFR为核心的药物靶点互作网络。本研究的核心发现不仅深化了对CRC遗传致病机制与肿瘤微环境调控规律的认知，更为CRC的精准诊疗提供了全新的生物标志物与潜在治疗靶点，为后续的机制研究与临床转化提供了坚实的理论依据。

参考文献

[1] Bray F， Laversanne M， Sung H， et al. Global cancer statistics 2022： GLOBOCAN estimates of incidence and mortality worldwide for 36 cancers in 185 countries［J］. CA Cancer J Clin， 2024， 74（3）： 229-263.

[2] Darmadi D， Mohammadian⁃Hafshejani A， Kheiri S. Global disparities in colorectal cancer： unveiling the present landscape of incidence and mortality rates， analyzing geographical variances， and assessing the human development index［J］. J Prev Med Hyg， 2024， 65（4）： E499-E514.

[3] Arnold M， Sierra MS， Laversanne M， et al. Global patterns and trends in colorectal cancer incidence and mortality［J］. Gut， 2017， 66（4）： 683-691.

[4] Van Cutsem E， Cervantes A， Adam R， et al. ESMO consensus guidelines for the management of patients with metastatic colorectal cancer［J］. Ann Oncol， 2016， 27（8）： 1386-1422.

[5] Ding S， Chen X， Shen K. Single⁃cell RNA sequencing in breast cancer： understanding tumor heterogeneity and paving roads to individualized therapy［J］. Cancer Commun， 2020， 40（8）： 329- 344.

[6] Lapuente⁃Santana Ó， Sturm G， Kant J， et al. Multimodal analysis unveils tumor microenvironment heterogeneity linked to immune activity and evasion［J］. iScience， 2024， 27（8）： 110529.

[7] Li H， Courtois ET， Sengupta D， et al. Reference component analysis of single⁃cell transcriptomes elucidates cellular heterogeneity in human colorectal tumors［J］. Nat Genet， 2017， 49（5）： 708-718.

[8] Huyghe JR， Bien SA， Harrison TA， et al. Discovery of common and rare genetic risk variants for colorectal cancer［J］. Nat Genet， 2019， 51（1）： 76-87.

[9] Ghaffar A， International Headache Genetics Consortium， Nyholt DR. Integrating eQTL and GWAS data characterises established and identifies novel migraine risk loci［J］. Hum Genet， 2023， 142（8）： 1113-1137.

[10] Wang P， Lin Z， Pan W. Unbiased causal inference with Mendelian randomization and covariate⁃adjusted GWAS data［J］. HGG Adv， 2025， 6（2）： 100412.

[11] Jiao M， Cui Y， Qiu X， et al. Causal relationship between complement C1QB and colorectal cancer： a drug target mendelian randomization study［J］. Front Genet， 2024， 15： 1403509.

[12] Park Y， Hauschild AC. The effect of data transformation on low⁃dimensional integration of single⁃cell RNA⁃seq［J］. BMC Bioinformatics， 2024， 25（1）： 171.

[13] Hänzelmann S， Castelo R， Guinney J. GSVA： gene set variation analysis for microarray and RNA⁃seq data［J］. BMC Bioinformatics， 2013， 14： 7.

[14] Zhang S， Zhang W， Ren H， et al. Mendelian randomization analysis revealed a gut microbiota⁃mammary axis in breast cancer［J］. Front Microbiol， 2023， 14： 1193725.

[15] Huang K， Huang H， Ma LL， et al. Identification of breast cancer susceptibility genes via trans⁃ethnic Mendelian randomization and colocalization analyses［J］. Medicine， 2025， 104（49）： e45751.

[16] Cao J， Spielmann M， Qiu X， et al. The single⁃cell transcriptional landscape of mammalian organogenesis［J］. Nature， 2019， 566（7745）： 496-502.

[17] Jin S， Guerrero⁃Juarez CF， Zhang L， et al. Inference and analysis of cell⁃cell communication using CellChat［J］. Nat Commun， 2021， 12（1）： 1088.

[18] Szklarczyk D， Kirsch R， Koutrouli M， et al. The STRING database in 2023： protein⁃protein association networks and functional enrichment analyses for any sequenced genome of interest［J］. Nucleic Acids Res， 2023， 51（D1）： D638-D646.

[19] Shen M， Ma Z， Cao Y， et al. SLC26A9 in triple⁃negative breast cancer stem cells： a network pharmacology and molecular mo⁃deling study［J］. Front Bioeng Biotechnol， 2025， 13： 1703343.

[20] Wen R， Zhou L， Peng Z， et al. Single⁃cell sequencing techno⁃logy in colorectal cancer： a new technology to disclose the tumor heterogeneity and target precise treatment［J］. Front Immunol， 2023， 14： 1175343.

[21] Zhu Q， Zhu C， Zhang X， et al. Comprehension of rectosigmoid junction cancer molecular features by comparison to the rectum or sigmoid colon cancer［J］. J Gastrointest Oncol， 2023， 14 （3）： 1307-1319.

[22] 马一诺，王群，巴微，等.肿瘤相关成纤维细胞在结直肠癌发生发展中的作用［J］. 生命的化学，2025， 45 （1）： 104-110.

[23] Geng Y， Feng J， Huang H， et al. Single⁃cell transcriptome analysis of tumor immune microenvironment characteristics in colorectal cancer liver metastasis［J］. Ann Transl Med， 2022， 10（21）： 1171-1170.

[24] Wu Y， Yu Y， Li D， et al. CDH1 genetic variants and its aberrant expression are the risk factors for colorectal cancer metastasis［J］. BMC Gastroenterol， 2025， 25（1）： 214.

[25] Moorman A， Benitez EK， Cambulli F， et al. Progressive plasticity during colorectal cancer metastasis［J］. Nature， 2025， 637（8047）： 947-954.

[26] Lee E， Lee TA， Yoo HJ， et al. CNBP controls tumor cell biology by regulating tumor⁃promoting gene expression［J］. Mol Carcinog， 2019， 58（8）： 1492-1501.

[27] David AP， Pipier A， Pascutti F， et al. CNBP controls transcription by unfolding DNA G⁃quadruplex structures［J］. Nucleic Acids Res， 2019， 47（15）： 7901-7913.

[28] Sicairos B， Alam S， Du Y. A comprehensive analysis of diffe⁃rent types of databases reveals that CDH1 mRNA and E⁃cadherin protein are not downregulated in most carcinoma tissues and carcinoma cell lines［J］. BMC Cancer， 2023， 23（1）： 441.

[29] White M， Mills ML， Millett LM， et al. MAPK⁃driven epithelial cell plasticity drives colorectal cancer therapeutic resistance［J］. Nature， 2026， 650（8102）： 748-758.

[30] Xiong H， Liao M， Zhang H， et al. DARS expression in BCR/ABL1⁃negative myeloproliferative neoplasms and its association with the immune microenvironment［J］. Sci Rep， 2024， 14（1）： 16711.

[31] Graham DB， Lefkovith A， Deelen P， et al. TMEM258 is a component of the oligosaccharyltransferase complex controlling ER stress and intestinal inflammation［J］. Cell Rep， 2016， 17（11）： 2955-2965.

[32] Sisinni L， Pietrafesa M， Lepore S， et al. Endoplasmic reticulum stress and unfolded protein response in breast cancer： the balance between apoptosis and autophagy and its role in drug resistance［J］. Int J Mol Sci， 2019， 20（4）： 857.

[33] Pernot S， Evrard S， Khatib AM. The give⁃and⁃take interaction between the tumor microenvironment and immune cells regulating tumor progression and repression［J］. Front Immunol， 2022， 13： 850856.

[34] Cao L， Meng X， Zhang Z， et al. Macrophage heterogeneity and its interactions with stromal cells in tumour microenvironment［J］. Cell Biosci， 2024， 14（1）： 16.

[35] Hicks SC， Townes FW， Teng M， et al. Missing data and technical variability in single⁃cell RNA⁃sequencing experiments［J］. Biostatistics， 2018， 19（4）： 562-578.

基金资助

国家自然科学基金项目(82374539)

四大慢病重大专项(2024ZD0521300/2024ZD0521304)

江苏省中医药科技发展计划重点项目(ZD202214)

江苏省中医药领军人才培养对象项目(SLJ0327)

南京中医药大学胃癌临床专病研究院项目(LCZBYJYZZ2024-001)