导语

今天给同学们分享一篇结肠癌+分型+预后的生信文章“Comprehensive genomics analysis of aging related gene signature to predict the prognosis and drug resistance of colon adenocarcinoma”，这篇文章于2023年2月8日发表在Front Pharmacol期刊上，影响因子为5.6。

结肠癌（COAD）是一种异质性肿瘤，衰老在癌症发生中至关重要。本研究旨在确定基于衰老的亚型，并构建一个预后标志来预测COAD患者的预后，并指导免疫治疗或化疗决策。

(资料图)

1.结肠癌（COAD）的单细胞RNA测序分析和标记基因识别

作者基于dim = 35对细胞进行了聚类，并获得了17个细胞亚群。图1A展示了腺瘤、血液、癌症、正常和癌旁样本的t-SNE图。图1B描绘了聚类后的17个细胞亚群。然后，作者使用一些典型的免疫细胞标记物对细胞亚群进行了注释。补充图S2提供了标记基因表达的TSNE图。图1C展示了注释细胞亚群的聚类特征。亚群2、3、5、7、9、11和12是表达EPCAM的上皮细胞；亚群0和6是表达CD3D、KLRD1和CD8A的自然杀伤（NK）T细胞；亚群1和10是表达MS4A1的滤泡性B细胞；亚群8是表达MZB1的浆细胞B细胞；亚群4是表达CD68、CD14和FCGR3A的单核细胞源性巨噬细胞（MDMD）；亚群13是表达DCN的成纤维细胞；亚群15和16是表达KIT的肥大细胞。因此，作者统计了每个亚群的细胞数量，并计算了其在不同样本中的比例（图1D）。上皮细胞和NK T细胞的细胞数量比其他亚群更多。图1E显示了每个亚群中的前5个标记基因，富集分析显示标记基因与人类T细胞白血病病毒1感染、Th17细胞分化、造血细胞系谱和Th1和Th2细胞分化密切相关。

图1 细胞亚群和标记基因的筛选

2.肿瘤微环境中的细胞衰老特征

为了表征单个细胞的肿瘤微环境中的细胞衰老，作者区分了细胞亚群中的非整倍体和二倍体。结果显示，在恶性细胞中有12,362个非整倍体细胞，30,833个二倍体细胞，它们的TNSE图在图2A、B中显示，表明癌组织中恶性细胞更多，而癌旁组织中恶性细胞较少。进一步，作者使用ssGSEA方法计算了恶性细胞和非恶性细胞中与细胞衰老相关的通路得分。结果显示，恶性细胞中的细胞衰老相关通路得分高于非恶性细胞（p < 0.0001）（图2C）。

图2 单细胞肿瘤微环境中的细胞衰老特征

3.基于批量RNA测序数据验证了细胞衰老

为了进一步验证细胞衰老的特征，作者基于批量RNA测序数据评估了肿瘤和癌旁组织中与衰老相关的通路。GOBP_REPLICATIVE_SENESCENCE、REACTOME_CELLULAR_SENESCENCE、REACTOME_DNA_DAMAGE_TELOMERE_STRESS_INDUCED_SENESCENCE和KEGG_P53_SIGNALING_PATHWAY在TCGA队列的肿瘤组织中显著富集（图3A、B）。通过使用ssGSEA方法计算衰老相关通路得分，发现RGOBP_REPLICATIVE_SENESCENCE、REACTOME_CELLULAR_SENESCENCE、REACTOME_DNA_DAMAGE_TELOMERE_STRESS_INDUCED_SENESCENCE和KEGG_P53_SIGNALING_PATHWAY在肿瘤组织中的衰老得分高于癌旁组织（p < 0.001）（图3C）。

图3 基于批量RNA-seq数据的细胞衰老验证

4.鉴定出了三种老化亚型

