多组学分享--微卫星稳定型结直肠癌的多组学表征揭示独特生物学和可靶向的脆弱性

作者，Evil Genius

这世界有太多不如意

但你的生活还是要继续

太阳每天依旧要升起

希望永远种在你心里

喜欢做的事，能坚持一天是一天。

今日我们分享多组学文献

核心问题：微卫星稳定（MSS）结直肠癌（CRC）是一种分子复杂且临床具有挑战性的疾病亚型，缺乏有效的靶向治疗。

Graphical Abstract

研究的整体框架：

多组学与表型数据：包括代谢组学、临床数据、基因组学、转录组学、蛋白质组学/磷酸化蛋白质组学、表观基因组学。
单组学和成对组学分析：展示突变、拷贝数变异（CNV）、代谢组学、甲基化、转录组和蛋白质组的单组学和成对分析。
通路、网络和细胞类型分析。
数据融合、多组学和表型：通过iPhenMap机器学习方法进行。
机制见解和靶点识别：在PDX和肿瘤类器官模型中进行体内外药物反应验证，展示生存曲线。
潜在的新重用药物和未来试验组合。

In Brief

整合了来自约2,000例CRC肿瘤的多组学、单细胞和临床数据（约4,500个分子谱，包括195,246个细胞/细胞核），以定义MSS CRC的复杂性。

分析识别出复制应激、右侧性、谷氨酰胺分解和MYC-E2F相关表观遗传程序作为关键特征。
FAK、PDK1/2和蛋白酶体抑制，联合trametinib，被确定为MSS CRC的潜在治疗策略。

Highlights

大规模单组学到多组学的MSS CRC遗传和表观遗传研究。
复制应激、侧别、代谢通量和表观遗传学塑造生物学。
多组学-表型整合/映射揭示疾病复杂性。
FAK、PDK1/2和蛋白酶体成为脆弱性靶点。

摘要

微卫星稳定（MSS）结直肠癌（CRC）是一种分子复杂且临床具有挑战性的疾病亚型，缺乏有效的靶向治疗。利用来自约2,000个肿瘤样本（包括316例先前未分析的MSS CRCs）的多组学数据，代表约4,500个分子谱和195,246个单细胞/细胞核，揭示了复杂的生物学和治疗脆弱性。

应用多组学-表型整合框架（包括iPhenMap，一种机器学习工具），揭示了复杂的分子程序，包括：
- 复制应激
- 右侧生物学
- 谷氨酰胺分解驱动的代谢重编程
- MYC-E2F转录网络调控表观遗传重塑和神经周围侵犯。
聚焦RAS突变MSS CRC，在患者来源肿瘤类器官和匹配的异种移植模型中进行功能验证，确认黏着斑激酶（FAK）、磷酸肌醇依赖性激酶（PDK）1/2和蛋白酶体抑制作为与丝裂原活化蛋白激酶（MEK）阻断联合的潜在治疗策略。
这项工作强调了整合多组学分析和功能建模在完善MSS CRC患者群体潜在治疗方法方面的力量。

Introduction

1. 结直肠癌的全球疾病负担

2022年，CRC是欧洲第三大最常诊断的癌症和第二大癌症死亡原因，约有190万新病例和超过90万死亡。

2. MSS CRC的临床挑战

大多数（>85%）CRC肿瘤表现为熟练的复制后DNA错配修复（pMMR）状态，导致微卫星稳定（MSS）。
大型临床CRC数据集的事后分析表明，携带MSS CRC肿瘤的患者表现出显著更高的疾病特异性死亡率和更差的预后。
这些病例中40%-50% 携带逆转录病毒相关序列（RAS）基因超家族的突变，RAS突变（RAS-mt）患者与RAS野生型（RAS-wt）对应者相比持续显示较差结局。
虽然微卫星不稳定（MSI）CRC患者现在通过免疫疗法成功治疗，但MSS转移性CRC（mCRC）的临床试验测试化疗与新型靶向药物（包括KRAS G12C突变抑制剂）和/或免疫疗法的新组合仅显示出相对于标准治疗的适度改善。
目前，MMR、RAS和BRAF突变状态 largely 指导临床决策。
由于MSS CRC患者生存率在最近几十年没有显著改善，需要协调努力揭示疾病复杂性。
在这种情况下，多组学数据集的稳健整合分析将继续提供超越单变量研究的重要见解。

3. 现有研究局限与本研究目标

据作者所知，迄今为止没有研究在MSS CRC背景下对bulk基因组学、表观基因组学、蛋白质组学和代谢组学与单核/单细胞RNA测序（sn/scRNA-seq）和临床参数（表型）进行全面整合。
为填补这一空白并利用完全整合的"床旁到床旁到患者衍生模型"方法，多学科和多中心CRC研究网络"COLOSSUS"（www.colossusproject.eu）对316例MSS RAS-mt CRC样本进行了详尽分子表征，对应于1,070个样本分子谱，包括119,489个单细胞核。
虽然这些队列设计围绕MSS RAS-mt CRC，但作者试图确定这些发现是否反映在更广泛的MSS CRC生物学中。为此，网络使用了来自9个独立患者和肿瘤模型队列的额外1,650个肿瘤样本的数据，代表RAS-mt和RAS-wt疾病，对应于约3,450个bulk分子谱和75,757个单细胞/细胞核。个体和整合的多组学和临床协变量被融合以揭示MSS CRC分子景观的复杂性。
关键的是，无论整合的具体组学层如何，都识别出潜在的功能和治疗依赖性。基于未满足的临床需求并与COLOSSUS研究的原始设计一致，MSS RAS-mt CRC被优先用于患者衍生模型中的靶点验证研究。

Results

结果1、患者队列和多组学分析

样本收集：

在COLOSSUS项目中收集的316例MSS RAS-mt样本进行了多组学分析，包括：
- RNA测序（RNA-seq）
- 低覆盖度全基因组测序（LC-WGS）
- 全外显子测序（WES）
- 代谢组学
- 总/磷酸化蛋白质组学（液相色谱-质谱/质谱，LC-MS/MS和反相蛋白阵列，RPPA）
- snRNA-seq
来自9个独立患者和肿瘤模型队列的额外1,650个肿瘤样本（包括RAS-mt和RAS-wt CRC样本）用于扩展和验证。

结果2、MSS RAS-mt和RAS-wt CRC的共享基因组架构

突变景观：

COLOSSUS MSS RAS-mt St2/3队列（n=88）的WES，以COSMIC Cancer Gene Census（CGC）Tier 1/2注释和TCGA MSS RAS-mt CRC队列（n=216）为基准，识别出一致的复发性驱动基因景观。
20个最频繁突变的驱动基因汇聚于主要致癌通路，包括RAS、TP53、PI3K、MYC和细胞周期信号。
大多数改变是无义和错义突变。
比较TCGA MSS RAS-mt（n=216）和RAS-wt（n=226）队列显示，尽管MSS RAS-mt肿瘤携带略高的突变频率，但驱动基因谱、突变类别和受影响通路在两组之间保持广泛保守。

拷贝数变异（CNA）：

使用GISTIC分析，结合整理的Angiopredict（APD）扩展队列（n=880）和COLOSSUS St2/3 MSS RAS-mt队列（n=182；总计n=1,062）的CNA数据。
在COLOSSUS队列中，高水平扩增和低水平增益占主导，而缺失仅限于低水平事件。
选定的CNAs与全基因组不稳定性测量密切相关，包括基因组异常分数（FGAs）和非整倍性。
大多数频繁CNAs与COLOSSUS St2/3 MSS RAS-mt队列中相应的基因表达变化显著相关（logistic回归；FDR<0.2），表明CNAs是这些基因表达水平的主要决定因素。
CNA影响基因的通路富集分析识别出PI3K和RAS信号（与突变富集通路相似），以及细胞衰老和代谢通路。
这些分析表明，MSS RAS-mt和RAS-wt CRC共享广泛相似的基因组架构，尽管突变频率存在定量差异，支持研究MSS CRC中常见的基因组机制和治疗脆弱性。

突变-CNA关联：

