简单易学5分+单细胞挖掘套路！

Integrated analysis of single‑cell RNA‑seq and bulk RNA‑seq unravels tumour heterogeneity plus M2‑like tumour‑associated macrophage infltration and aggressiveness in TNBC单细胞RNA-seq和bulk RNA-seq整合分析揭示了TNBC中肿瘤异质性以及M2样肿瘤相关的巨噬细胞浸润和侵袭

一、 研究背景

三阴性乳腺癌（TNBC）的特征是缺乏孕酮受体 、雌激素受体和人表皮生长因子受体2(HER2)，与其他乳腺癌类型相比，TNBC表现出更多的血管生成和上皮-间充质转换(EMT)，且具有广泛的肿瘤间和肿瘤内异质性，在单细胞水平上鉴定TNBC细胞的异质性和高转移性肿瘤亚群将有助于TNBC的精确治疗。

二、 分析流程

三、 结果解读

1.TNBC细胞中肿瘤细胞的异质性

作者整合了两个单细胞转录组（GSE118389和GSE75688），用Seurat包基于regularized negative binomial regression的方法处理批次效应以及进行单细胞RNA-seq分析，采用UMAP方法进行非线性降维（图1.B）,显示出每个患者的肿瘤细胞高度聚集。而免疫细胞则根据其细胞类型而不是来源表现出较高的聚类（图2.A）。

TNBC的肿瘤复发与几个侵袭性癌症特征密切相关，如细胞干性、血管生成和EMT，因此作者在单细胞水平上分析了这三个特征。在TNBC细胞中，细胞干性、血管生成和EMT互相呈正相关（图1.C-E），这三个特征的相关性表面TNBC的转移性更强，复发率更高。作者还鉴定出三者分数都很高的细胞(散点图中以红点标出)，提示它们在肿瘤进展、转移和复发中可能具有重要作用。

表1.1TNBC患者信息(单细胞RNA-seq)

图1.在单细胞水平上TNBC肿瘤内的异质性

2.TNBC细胞中免疫细胞的异质性

免疫细胞经UMAP降维被分为三类：T细胞、B 细胞和巨噬细胞(图.2A)，热图(图.2B)显示了三个子类的代表基因。接着作者分别对三类免疫细胞进行分析。

对于T细胞，作者根据细胞在每个样本中的GSVA富集分数进行聚类分析，将其与T细胞功能状态如调节性、共刺激性、初始性、细胞毒性和耗竭性等联系起来，分别表现出初始T，共刺激T，调节T和耗竭性T细胞的基因表达特征。与其他患者相比，BC09患者具有更多的细胞毒性T细胞。相反，PT058和 PT084患者的细胞毒性T细胞活性水平较低。高细胞毒活性的T细胞表达高水平的CST7、GZMA、GZMB、IFNG、NKG7、PRF1 和TNFSF10。另外，作者发现TNBC组织中的T细胞即使在同一患者中也表现出混合的初始T细胞状态，提示初始T细胞状态可能与TNBC中T细胞细胞毒性活性的效率有关。(图.3)

对于肿瘤浸润的B细胞，作者分析了B细胞的功能状态， 如抗凋亡、初始记忆、细胞因子、促凋亡、增殖、生发中心和相关基因表达特征。在TNBC患者之间和同一TNBC患者中，B免疫细胞及相关基因表达均表现出异质性分布。(图.4)

对于肿瘤浸润的巨噬细胞，作者在巨噬细胞中发现了较高的M2活性。在许多癌症类型中，M2样TAMS通过诱导血管生成、EMT和细胞干性以及抑制免疫监测来促进肿瘤的发生。在TNBC的大多数巨噬细胞中，M2型基因如TGFBI、MS4A6A和CD163被上调。此外，已知促进血管生成、EMT和肿瘤进展的基因，如IL8和PLAUR，在TNBC的巨噬细胞中过表达。总的来说，M2样TAMs是TNBC巨噬细胞的主要类型。

图2.肿瘤浸润免疫细胞群的聚类

图3.在单细胞水平上TNBC中的T细胞特征

图4.在单细胞水平上TNBC中的B细胞特征

图5.在单细胞水平上TNBC中的免疫浸润细胞整合分析

3.M2样TAMS与TNBC患者生存不良有关

作者用CIBERSORT计算了TCGA-TNBC患者bulk RNA-seq转录组图谱中M2样TAMS的丰度。表2总结了TCGA-TNBC患者信息。生存分析表明，高丰度的TAMs与 TNBC患者较差的无复发生存(RFS)有关(图.6.A)。为了进一步了解TAMs在TNBC中的潜在作用，作者进行了WGCNA分析，以软阈值=5，拟合指数R2=0.95构建无标度的共表达网络来识别与TAM相关的基因特征。一共生成7个模块，其中蓝色模块与与M2样TAM得分相关性最高(r=−0.48，P=1e−6，图.6.B)。图.6.C中的每一个点代表一个基因，横坐标Module Membership表示基因与模块（module eigengene）的相关性，纵坐标基因显著性表示基因与TAM的相关性，这里可以看出与TAM高度显著相关的基因往往是与TAM显著相关的模块中的重要元素。接着，作者提取了蓝色模块的hub genes并进行单因素COX回归分析，以p=0.01为阈值，火山图显示筛选得到146个基因(图.6.D)，并使用TAM相关基因计算预后风险评分。