作者使用上述四个富集通路中的基因进行了单变量Cox回归分析。共鉴定出与预后相关的16个基因（图4A）。为了进一步确定亚型，作者对432个TCGA-COAD样本的表达谱进行了一致性聚类分析。根据16个衰老相关基因的表达谱，聚类结果显示，聚类=3具有相对稳定的聚类效果（图4B）。考虑到一致性矩阵k=3是一个较好的选择，作者将整个队列分为三个亚型（图4C）。接下来，Kaplan-Meier曲线显示三个亚型之间存在显著差异，其中clust3的生存概率最低，而clust1在TCGA中具有最佳预后（p=0.0021）（图4D）。类似的结果也在GSE39582队列中观察到（p=0.0027）（图4E）。

图4 鉴定出了三种衰老亚型

5.Clust3的预后最差，T分期较高

作者随后在TCGA队列中比较了三个亚型的临床病理特征分布（性别、T分期、N分期、M分期、分期、年龄和生存状态）。结果发现三个亚型在T分期和生存状态上存在显著差异（图5A）。此外，桑基图详细展示了三个亚型在T分期和生存状态上的分布情况（图5B）。clust3患者预后最差，T分期较高（主要为T3和T4期）。

图5 Clust3的预后最差且T分期较高

6.Clust3展示了更高的TMB和突变

此外，突变特征在衰老亚型中进一步评估。拷贝数变异在三个亚型中发生了显著变化（图6A）。同时，SNV结果显示APC（69%），TP53（51%），TTIN（49%）和KRAS（39%）在前15个突变基因中表现出较高的突变频率（图6B）。此外，与clust1相比，clust3中的TMB和突变基因数量均有所增加（图6C，D）。

图6 Clust3展示了更高的TMB和突变

7.衰老亚型与上皮-间质转化、缺氧、血管生成和与肿瘤相关的途径相关联

在肿瘤中，细胞衰老促进细胞外基质的裂解，导致生长因子的释放，从而促进上皮-间质转化（EMT），进而导致肿瘤转移。因此，作者阐明了衰老亚型与EMT评分之间的关系。EMT评分在衰老亚型之间有明显的差异，其中clust3的EMT评分最高，相比clust1和clust2（图7A）。同时，缺氧和血管生成评分在clust3中高于clust1和clust2（图7B，C）。图7D发现，在这三个亚型中，有9个与肿瘤相关的通路发生了显著改变，包括细胞周期、HIPPO、MYC、NOTCH、NRF1、PI3K、TGF-beta、RAS、TP53和WNT。

图7 EMT、缺氧、血管生成和与肿瘤相关的通路之间的衰老亚型关联

8.衰老亚型与免疫特征之间的关系

免疫浸润评分，包括基质评分（StromalScore）、免疫评分（ImmuneScore）和ESTIMATEScore，在三个亚型之间存在显著差异（图8A），作者发现clust3的免疫浸润程度更高。三个亚型之间的免疫细胞浸润也有显著变化（图8B）。此外，clust3在几个与炎症相关的通路上得分更高，包括JAK-STAT信号通路、NF-Kappa B信号通路、Toll样受体信号通路、B细胞受体信号通路、T细胞受体信号通路和炎症反应（图8C-H）。此外，clust3的TIDE评分更高（图8I），表明clust3更容易发生免疫逃逸。

图8 衰老亚型与免疫特征的关联

9.衰老基于风险模型的构建和验证

通过对三个亚型进行差异分析，鉴定出2,085个差异表达基因（DEGs），并用于单变量Cox回归分析。选择了对预后影响更大的194个基因，其中包括180个风险基因和14个保护基因（图9A）。为了减少基因数量，进行了LASSO Cox回归分析。随着lambda的逐渐增加，独立变量系数趋近于零的基因数量逐渐增加（图9B）。采用10折交叉验证，并在每个lambda下显示置信区间（图9C）。当lambda = 0.0347时，选择了25个基因进行进一步分析。基于AIC的逐步多元回归分析最终确定了11个基因（MFNG、GPRC5B、TNNT1、CCL22、NOXA1、PABPC1L、PCOLCE2、MID2、CPA3、HSPA1A和CALB1）。随后，在TCGA队列中进行的生存分析显示，高风险患者的预后较低风险患者差（p < 0.0001），1年AUC为0.81，3年AUC为0.77，5年AUC为0.75（图9D）。为了验证其稳健性，还在GSE39582和GSE17537队列中进行了生存分析。高风险患者在GSE39582（p < 0.0001）和GSE17537（p = 0.006）中的预后较低风险患者差，预后预测表现良好（图9E，F）。