应用基于稳健回归的框架与广泛重采样，识别可重现的突变-CNA关联。
PIK3CA和TP53突变与较高的FGAs一致相关（>10%迭代/0.1频率），跨越COLOSSUS和TCGA MSS CRC队列，无论RAS状态如何。
在COLOSSUS MSS RAS-mt肿瘤中，PIK3CA突变与臂水平CNAs相关，包括8q扩增（包含PTK2基因；其蛋白表示为黏着斑激酶[FAK]）和8p和18q缺失。
TP53突变同样与7p11.2（EGFR）扩增和所有MSS CRC队列中的多个复发性CNA事件相关，与其在染色体不稳定性中的作用一致。
这些发现突出MSS CRC中可操作靶点/通路——PI3K-AKT、EGFR和FAK——的潜在协调改变。

CRISPR依赖性分析：

分析MSS CRC细胞系（n=30）中来自Cancer Dependency Map（DepMap）的CRISPR-Cas9敲除依赖性数据。
在COLOSSUS St2/3和整理的APD扩展队列中识别的114个显著（GISTIC，FDR调整q<0.2）CNA基因中，22个也在CRISPR分析中显示依赖性（基因效应评分阈值<-0.25）。
snRNA-seq/scRNA-seq分析MSS CRC样本（n=38；195,246个单核/细胞）进一步识别出EGFR和PTK2（FAK），它们在上皮细胞（以及成纤维细胞）中高表达，且已有药物获批或处于癌症临床试验中。
鉴于抗EGFR治疗已临床批准用于RAS-wt CRC，FAK被优先用于功能验证，作为MSS CRC中潜在可操作的治疗靶点。

结果3、代谢组学分析通过多组学整合识别MSS CRC中保守的代谢程序

代谢组学分析：

对COLOSSUS St2/3 MSS CRC队列（n=177）通过LC-MS/MS进行1,478种非靶向和靶向代谢物分析，识别出207种靶向代谢物（129种独特）。
基于相关性的聚类揭示两个主要相互连接的代谢物簇，具有实质代谢异质性。
使用MetaboAnalyst 6.0的通路富集识别出8个显著通路（FDR<0.05），包括丙氨酸、天冬氨酸和谷氨酸代谢。
使用基于GSVA的ssGSEA的跨组学整合显示，其中6个通路在COLOSSUS St2/3 CRC队列的代谢组学、转录组学和蛋白质组学数据集中显示一致富集。
这些富集模式进一步在同一队列中使用singscore（一种独立的基因集评分方法）以及一个具有匹配RNA-seq和代谢组学数据的独立验证队列（GSE89076；n=37）中得到支持。
在MSS肿瘤类器官（n=12）中，具有匹配转录组和代谢组数据，丙氨酸、天冬氨酸和谷氨酸代谢保持持续富集。

结果4、区室化谷氨酸-谷氨酸代谢程序的特征识别治疗脆弱性

谷氨酸代谢详细研究：

在该通路中，7种代谢物被MetaboAnalyst显著富集。
线性回归分析整合代谢物与基因和蛋白质表达（GLS、GLUL和GOT2）显示谷氨酰胺代谢的分子层特异性调控：
- 在转录水平，GLUL（催化谷氨酸合成谷氨酰胺）与谷氨酰胺显示显著正相关（p<0.05）。
- 虽然GLS和GOT2未显示显著的转录水平关联，但GLS（将谷氨酰胺转化为谷氨酸）表现出与其在谷氨酰胺利用中作用一致的负趋势。
- 然而，GLS蛋白丰度与谷氨酰胺水平显著相关，表明谷氨酰胺利用在蛋白质水平比转录水平更密切反映。
谷氨酸丰度在肿瘤中超过谷氨酰胺，两种代谢物强相关（Pearson r=0.84，p<0.001），表明协调的谷氨酰胺-谷氨酸代谢。
这一关系在GSE89076队列中独立重现。

单核RNA-seq验证区室特异性组织：

分析COLOSSUS MSS KRAS-mt CRCs的snRNA-seq（n=14；84,935个单核）。
丙氨酸、天冬氨酸和谷氨酸代谢在上皮区室中富集，而单个酶表达显示细胞特异性：
- GLUL优先在髓系细胞中富集。
- GLS和GOT2在恶性上皮细胞中明显。
与这些发现一致，TCA循环基因在恶性和正常上皮细胞中均相对富集，表明上皮肿瘤细胞在谷氨酸连接的氧化代谢中发挥作用。

细胞组成调整：

使用MuSiC工具和COLOSSUS snRNA-seq数据的bulk转录组反卷积显示，调整上皮、髓系和成纤维细胞组成后保留了正谷氨酰胺-GLUL关联。
尽管GLS和GOT2显示类似方向趋势，但其关联在转录水平保持较弱或不显著，支持GLS调控超越转录丰度。
在COLOSSUS MSS RAS-wt snRNA-seq队列和两个独立公共MSS CRC scRNA-seq数据集中观察到可比的区室特异性代谢组织，支持这些发现在MSS CRC队列中的稳健性。

代谢组学与RPPA整合：

蛋白质丰度（总和磷酸化形式）与连接到mTOR和AKT信号、凋亡（CASP9）和缺氧（HIF1A和CA9）的代谢物正相关。
总FAK也与十七烷酸、2-甲基苯甲酸和HMMNA相关，将脂质代谢与FAK信号联系起来，进一步暗示这些蛋白质及其相关通路作为潜在治疗机会。

结果5、增强子甲基化突出与MSS CRC治疗通路相关的汇聚基因组和代谢改变

增强子特异性CpG甲基化分析：

分析MSS CRC中增强子特异性CpG甲基化模式与转录组和CNA谱的关系。
使用ELMER在原发性肿瘤（TCGA MSS RAS-mt，n=153）和匹配的正常邻近组织（n=8）之间进行差异增强子特异性甲基化分析。
识别出2,995个低甲基化和776个高甲基化增强子特异性CpG探针。
与邻近基因表达整合识别出4,140个低甲基化和893个高甲基化CpG-基因对，显示显著的逆增强子-基因关联（ELMER FDR<0.05；Pearson r<-0.2）。

PTK2（FAK）和MYC的增强子甲基化：

在增强子低甲基化相关基因中，7个CpG增强子探针与PTK2（FAK）基因表达显著负相关。
此外，37个扩增和27个缺失基因（来自合并的COLOSSUS St2/3和整理的APD扩展队列；n=1,062）在这些肿瘤中的增强子低/高甲基化基因集中富集。
这些发现支持CNAs子集的汇聚基因组和表观遗传调控，包括PTK2（FAK），并表明增强子甲基化在基因组改变之外对转录调控有贡献。

基序富集分析：

低甲基化增强子区域的基序富集分析识别出几个转录因子基序，最突出的是AP-1（FOS/JUN家族）。
鉴于AP-1、MYC信号和AKT-mTOR通路激活之间的既定联系，评估了MYC通过增强子甲基化的调控。
多个增强子低甲基化CpG位点显示与MYC表达的负相关，并与这些CRC样本中邻近基因的转录相关。
一致地，MSigDB Hallmark分析揭示MYC靶点、mTOR信号和E2F通路在增强子低甲基化基因中富集。
这些结果支持增强子相关甲基化和MYC驱动转录程序之间的联系。

谷氨酰胺分解相关基因：

代谢组学工作台基于基因集分析进一步显示，增强子低甲基化基因谱与谷氨酰胺和谷氨酸相关代谢通路相关，表明在该队列中谷氨酰胺/谷氨酸代谢富集。
与先前观察到的谷氨酰胺代谢的多分子层调控一致，GLS表达与两个增强子甲基化探针显著负相关，这两个探针在肿瘤中相比正常组织显著差异低甲基化（p<0.001），表明通过增强子甲基化调控。

结果6、多组学数据融合识别MSS CRC中临床相关分子程序

iPhenMap分析：