作者基于PCA以146个基因的第一主成分计算了保护得分Anti-TAMPCA(图.6.E) ，具体计算过程见下图，其中Eij表示肿瘤样本j中第i个hub gene的 log2(RSEM + 1) 表达量，Ci表示表示TAM细胞相关基因对应的coefficient。在TCGA-TNBC队列中，高保护得分的患者比低保护得分的患者有更好的预后 (HR=0.11，P<0.001，图.6.F)。作者还在GEO整合数据集(GSE19615, GSE21653, GSE31519)中进一步验证了TAM相关基因预测TNBC患者RFS的能力(HR=0.56，P=0.033，图.6.G)。接着，作者使用ssGSEA鉴定了TAM相关基因特征的潜在机制，图.6.H显示胆汁酸代谢、雌激素反应和P53通路与TAM相关的保护评分呈较高的负相关。

表2.TCGA-TNBC患者信息

图6.bulk RNA-seq分析显示TNBC患者高TAMs与预后不良相关

4.根据与TAM相关的基因进行分子分型

根据TAM相关基因的表达，使用ConsensusClusterPlus包基于K均值的Consensus Clustering(一致性聚类)将TNBC分为不同的亚组。一致性聚类通过基于重采样的方法来验证聚类合理性(即找到一个适合的K值)。常见标准是选择CDF(Cumulative distribution function)(图7.A)下降坡度小的K值，作者在这里选择了K=2作为最优参数，将训练队列划分为不同的亚组。然后图7.B建立一致性矩阵(consensus matrix)，一致性矩阵的数值在[0,1]，等于1代表多次聚类两个数据点全部在同一个类里面，等于0代表多次聚类全部不在同一个类里，图中显示出K=2时聚类效果较好。图7.C热图显示了作者使用limma包筛选得到的两个cluster之间的差异表达基因。图7.D显示cluster 1患者的预后优于cluster 2患者，图7.E的小提琴图显示了与cluster 1相比，cluster 2的TAM评分显著升高(P<0.005)，总的来说，低TAM评分的cluster 1患者预后优于高TAM评分cluster 2患者。

作者还进一步进行了GSEA分析，与cluster 2相比，cluster 1上调的通路包括与雌激素早期反应相关的通路、与雌激素晚期反应相关的通路、M2巨噬细胞的上调以及p53通路；cluster 1下调的通路包括E2F靶点、G2M检查点和MYC靶点相关通路。(图7.F-G)

图7.根据TAMs相关基因特征进行分子亚型识别

5.通过随机森林算法和SVM算法对特征重要性进行排序， 揭示与M2样TAM评分相关的最重要基因

首先，作者通过ranger包找到回归过程中的最佳超参数（这里的超参数即为机器学习算法中的调优参数，如正则化系数λ等）。然后，使用randomforest包构建随机森林回归模型，使用e1071包构建SVM模型并进行特征排序，提取了前50个基因，随机森林算法和SVM算法有19个基因重叠，作者发现其中的VIPR1、ABCD4、SLC40A1和SCNN1A是不良的预后因素，这四个基因在之前的研究中显示与代谢变化有关，这提示了TAM与TNBC代谢变化之间的潜在联系。而其他15个基因与较好的RFS有关。

6.与TAM相关的基因特征预测免疫治疗反应

作者接着在Python中采用PyTorch构建神经网络，通过TAM相关基因来预测免疫治疗的反应。神经网络的构造过程如图8.A所示。神经网络由输入值x，权重w，偏置b，激活函数构成。作者使用了随机梯度下降法求得最合适的权重，使得预测值与真实值之间的误差最小，梯度下降的learning rate或步长设置为0.001。将TCGA-TNBC队列作为训练集，将整合的GEO队列作为测试集来评估准确性。在训练过程中，dropout rate设置为0.2，dropout是指构造神经网络前向传播的时候，让某个神经元节点以一定的概率p隐藏，这样可以使模型泛化性更强，不会过于依赖某些局部特征，减少过拟合。激活函数为ReLU。训练集的损失值随迭代次数Epoch的增加而减小，表明训练过程的效率很高(图8.B)。当对训练集迭代次数达到1000时，值得注意的是，使用混淆矩阵在测试集中正确地识别了所有样本(图8.D)。在ROC图中，曲线下面积(AUC)达到100%(图8.C)。

图8.基于神经网络的深度学习框架构建与验证

小结

本篇文章中使用公共数据库的scRNA-seq和bulk RNA-sec数据，分别分析了三阴乳腺癌肿瘤细胞和肿瘤免疫浸润细胞的肿瘤间和肿瘤内异质性。在肿瘤细胞中，在单细胞水平上发现了血管生成、stemness和EMT之前的高度相关性，并发现了三者同时高表达的细胞子集。在肿瘤免疫浸润细胞中，比较了T细胞、B细胞、巨噬细胞在单细胞水平上的肿瘤异质性，还发现了巨噬细胞中M2样TAMs与预后不良有关，并基于已识别的M2样TAMs相关特征建立神经网络，以预测TNBC对免疫治疗的反应。

值得注意的是，神经网络在测试集中达到了100%的准确率，这一结果强调了M2样TAMs在TNBC治疗中对预测免疫治疗的重要性，并提示M2样TAMs的消融可能对治疗肿瘤有效。然而，巨噬细胞靶向方法具有全身性毒性，因为靶体作用于所有类型的巨噬细胞，因此还需要大量的临床研究与实验。