图9 构建和验证基于衰老的风险模型

10.风险评分与临床病理特征和生物特征的关联

为了进一步阐明风险评分与临床病理特征之间的关系，作者比较了TCGA队列中临床病理特征在风险评分上的差异，并发现临床分期较高（T分期、N分期、M分期和分期）的患者风险评分较高（图10）。此外，clust3患者的风险评分较高（图10）。为了评估风险评分与生物学特征之间的关系，作者进行了风险评分与缺氧、血管生成和EMT评分的相关性分析。图11A-C显示，缺氧、血管生成和EMT评分与风险评分均呈正相关。接下来，作者进行了基于衰老的风险评分与底层调控的KEGG通路之间的相关性分析，以找到与风险评分相关的通路（图11D），进一步的统计结果显示，只有一些通路在高风险组和低风险组之间存在显著差异（图11E）。此外，作者比较了高风险组和低风险组之间与衰老相关通路的评分。结果显示，高风险患者在GOBP_AGING、GOBP_MUSCLE_ATROPHY和GOBP_NEGATIVE_REGULATION_OF_CELL_AGING方面得分较高；而低风险患者在GOBP_STRESS_INDUCED_PREMATURE_SENESCENCE、KEGG_P53_SIGNALING_PATHWAY和REACTOME_DNA_DAMAGE_TELOMERE_STRESS_INDUCED_SENESCENCE方面得分较高。

图10 不同临床病理特征中风险评分的分布

图11 风险评分与生物特征的关联

11.化疗反应预测

此外，作者评估了高风险组和低风险组对传统化疗药物的反应性。如图12所示，低风险患者的厄洛替尼（p < 0.001）、舒尼替尼（p < 0.001）、MG-132（p < 0.001）、CGP-082996（p < 0.01）、AZ628（p < 0.01）、索拉非尼（p < 0.001）、VX-680（p < 0.01）和Z-LLNle-CHO（p < 0.01）的IC50值显著低于高风险患者，这表明低风险患者对厄洛替尼、舒尼替尼、MG-132、CGP-082996、AZ628、索拉非尼、VX-680和Z-LLNle-CHO更敏感。

图12 化疗反应性预测

12.风险基因表达与甲基化水平之间的关系

此外，作者在TCGA中分析了11个风险基因的甲基化水平，并对风险基因表达和甲基化水平进行了相关性分析。图13显示，CALB1、CPA3、NOXA1和TNNT1的表达与它们的甲基化水平呈显著正相关；而CCL22、GPRC5B、HSPA1A、MFNG、PABPC1L和PCOLCE2的表达与它们的甲基化水平呈显著负相关。

图13 对风险基因表达和甲基化水平进行了相关性分析

13.Nomogram构建

基于单变量和多变量Cox回归分析的基础上，结果显示风险评分是一个独立的预后因素（图14A）。为了更好地量化COAD患者的风险评估和生存概率，作者构建了一个用于估计1年、3年和5年生存率的Nomogram，使用了风险评分和临床病理特征。图14B显示风险评分对COAD患者的生存率影响最大。该Nomogram的校准曲线显示观察值和预测值之间具有高度一致性（图14C）。为了评估该模型的可靠性，进行了DCA分析，并确认Nomogram和风险评分在预测预后方面具有最强大的能力（图14D）。

图14 COAD患者预测预后的Nomogram构建

总结

总之，作者开发了能够区分结局、T分期、突变和免疫特征的衰老基础亚型，这可能指导进一步研究结肠癌异质性的机制。此外，作者构建并验证了一个基于衰老的标志，并为结肠癌患者的预后预测提供了可靠的工具。