同时整合多组学数据与临床协变量，应用整合表型映射（iPhenMap），一种与先前发表的PhenMap框架相关的内部半监督机器学习工具，应用于COLOSSUS St2/3队列的98个样本。
iPhenMap识别出6个最优映射变量（MVs），代表跨多组学数据类型及其与临床特征关系的潜在模式。
MV1主要与转录组和蛋白质组特征相关，而MV4主要与转录组和CNA特征相关。相反，MV2、3、5和6主要由代谢组学信号驱动，特别是非靶向代谢组学。
这些MVs捕获了跨组学可变的稳健生物学。

独立验证：

使用多组学因子分析（MOFA）进行独立验证，一种无监督多组学因子分析框架，识别出10个由转录组、蛋白质组、代谢组和CNA变异驱动的因子，与iPhenMap MVs显著相关。
此外，将iPhenMap应用于独立TCGA多组学CRC队列（n=64），该队列仅包含转录组和蛋白质组数据及有限临床协变量，产生广泛一致的分子程序，特别是MV1。
总体而言，这些发现支持iPhenMap捕获一致、生物学上有意义的多组学景观的能力，跨分析和独立队列观察到的协调模式一致。

临床关联：

在iPhenMap工具中实施的整合多变量回归分析识别出个体MVs与多个临床病理特征之间的显著关联，包括神经周围侵犯（PNI）、诊断时ECOG状态、淋巴血管侵犯、组织学类型、肿瘤位置、分组分期类别、性别和复发。
值得注意的是，右侧肿瘤与MV1和MV4相关，表明右侧性在MSS RAS-mt CRC中可能发挥潜在作用。
一致地，MOFA因子和来自TCGA队列的iPhenMap MVs与临床协变量的关联也与COLOSSUS队列中的iPhenMap MVs一致。

结果7、MV1定义与基因组不稳定性、免疫排斥和治疗脆弱性相关的复制应激程序

MV1的生物学基础：

为定义MV1的生物学基础，检查了显著载荷（>2或<-2）的基因和蛋白质。
通路富集（MSigDB-based、KEGG、REACTOME和BioCARTA）揭示MV1阳性基因和蛋白质中的先天免疫和炎症信号程序，提示肿瘤内在程序与右侧CRC生物学相关。
然而，MV1阳性特征也与DNA复制和复制应激通路相关，主要由微型染色体维持复合物（MCMs）3、4和7驱动。
有趣的是，MCM3、MCM4和MCM7表达在转录和蛋白质水平均与多个独立队列中的复制应激签名评分显著相关（p<0.05），包括MSS RAS-mt队列（COLOSSUS St2/3；n=134，COLOSSUS mCRC；n=27和TCGA；n=159），MSS RAS-wt TCGA（n=136）队列和混合RAS-mt/wt肿瘤类器官队列（n=131）。
这些MCM基因也与AURKB表达适度正相关，AURKB是有丝分裂期间染色体分离和胞质分裂的关键调节因子，在扩展bulk转录组数据中。
一致地，MCM3、MCM4和MCM7蛋白表达与COLOSSUS St2/3队列中的FGAs中度相关，将复制应激与基因组不稳定性联系起来。
这些观察支持一个模型，其中选定的MCM相关复制应激可能与基因组不稳定性（FGAs）相关，可能通过涉及AURKB的失调有丝分裂控制。这可能代表复制应激丰富的MSS CRC中的治疗靶点。

MV1与右侧性：

尽管MV1阳性肿瘤显示右侧CRC的适度富集，但MCM3、MCM4和MCM7的基因表达在独立混合CRC队列中未能区分MSS与MSI肿瘤，表明该程序反映增殖而非微卫星状态。
尽管如此，选定的MCM基因和蛋白质表达与MMR基因/蛋白质表达强烈相关，连接DNA复制和修复。

免疫景观：

尽管肿瘤内在先天免疫和炎症转录程序富集，MV1-high肿瘤未表现出广泛免疫炎症的肿瘤微环境。
免疫反卷积（使用MCP-counter）在COLOSSUS St2/3队列和已发表TCGA免疫谱中广泛揭示MV1评分（包括选定的MCM表达）、FGAs和巨噬细胞/髓系群体之间的 predominantly 负相关，仅与选定淋巴样签名适度正相关。

单细胞分析：

单细胞分析进一步阐明了与MV1相关程序相关的细胞类型。COLOSSUS snRNA-seq和独立MSS scRNA-seq数据集（GSE200997和GSE178318）显示MV1阳性基因 predominantly 在上皮细胞中表达，在免疫和基质区室中表达极少。
选定的MCM表达在S期细胞周期群体中富集，支持增殖驱动程序。
与MSS CRC中免疫排斥一致，这些发现表明复制应激和基因组不稳定性可能有助于免疫逃逸，尽管与MSI相关转录程序部分重叠。

蛋白酶体脆弱性：

鉴于MSS CRC中免疫治疗脆弱性的有限证据，探索了与MV1相关的替代靶点。泛素-蛋白酶体系统与DNA复制密切相关；蛋白酶体签名正向与多个队列中的MCM3、MCM4和MCM7基因表达相关，提示蛋白酶体相关治疗脆弱性涉及bortezomib。
相反，MV1阴性基因富集于基质和细胞外基质程序，而MV4与黏液组织学、右侧肿瘤位置和染色体8q24.3（PTK2/FAK）扩增相关。
与这些观察一致，MV1-low肿瘤显示PTK2（FAK）信号基因集的富集，进一步支持FAK作为MSS CRC中的额外治疗候选。
总之，这些发现突出MV1相关生物学，特征为部分由MCM驱动的复制应激和MSS CRC中的基因组不稳定性，识别超越免疫治疗的潜在治疗脆弱性。

结果8、MSS RAS-mt原发性CRC中与PNI相关的代谢重编程揭示信号重连和治疗脆弱性

PNI相关代谢特征：

研究MSS CRC中与临床特征相关的代谢变化，聚焦于两个关键代谢模式，包含23个已知代谢物（20个独特），通过iPhenMap在COLOSSUS St2/3 MSS RAS-mt CRC（n=98）中贡献于MV5（载荷>2或<-2）。
MV5评分在PNI阳性肿瘤中适度更高（p=0.024），表明PNI状态与MV5潜在变量中共享的代谢特征之间的关联。

代谢物-PNI关联：

使用重采样线性回归（70%训练和30%测试集，1,000次迭代）分析161个肿瘤样本，包括98个来自原始iPhenMap数据集和63个额外COLOSSUS St2/3样本。
23个MV5相关代谢物中的10个在重采样迭代中显示与PNI的可重现关联（出现频率>10%）。
这些代谢物包括脂肪酸-棕榈酸、氨基酸相关代谢物（如酮亮氨酸、2-乙基羟基乙酸）以及与线粒体和氧化过程相关的代谢物（如丁酸和氢化肉桂酸），与先前在肿瘤-神经相互作用背景下报道的代谢特征一致。
富马酸，与线粒体代谢相关，是唯一与PNI负相关的代谢物。

代谢物-PNI-蛋白质相互作用：

检查这5个正PNI相关代谢物与总蛋白质表达的关系。
使用2,693个高变LC-MS/MS蛋白质（SD>0.5），通过重复训练-测试重采样（100次迭代）建模代谢物-PNI相互作用效应。
相互作用估计显示训练和测试数据集之间一致效应方向的趋势。
为提高稳健性 despite 有限统计功效，优先考虑高训练显著性频率（≥0.6）且在匹配测试集中效应方向一致的相互作用，识别出81个具有一致代谢物-PNI相互作用效应的蛋白质。
在47个显示正相互作用信号的蛋白质中，通路分析表明富集MYC和E2F转录因子靶点、G2-M检查点、凋亡、mTOR通路以及氧化磷酸化和血红素代谢。
转录因子分析进一步突出MYC、E2F靶点、BRCA1和相关调节因子，而毒理基因组富集识别出bortezomib反应签名在顶部类别中，将这一代谢程序与蛋白质毒性应激反应和潜在治疗脆弱性联系起来。

信号水平整合：

将代谢组学数据与COLOSSUS St2/3队列的RPPA谱（子队列n=121；67种蛋白质）整合，使用相同的重复训练-测试重采样策略。
在选定代谢物与总和磷酸化蛋白质之间观察到可重现的正相关，包括磷酸化AKT、SRC、HIF1α和磷酸化STAT3。
这些发现暗示AKT、SRC和缺氧相关信号通路，表明PNI相关代谢程序在磷酸化蛋白质信号水平上得到反映。
这些发现值得在MSS RAS-wt CRC队列中进行独立验证。总之，这些发现识别出MSS RAS-mt CRC中PNI相关代谢程序，特征为协调的代谢、蛋白质组和磷酸化蛋白质组信号变化，汇聚于MYC-E2F、AKT、SRC和缺氧通路，突出潜在治疗脆弱性。

结果9、多组学指导的功能验证识别MSS RAS-mt CRC肿瘤类器官中的治疗脆弱性

候选通路和脆弱性：

汇总通过MSS CRC广泛多组学分析识别的候选通路和治疗脆弱性，按跨分析频率排名。
最高排名的脆弱性包括：
- PDK1/2激酶信号（与PI3K-AKT-mTOR相关）
- FAK激酶活性
- 蛋白酶体功能
- 线粒体谷氨酰胺分解驱动的代谢重编程
- 缺氧应激反应
- SRC通路激活
- 有丝分裂复制应激（与AURKB相关）
尽管HIF1A作为潜在治疗靶点出现，但选择临床先进的HIF2A抑制剂PT-2977进行离体功能测试，因为其能够抑制重叠的缺氧驱动转录程序。
鉴于未满足的临床需求，MSS RAS-mt CRC被优先用于患者衍生模型中的功能验证。

药物测试：

在来自转移性MSS RAS-mt PDX肿瘤的CRC肿瘤类器官中进行药理学测试，这些肿瘤在突变和CNA方面与患者肿瘤相对相似。
从131个模型面板中选择肿瘤类器官，以捕获所选治疗靶点识别的分子异质性，同时保持稳健的离体生长特征。
所有化合物均与MEK1/2抑制剂trametinib联合评估，基于MEK抑制增强RAS-mt肿瘤对其他治疗干预敏感性的证据。

结果：

在最高排名靶点中，GSK-2334470抑制PDK1/2、defactinib抑制FAK和bortezomib抑制蛋白酶体，各自与trametinib联合，在大多数测试模型中显著降低肿瘤类器官活力：
- GSK-2334470显著抑制所有9个检查模型的生长。
- Defactinib在9个模型中的8个中增强trametinib的抗增殖活性。
- Bortezomib在9个模型中显著增强trametinib的效果，显著降低细胞活力，在几个模型中甚至最低剂量33 nM就接近完全丧失活力细胞。
这些发现表明，靶向特定信号通路或功能特征可以使MSS RAS-mt CRC肿瘤类器官对MEK阻断敏感。
相反，SRC抑制（bosutinib）、谷氨酰胺酶抑制（telagenastat）、HIF-2α抑制（PT-2977）或Aurora B激酶抑制（barasertib） 与trametinib联合提供有限或无额外益处，表明在多组学分析中一致性较低的识别的通路在患者衍生肿瘤类器官模型中具有较低的功能相关性。

结果10、体内验证研究

PDX模型验证：

为验证最有前景的与trametinib联合的治疗组合，匹配的PDX模型用单药或与trametinib联合治疗5周。
短期疗效通过肿瘤生长抑制评估，长期疗效通过生存直至预定义研究终点评估。
与离体发现一致，GSK-2334470联合trametinib与安慰剂相比在两个（共三个）PDX模型中显著降低肿瘤生长，并在两个（共三个）模型中相比安慰剂或GSK-2334470单药延长生存。
为靶向FAK信号，使用VS-4718（defactinib的结构类似物，具有改善的口服生物利用度和药效学特性）。VS-4718与安慰剂相比在两个（共三个）模型中显著抑制肿瘤生长，无论是单药还是与trametinib联合。类似地，与trametinib联合给药时在两个（共三个）模型中相比安慰剂延长生存。
相反，bortezomib在检查的两个模型中对肿瘤生长影响有限，尽管在一个模型中显著延长生存。
总之，这些离体和体内研究验证了PDK1、FAK和蛋白酶体抑制作为通过MSS RAS-mt CRC整合多组学分析识别的治疗可操作脆弱性。
尽管个体患者衍生模型的治疗反应各异，但PDK1和FAK抑制在两个（共三个）独立肿瘤类器官和PDX模型中的一致活性支持优先考虑这些通路进行进一步临床前和临床研究。

讨论

主要贡献

通过整合来自约2,000个样本和195,246个单细胞/细胞核的多组学和临床数据，来自12个内部和公开可用的患者和患者衍生肿瘤类器官/PDX队列，提供了MSS CRC的全面分子表征。
重要的是，大幅扩展了可用于更广泛研究社区的RAS-mt CRC资源，增加了来自316个样本的跨平台分子表征数据，包括>1,000个分子谱以及来自119,489个细胞核的单核数据。
尽管先前大规模研究已定义CRC的分子景观，但MSS CRC此前未被系统地使用基因组（突变和CNA）、代谢组学、转录组学（包括单细胞）、蛋白质组学（总和磷酸化蛋白质组学）、表观基因组学和机器学习驱动的多组学数据融合方法进行研究。
将这些互补分子层与临床协变量结合，使得能够识别协调的生物学程序、机制见解和治疗可操作脆弱性，随后在PDX及其肿瘤类器官模型中进行评估。

基因组分析

基因组分析确认影响PI3K、WNT和TP53通路的复发性改变，同时揭示MSS RAS-mt和RAS-wt肿瘤之间的背景依赖性差异。
值得注意的是，PIK3CA突变在MSS RAS-mt中频率约高2倍，而PTEN丢失相关CNAs在PIK3CA突变MSS RAS-wt肿瘤中富集，提示PI3K通路激活的额外机制。
在MSS CRC中，PIK3CA突变与CNAs的复发性共现提示对持续PI3K信号的协调选择，而非完全独立的基因组事件。
鉴于PI3K-AKT-PDK1和黏着斑信号之间的既定串扰，这些改变可能共同促进肿瘤表型，突出PI3K轴作为基因组定义MSS CRC中的背景依赖性治疗脆弱性。

FAK作为中心节点

与此一致，拷贝数与转录组和表观基因组数据的整合表明，许多复发性CNAs，包括PTK2（FAK）扩增，伴随着基因表达的相应变化，而增强子低甲基化提供了独立于基因剂量的额外调控PTK2的机制。
磷酸化FAK与脂质代谢相关代谢物之间的正相关进一步表明，黏着斑信号与更广泛的代谢重编程功能整合。
此外，PTK2在DepMap CRISPR数据集中成为复发性CNA相关基因中的顶级依赖性候选，单细胞分析显示 preferential 上皮表达，尽管基质FAK分布更广泛。
总之，这些互补观察支持FAK作为连接MSS CRC中基因组、表观遗传和代谢调控的中心节点。
这些发现突出整合分析超越CNAs以完全捕获MSS CRC基因调控复杂性的重要性。
此外，FDA最近批准了avutomectinib（一种RAS/MEK抑制剂）与FAK抑制剂defactinib联合用于KRAS-mt复发性低级别浆液性卵巢癌。与这一临床先例一致，作者的MSS RAS-mt CRC离体和体内研究显示对FAK抑制的模型依赖性反应，包括与MEK抑制联合，进一步支持FAK作为这一临床重要亚组中的治疗可操作脆弱性。

代谢分析

MSS CRC特异性代谢组学分析仍然有限。
整合代谢分析表明，谷氨酰胺代谢在 distinct 调控层中组织，而非单一转录程序。
GLUL在转录水平和GLS在蛋白质水平的差异调控支持肿瘤生态系统内的区室化代谢组织。
尽管增强子低甲基化提示GLS表达的 epigenetically permissive 状态，但其功能活性似乎在 downstream 调控。
在缺乏免疫/基质细胞的离体模型中对GLS抑制的有限反应进一步突出细胞背景和代谢异质性在塑造谷氨酰胺酶依赖性中的重要性。
总之，这些观察表明，在MSS CRC中有效靶向谷氨酰胺代谢需要考虑调控分层和肿瘤-微环境相互作用。
相反，FAK和PI3K-PDK1-AKT-mTOR轴在整合分析中作为复发性脆弱性出现，随后在多个肿瘤类器官和PDX模型中得到验证。

表观遗传分析

与这些观察一致，表观遗传分析进一步显示广泛的增强子低甲基化涉及AP-1、MYC和E2F转录网络，连同mTOR通路组分的富集，为MSS CRC中靶向PI3K-PDK1信号提供额外支持。
复发性CNAs和增强子相关基因之间的显著重叠，特别是涉及PTK2（FAK），表明结构基因组改变和增强子失调合作增强致癌转录程序。
总之，这些发现说明整合基因组和表观遗传调控如何促进MSS CRC中生物学和治疗相关的脆弱性。

iPhenMap与MV1

本研究的一个关键方面是使用iPhenMap将多组学数据与COLOSSUS St2/3队列的临床协变量整合。
这种半监督多组学和表型整合框架定义了六个与转录组学、蛋白质组学和代谢组学程序相关的MVs，具有临床相关性。
其中，MV1出现为最稳健的，特征为MCM3、MCM4和MCM7，DNA复制机器的关键组分。
它们与FGAs的关联支持MV1、复制应激和基因组不稳定性之间的联系。
持续的复制应激可能增加对蛋白质稳态和蛋白酶体通路的依赖，这些通路在致癌应激下维持肿瘤细胞存活，与作者整合分析中观察到的bortezomib相关基因签名的富集一致。
MV1相关程序也优先在右侧MSS CRC中富集，并 predominantly 定位于上皮区室。
尽管右侧CRC传统上与MSI和免疫丰富表型相关，但作者的发现表明，一部分MSS右侧肿瘤表现出肿瘤内在程序，特征为复制应激连同先天免疫和炎症信号，尽管免疫细胞浸润有限。
这些观察扩展了当前对右侧CRC生物学的理解超越解剖位置，值得进一步研究其生物学和治疗意义。

AURKB

尽管染色体分离机器的多个组分被涉及，但AURKB因其与MV1、复制应激的关联及其在该网络中的中心调控作用而被优先考虑。
然而，AURKB抑制未能在患者衍生肿瘤类器官中显示活性。这表明AURKB可能主要作为有丝分裂失调的标志物而非必然依赖性，可能由于其他有丝分裂激酶的补偿。
这一结果进一步说明区分分子关联与功能治疗依赖性的重要性。

PNI

多组学融合与临床协变量突出PNI作为与原发性肿瘤中协调代谢和蛋白质组重塑显著相关的关键临床特征。
PNI相关代谢签名涉及10种代谢物连同缺氧、AKT、SRC和STAT3蛋白。
与PNI和代谢物相互作用相关的蛋白质额外富集bortezomib相关签名和与复制应激连接的MYC-E2F程序。
这些发现表明，PNI阳性MSS原发性CRC可能代表一个生物学上 distinct 的亚组，特征为代谢适应、应激反应信号和蛋白质稳态依赖。
尽管这一假设需要在独立队列和适当的功能模型中进行验证，但它为评估破坏PNI阳性疾病中这些汇聚通路的治疗提供了理由。

MEK抑制作为骨干

尽管直接KRAS抑制正在快速发展，但其目前在CRC中的适用性仍受突变特异性覆盖和适应性MAPK再激活的限制。
因此，作者选择MEK抑制（通过trametinib）作为通路水平策略，以在遗传异质性RAS-mt MSS CRC模型中抑制MAPK信号，同时为评估联合治疗提供共同骨干。

体内验证

在CRC肿瘤类器官中显示活性的组合中，匹配的MSS RAS-mt mCRC PDX模型中的体内测试识别出与VS-4718（FAK抑制剂）、GSK-2334470（PDK1抑制剂）和bortezomib（蛋白酶体抑制剂）的可变但令人鼓舞的反应。
VS-4718，单独和与MEK阻断联合，显示对肿瘤生长和生存的模型依赖性效应。
总之，这些发现支持在MSS RAS-mt CRC中进一步研究FAK导向的治疗策略，同时也强调治疗反应的异质性。未来研究应定义对FAK、MEK、PDK1和蛋白酶体靶向治疗反应的生物标志物，并确定这些脆弱性是否延伸到MSS RAS-mt之外到MSS RAS-wt CRC。

结论

总之，MSS CRC的深入分子表征揭示了生物学见解和有前景的治疗机会。
通过无缝整合基因组、表观遗传、蛋白质组、代谢组和snRNA-seq/scRNA-seq表达分析，大幅推进了对该疾病复杂和异质性生物学的理解，同时识别出有前景的新脆弱性背景，具有进一步转化潜力。
鉴于MSS CRC患者缺乏有效治疗，这些脆弱性现在可能为精准肿瘤学治疗框架提供基础，旨在改善这一难治性疾病背景下的结局。

局限性

外部验证数据集并非在所有COLOSSUS队列中代表的分子层中均统一可用，限制了选定发现的直接复制。
在直接正交验证不可行的情况下，使用互补分子数据类型和重采样分析支持生物学一致性，包括增强子甲基化相关转录组签名和PNI相关代谢物-蛋白质相互作用。
重要的是，作者的逐步靶点优先排序策略，从整合多组学和表型分析到离体和体内药物验证，并未确认所有计算优先排序的靶点的治疗相关性。

方法总结

全外显子测序（WES）：

体细胞突变检测：

WES数据使用现有内部生物信息学管道分析。原始测序reads使用Burrows-Wheeler Aligner映射到人类参考基因组（GRCh38）。使用Picard去除重复，使用GATK进行indel重比对和碱基质量重校准后，使用MuTect2识别体细胞变异，并用Annovar注释。
在缺乏匹配种系样本的情况下，使用公开数据库应用严格过滤策略以排除常见SNPs。Annovar用于去除gnomAD中发生等位基因频率>0.001的变异。基于以下标准保留高置信度变异：最小覆盖深度≥10x，变异等位基因频率≥0.1，≥2个支持替代reads。所得VCF文件使用vcf2maf包转换为MAF文件。排除一个>1,000个突变变异（超突变）的样本，留下88个WES样本用于下游分析。

TCGA突变数据处理和队列分层：

从TCGA PanCancer Atlas项目通过cBioPortal获取结肠腺癌（COAD）和直肠腺癌（READ）的体细胞突变和临床注释数据。MAF文件和相应临床样本数据使用maftools包导入R。
微卫星不稳定性（MSI）状态使用cBioPortal提供的MSIsensor评分确定。MSIsensor评分≥20的肿瘤分类为MSI，而<20的认为MSS。在MSS组内，使用非同义肿瘤突变负荷（TMB）截止值≥10突变/Mb排除超突变肿瘤。
在MSS队列内，肿瘤进一步按RAS突变状态分层。携带KRAS或NRAS突变的样本分类为RAS突变（RAS-mt），而两个基因均无改变的指定为RAS野生型（RAS-wt）。

低覆盖度全基因组测序（LC-WGS）：

使用QIAamp DNA FFPE Tissue Kit从FFPE组织切片提取DNA。使用NanoDrop 2000分光光度计定量DNA样本。
对于CNA profiling，DNA样本使用先前优化的LC-WGS协议测序。简言之，使用KAPA HyperPlus Library Preparation Kit从DNA制备shotgun全基因组DNA文库。全基因组文库在Illumina测序仪上以低覆盖度（约0.1x）测序（HiSeq 4000单端1×51 bp reads或NovaSeq单端1×100 bp reads）。
CNAs使用现有内部生物信息学管道识别。原始测序reads使用Burrows-Wheeler Aligner（v0.7.15）比对到人类参考基因组（hg19），生成SAM文件。Samtools（v1.3.1）将SAM转换为BAM。BAM文件使用Picard tools（v2.5.0）进行重复标记。

DNA拷贝数和GISTIC分析：

使用QDNAseq包（v1.22）在R（v3.6.3）中创建片段，使用ASCAT方法（v2.0.7）推断等位基因特异性拷贝数谱。
使用GISTIC-2.0（v2.0）识别最频繁和过度代表的CNAs，置信阈值0.05（5%）。
增益/扩增跨阈值离散GISTIC调用表示（0-二倍体，1-低水平增益，2-高水平扩增），而丢失/缺失使用较低评分（0-二倍体，1-低水平缺失，2-高水平丢失）。
对于GISTIC分析，COLOSSUS St2/3（n=182）和三个不同批次的COLOSSUS mCRC MSS RAS-mt队列（总n=121）独立与整理的APD扩展队列（n=880）的segmented数据合并。来自整理队列的COLOSSUS St2/3 MSS RAS-mt子集的GISTIC评分在下游分析前使用sva R包中的ComBat工具校正批次效应。

CNA与基因表达的关联：

GISTIC2.0识别的CNAs跨患者分析以确定每个区域内事件的频率。对于每个GISTIC2.0峰，区域内所有基因的平均RNA表达使用logistic回归评估，并关联统计检验以确定与CNA状态的相关性。

突变和CNA相互作用分析：

使用logistic回归测试顶部突变基因与基因组异常分数（FGA）的关系，FGA作为预测因子，二元突变状态作为结果。
使用重复75/25训练-测试分割跨多个（100次迭代）随机种子和每个种子1,000次迭代实现基于重采样的框架。对于每个基因，从每个迭代的二项广义线性模型中提取回归系数和相关p值。稳健性量化为重采样迭代中突变状态与FGA之间关联达到名义显著性（p<0.05）的比例。
CNA事件与显著突变（PIK3CA/TP53）之间的关联使用Fisher精确检验进一步评估，比较突变（0-野生型 vs 1-突变）与CNA存在（GISTIC2评分0 vs 1或2）。

细胞系DepMap衍生的CRISPR依赖性和基因组数据：

从DepMap门户收集MSS CRC细胞系（n=35；21 RAS-mt和14 RAS-wt）的CRISPR依赖性评分（DepMap Public_24Q2+Score_Chronos），用于显著扩增区域内的所有基因。匹配DepMap的细胞系突变和CNA数据来自Sanger's Cell Model Passports。

单核RNA测序（snRNA-seq）：

使用Chromium Single Cell 3'试剂盒v3.1生成单核RNA测序文库。核加载到Chromium Single Cell G Chip以获得约25,000个单核转录组回收。文库在NovaSeq 6000平台上测序。

snRNA-seq分析：

原始测序reads比对到GRCh38人类参考基因组（版本3.0.0），使用Cell Ranger（v5.0.0）处理为基因表达矩阵。个体样本（共18个）使用CellBender（版本0.2.2）校正环境RNA污染。'expected cells'参数设置为10,000，'total droplets included'参数设置为25,000，所有其他参数默认。
处理后的矩阵使用Seurat（v4.3.0）在R中分析。核过滤保留每个核表达>1,000和<5,000基因且线粒体RNA<5%的核。使用DoubleFinder（v2.0.3）检测双细胞。数据使用Seurat的NormalizeData函数归一化，并在回归转录本数、线粒体RNA百分比和细胞周期评分时缩放。使用PCA对top 2,000可变基因进行降维。样本使用Harmony（v2.0.3）默认设置整合，考虑样本特异性效应。细胞聚类后，共识别出119,489个核。基于经典标记基因，使用表S10提供的标记注释细胞：髓系细胞、内皮细胞、成纤维细胞、上皮细胞、B细胞和T细胞。

公共数据的单细胞RNA-seq（scRNA-seq）分析：

来自GEO Omnibus的scRNA-seq数据在Seurat（v4.1.1）中处理。对于先前发表和公开可用的GSE200997数据集，仅保留MSS肿瘤衍生细胞，包括14个MSS样本的24,764个单细胞。基因在少于三个细胞中检测到的和细胞检测基因少于200个的被排除，线粒体转录本>25%的细胞被移除。使用scDblFinder（v1.0）识别双细胞。
对于公开可用的GSE178318数据集，仅保留治疗初治样本，包括原发性肿瘤（n=3）和匹配的肝转移（n=3）的50,993个单细胞。细胞基于基因复杂性（≥500检测基因）、线粒体含量（≤15%）和原始UMI计数和基因特征的批次特异性阈值（均值±3 SD内）过滤。在少于三个细胞中检测到的基因和非蛋白质编码基因被排除。
对于上述数据集，数据对数归一化，使用方差稳定变换方法识别高变基因，通过PCA降维，然后跨样本Harmony批次校正。邻居图构建、Louvain聚类和UMAP嵌入使用Harmony校正维度进行。广泛细胞类型注释（表S10）使用SingleR和Human Primary Cell Atlas参考分配，并合并为六个主要区室。

基于scRNA-seq/snRNA-seq的恶性上皮细胞定义：

使用scATOMIC R包（默认转录组签名）将上皮状态分类为恶性/癌性 vs 正常。连续scATOMIC恶性评分仅用于排序和可视化，因为其分布在这些状态之间重叠很大。

基于scRNA-seq/snRNA-seq的基因富集分析：

使用escape R包中实现的ssGSEA方法进行单细胞基因集富集分析。自定义基因集应用于Seurat对象，使用runEscape（groups=1,000，最小基因集大小=5）计算每个细胞的富集评分。富集分布使用escape R包的geyser图可视化，而降维和基因表达模式使用Seurat R包生成的UMAP投影和特征表达图显示。

非靶向和靶向代谢组学：

组织切片在650 μL 60%甲醇中在RiboLyser管中匀浆。样本在-80°C超低温冰箱中储存过夜，第二天在4°C、15,000×g离心15分钟。蛋白质沉淀用于归一化代谢物数据，使用BCA蛋白 assay 定量。
将500 μL上清转移到新鲜微量离心管中，加入500 μL氯仿。混合物在热混匀器（28°C，1,000 rpm）上摇动2小时。随后，将200 μL上层极性相转移到LC-MS小瓶。下层非极性相（450 μL）在SpeedVac浓缩器中干燥，重悬于200 μL 100%乙醇中，转移到LC-MS小瓶。对于每个相（极性/非极性），通过从每个样本收集10 μL生成QC样本。
肿瘤类器官沉淀使用类似的提取工作流程。
FIA-MS：使用Infinity II自动进样器和泵耦合到scimaX磁共振质谱仪，正负离子模式运行。极性和非极性组分随机顺序分析，每个样本测量三次。QC样本每15次运行注射一次以监测批次效应和仪器稳定性。样品通过80% 2-丙醇连续流引入MS。MS在全扫描模式下运行，m/z范围107-800。

代谢组学数据分析和代谢物注释：

数据采集使用FTMS Control和Compass Hystar软件。数据分析（峰提取、注释和批次效应强度校正）使用MetaboScape 2022b进行。代谢物注释通过比较精确质量和同位素配置进行。注释数据库为Human Metabolome Database。
COLOSSUS患者队列的代谢组学数据进行log2变换和批次校正（组织重量、样本收集中心和5次不同 assay 运行），使用sva R库中的ComBat工具。在每个步骤进行PCA以评估和校正批次效应。高变代谢物（SD>1）用于进一步分析。
对于肿瘤类器官代谢组学，小样本量12，去除未通过缺失性、mock背景或合并QC可重复性过滤器的特征。数据对数变换，归一化到合并QC中位信号，使用minimum-minus-one方法插补。鉴于小样本量，代谢物额外按变异性（SD>0.5）优先排序。

Bulk RNA测序（RNA-seq）：

从FFPE组织切片提取RNA，使用RNeasy FFPE试剂盒。RNA浓度使用NanoDrop 2000分光光度计测量，记录A260/280吸光度比。RNA质量使用Agilent 2100 Bioanalyzer系统评估。
使用QuantSeq 3' mRNA-Seq Library Prep Kit（FWD；Lexogen）生成RNA测序文库。简言之，通过oligo-dT引物在第一链合成期间从总RNA富集poly(A)+ RNA。第一链合成和RNA去除后，使用含5' Illumina兼容接头序列的随机引物进行第二链合成。所得双链cDNA使用磁珠纯化。
文库扩增通过PCR进行以掺入簇生成所需的完整接头序列。所得cDNA文库在Illumina平台上测序（HiSeq 4000单端1×51 bp reads或NovaSeq单端1×100 bp reads）。

Bulk RNA-seq分析：

使用内部生物信息学管道生成转录组谱。简言之，HISAT2用于将FASTQ测序reads比对到人类参考基因组，featureCounts用于生成基因水平原始计数矩阵。为最小化低转录本覆盖相关的潜在偏差，仅保留>100万有用映射reads的样本用于下游分析。
原始计数数据在DESeq2中使用设计公式~1处理。基因过滤保留在至少3个样本中至少5个计数的基因。应用DESeq2归一化，然后进行方差稳定变换（VST）用于下游分析。
来自多中心、文库制备时间点和/或测序泳道的批次效应使用sva R包中实现的经验贝叶斯ComBat方法校正。批次结构和校正效果最初通过PCA评估。下游分析使用批次校正的VST表达值进行，关键发现使用独立数据集进一步验证。

蛋白质组学和质谱数据分析：

原始质谱数据文件使用MaxQuant软件（版本2.0.3.0）分析LC-MS/MS蛋白质定量。默认设置用于分析，除非另有说明。包括以下参数：Trypsin/P特异性；可变修饰（甲硫氨酸氧化和蛋白质N端乙酰化）；固定修饰（半胱氨酸羧甲基化）；多重性设置为1；第一次搜索容差20 ppm；主搜索容差前体离子6 ppm和碎片离子20 ppm；数据搜索UniProt人类参考蛋白质组数据库；最小肽长度7个氨基酸；未指定标记（无标记分析）；肽和蛋白质FDR 0.01；结果使用"re-quantify"选项精炼；"match between runs"启用，1分钟保留时间窗口；无标记定量（LFQ）最小比率计数1。
所得数据进一步预过滤保留≥80%有效（非零）值的条目，然后log2变换和全局中位数中心化（归一化）。为解决低丰度蛋白质中常见的缺失值，采用左移正态分布基础插补策略。批次效应（收集中心和裂解时间差异）在下游分析前使用sva R包中的ComBat函数校正。

甲基化数据和联合RNA-seq分析：

TCGA结肠和直肠癌甲基化微阵列数据从UCSC Xena Browser下载。从443个含甲基化数据的样本中，使用Phenotype文件从同一数据源选择395个原发性肿瘤和45个正常样本。三个代表"复发肿瘤"和"转移"的样本未被选择以保持原发肿瘤一致性。相应的RNA-seq数据也从同一来源下载。
匹配的TCGA甲基化和基因表达数据从肿瘤（MSS RAS-mt n=153和MSS RAS-wt n=136）和匹配正常（n=12）样本中进一步选择。保留97,187个已发表增强子探针用于分析。
ELMER R包（v.2）用于评估基于增强子的甲基化变化，包括CpG位点的高/低甲基化、显著甲基化-基因对、富集转录因子结合基序和相关可视化。显著探针-基因对使用ELMER基于置换衍生经验p值和FDR校正识别。对于这些对，post hoc计算Pearson相关，应用效应量阈值（r<-0.2）以优先考虑生物学有意义关联。
为总结全局趋势，为每个样本计算所有显著探针-基因对的平均甲基化和匹配表达值，标准化为Z评分。平均甲基化和表达之间的关系使用线性回归和Pearson相关评估。
为可视化全基因组甲基化模式，为每个保留CpG探针计算肿瘤平均甲基化β值，然后按染色体聚合到5-Mb基因组窗口。环状基因组图在R中使用BioCircos生成。

增强子低甲基化相关基因表达的评估：

使用bulk RNA-seq数据评估增强子低甲基化与基因表达之间的关系。低甲基化基因从使用ELMER识别的显著探针-基因对定义（Pearson相关≤-0.2和FDR<0.05）。仅保留转录组数据集中存在的基因用于下游分析。对于每个基因，表达总结为所有样本的中位表达。
为确定低甲基化基因是否表达高于偶然预期，应用基于置换的框架。简言之，从表达基因背景中采样1,000个与基因集大小匹配的随机基因集，排除低甲基化基因集。对于每个迭代，计算采样基因集的中位表达值，生成中位表达值的零分布。将低甲基化基因集的观察中位表达与零分布比较，经验p值计算为置换中位值大于或等于观察值的比例。

谷氨酰胺代谢的基因和蛋白质-代谢物分析：

研究bulk肿瘤样本中谷氨酰胺代谢和下游通路活性，分析稳态代谢物丰度与谷氨酰胺合成和利用关键酶的基因表达之间的关联。具体检查GLUL（催化谷氨酰胺合成）、GLS和GOT2（分别介导谷氨酰胺分解和线粒体转氨）。
初始分析使用线性回归模型：代谢物丰度~基因表达，其中代谢物丰度对应于谷氨酰胺或谷氨酸，分子特征对应于GLUL、GLS、GOT2的基因表达或蛋白质丰度。这些模型用于建立基线关联并基于底物-产物动力学解释关系。

使用反卷积调整肿瘤内细胞异质性进行基因-代谢物关联分析：

为考虑bulk转录组数据中的肿瘤内细胞异质性，使用MuSiC估计细胞类型比例。为每个bulk样本推导正常上皮细胞、恶性上皮细胞、髓系细胞和成纤维细胞的细胞类型分数，用作下游分析中的协变量。
为在考虑细胞组成的同时评估基因-代谢物关联，拟合线性回归模型：代谢物丰度~基因表达+细胞类型分数，其中基因表达对应于GLUL、GLS或GOT2，细胞类型分数对应于一个单独的反卷积区室。模型独立拟合谷氨酰胺和谷氨酸以避免细胞类型估计之间的多重共线性。对于每个基因-协变量组合，提取基因项的回归系数、95%置信区间和p值。结果使用forest图总结。

MYC靶点甲基化-表达关联分析：

使用手动策划的MYC靶基因集研究DNA甲基化与转录活性之间的关联。对应于MYC靶点的探针-基因对基于显著负相关（Pearson r<-0.2）过滤。对于每个探针-基因对，跨样本拟合线性回归模型（甲基化~基因表达）以估计效应量和统计显著性。所得p值使用Benjamini-Hochberg方法调整多重检验，识别显著关联。顶级关联按调整p值排序，并使用回归系数条形图可视化。

iPhenMap模型：

整合表型映射（iPhenMap），先前发表的PhenMap框架的扩展，用于单组学和表型数据的联合建模，是一种贝叶斯方法，捕获跨多组学数据集的共享和数据类型特异性变异，同时将潜在变量（映射变量，MVs）与表型协变量连接。
让y_i=(y_i^(1),...,y_i(M))^{T表示来自肿瘤i的测量，跨M个组学平台（如转录组学、蛋白质组学、拷贝数改变和代谢组学）分析。每个肿瘤额外由l个表型变量表征，表示为x_i=(x_i1,...,x_iJ+1)}T，其中初始项表示截距。对于n个肿瘤（i=1,...,n），iPhenMap使用一组低维潜在映射变量z_i=(z_i1,...,z_iQ)^T建模每个数据类型t的观察：
- P(y_i^(t)) = MVN_pi[W^(t)z_i, Σ_j=1^(t)]（方程1）
- P(z_i) = MVN_q[βx_i, Φ]（方程2）
在此公式下，每个数据类型通过载荷矩阵W^(t)∈R(p1×q)连接到潜在变量，将特征映射到q个潜在MVs。残差变异由对角协方差矩阵Σ^{(t)=diag(σ_1t}2,...,σ_p1t^{2)捕获。表型协变量与潜在变量之间的关联由回归系数β控制，而潜在空间中的残差变异性由Φ=diag(φ_1}2,...,φ_q^2)描述。
先验分布：为捕获跨和组学数据集内的结构化依赖性，iPhenMap对载荷矩阵W^(t)施加稀疏诱导先验。分层公式结合spike-and-slab先验以确定潜在变量是否在给定数据类型中活跃，以及自动相关性确定（ARD）先验以强制执行特征水平稀疏性。
为进一步调节表型协变量的影响，对回归系数β放置g先验。该先验稳定估计并降低对协变量缩放差异的敏感性。
由于这些回归系数是MV特异性的，它们捕获协变量对潜在结构的共享和数据集特异性效应。模型参数在贝叶斯框架内推断，非常适合基因组学典型的有限样本量的高维设置。
工具可在https://github.com/syspremed/iPhenMap获取。

基于MOFA的多组学分析和与iPhenMap的比较：

COLOSSUS St2/3 CRC（n=98）多组学数据使用与iPhenMap一致的框架处理。归一化转录组和蛋白质组使用MOFA整合以推导捕获协调分子变异的低维因子。
由于MOFA无法同时建模因子和临床变量，这些因子随后使用多变量分析在第二步中与临床变量关联。每个因子的贡献用于解释 underlying 分子程序。为与iPhenMap比较，MVs和MOFA衍生因子使用Pearson相关分析相关。

iPhenMap和MOFA作为互补模型：

iPhenMap采用完全贝叶斯Gibbs采样框架，具有混合（半监督）策略，其中潜在映射变量（MVs）从多组学数据推断，同时通过模型内回归由临床和生物学协变量引导，使分子模式能够在模型训练期间直接与表型特征连接。
MOFA是一种广泛使用的无监督变分贝叶斯框架，将多组学变异分解为仅从分子数据衍生的潜在因子，使得跨多样数据类型的高效建模成为可能。
这些框架提供对共享生物信号的互补视角。因此，应用两种方法整合COLOSSUS St2/3多组学数据集，在MOFA中通过潜在因子回归 post hoc 纳入协变量。值得注意的是，两种模型产生广泛一致的结果。

TCGA验证iPhenMap COLOSSUS结果的MV1：

为外部验证来自iPhenMap COLOSSUS St2/3队列的MV1，TCGA CRC匹配转录组、蛋白质组和临床数据（n=64）使用iPhenMap框架整合，类似于COLOSSUS St2/3多组学数据。
临床协变量被策划和数值编码，包括肿瘤侧别、病理分期、性别和淋巴管侵犯（该队列中PNI缺失）。RNA和蛋白质特征通过中位数绝对偏差过滤以保留最高变特征，然后TCGA数据集限制为与先前识别的COLOSSUS St2/3 CRC MV1相关载荷重叠的基因和蛋白质。
多组学数据然后与临床协变量联合建模（在单一步骤中）使用iPhenMap推导捕获协调分子和表型变异的MVs。所得元变量（包括MV1）与临床特征的关联通过模型内估计的回归系数评估，以验证COLOSSUS St2/3 CRC队列的结果。

代谢物-表型-蛋白质相互作用机器学习方法：

为表征患者样本中与表型（PNI）相关的代谢变化，对来自COLOSSUS St2/3 MSS RAS-mt CRC队列（n=161）的肿瘤样本应用基于重采样的线性回归分析。样本随机分为训练和测试集（70/30），分析超过1,000次迭代以解决肿瘤异质性。
代谢物特征首先从iPhenMap潜在变量中优先排序，选择在MV5中绝对缩放载荷系数大于2的注释代谢物，排除未知代谢物。与表型一致相关的代谢物使用20%（复发信号）和10%（更可变模式）频率阈值识别。
为评估表型相关代谢变异如何影响蛋白质表达，使用带交互项的线性回归（蛋白质丰度~PNI+代谢物+PNI×代谢物）为每种蛋白质建模表型状态（PNI）和代谢物水平之间的交互效应。
使用重复分层训练-测试重采样（100次迭代）评估稳健性，样本分为70%训练和30%测试集，同时保留类别平衡表型。在每个迭代中，交互项系数在训练和测试集中独立估计。如果任一分割中完整观察少于5个、缺乏两个表型类别、未能收敛、产生非有限估计或表现出退化拟合（残差自由度<0或R²>0.9999），则模型被排除。
结果通过计算以下指标跨迭代总结：（i）训练集中名义显著交互效应的频率，（ii）测试集中名义显著交互效应的频率，（iii）训练和测试集之间显示一致交互方向的迭代比例。方向验证的交互定义为训练中名义显著且测试中效应方向一致，无论测试集显著性如何。
主要候选交互基于可重现方向性优先排序，定义为训练显著性频率≥0.6和方向验证频率≥0.6，并按验证频率、训练频率和中位训练集p值排序。对于下游生物学解释，对优先交互的蛋白质进行富集分析。正方向交互进一步基于跨迭代中位交互估计为正值定义。
蛋白质测量来自COLOSSUS St2/3队列的LC-MS/MS（n=161）和RPPA（n=121）数据集。

通路富集分析：

来自不同实验的基因/蛋白质使用Enrichr、hypeR、ssGSEA、基于GSVA的ssGSEA、singscore、基因集平均评分或ClusterProfiler及下述相应基因集进行基因富集分析。
MetaboAnalyst 6.0（基于网络的工具；使用前将代谢物名称转换为HMDB身份）用于基于代谢物的通路和网络分析。
用于富集分析的基因集包括MSigDB的Hallmark和curated Wong_Proteosome_Gene_Module、KEGG（包括MMR签名）、REACTOME和Biocarta、Metabolomics Workbench、复制应激基因集、ENCODE TF ChIP-seq 2015、WikiPathways和TG-GATES。
基于snRNA-seq/scRNA-seq的富集分析使用Seurat和escape R包进行。

免疫反卷积：

使用MCP-counter工具对bulk RNA-seq数据进行免疫反卷积。

肿瘤类器官和PDX模型RNA-seq数据分析：

来自Leto等人的RAS-mt（n=48）和RAS-wt肿瘤类器官（n=83）的归一化RNA-seq表达数据跨重复模型平均用于进一步分析。
类似地，来自Leto等人的PDX样本的归一化RNA-seq表达数据过滤保留未治疗肿瘤类器官，共享相同模型标识符的技术/生物学重复样本通过平均表达值合并。基因符号在进一步分析前清理，最终PDX队列（n=441）。
由于肿瘤类器官在测试特征中显示有限差异，PDX未基于RAS或微卫星状态选择；MSI病例不常见。这些肿瘤类器官的突变也从Leto等人获得。

使用靶基因面板表达指标选择肿瘤类器官模型：

分析来自MSS RAS-mt肿瘤类器官的RNA-seq表达数据以指导模型选择。来自n=131模型（合并重复后）的可用数据用于评估预定义靶基因面板（PTK2、PDK1、HIF1A、GLS、SRC和AURKB）和蛋白酶体基因签名的变异性。
对于每个样本，跨6个选定靶基因计算面板水平汇总指标，包括表达的均值、总和、SD、变异系数（CV；SD/均值）和范围（最大值-最小值）。重复合并为每个病例的单一代表值，取每个指标的中位数，同时保留重复变异性作为最小-最大范围。
为评估生物学异质性，CRC病例对每个指标独立分层，使用数据驱动的百分位数阈值。具体而言，CRC水平中位值分布的第33和第67百分位数用于定义低、中和高指标层。选择9个MSS RAS-mt肿瘤类器官，使其靶面板指标跨多个指标跨越这些层，而非代表单一极端表达状态。
所有样本保留在分析中，明确评估重复一致性以确保观察到的差异反映生物学而非技术变异。最终模型选择还考虑实验可管理性，包括生长动力学和反复冻融循环的耐受性。

量化和统计分析

Wilcoxon秩和检验用于两组之间的非参数比较，包括代谢物丰度比较，而Kruskal-Wallis检验用于多组比较，包括单核细胞类型富集分析。
多重检验校正使用FDR和Benjamini-Hochberg（BH）程序（如适用）。应用于个体分析的统计检验和校正在图例中指定。
对于肿瘤类器官研究，显著性使用单因素方差分析与Sidak多重比较校正确定（使用原始数据）。
对于PDX生存研究（Kaplan-Meier曲线；长期效应），显著性使用整体和成对log rank检验确定。
对于PDX肿瘤体积研究（短期效应），显著性使用单因素方差分析与Tukey多重比较校正确定。