今天给大家分享的是一篇今年9月份发表在《Frontiers in Oncology》（IF:6.244）上的文章。目前大量证据表明，卵巢癌 (OC) 的异质性是其治疗失败的主要原因。单细胞RNA测序技术 (scRNA-seq) 是在单细胞水平上分析肿瘤异质性的有力工具，有助于在遗传和细胞水平上更好地了解细胞功能。本研究中，作者将scRNA-seq与bulk RNA-seq相结合，最后确定了两个标记基因，在OC中具有重要的预后意义和免疫治疗价值。

scRNA-Seq与Bulk RNA-Seq整合分析卵巢癌免疫细胞异质性并建立分子风险模型

背景

卵巢癌（OC）是一种常见的妇科恶性肿瘤，死亡率很高，并具有异质性。目前主要的治疗手段是手术和化疗，且其异质性是其治疗效果不佳的主要因素之一。超过90％的卵巢恶性肿瘤被归类为上皮性卵巢癌，包括五种主要类型：高级别浆液性（70~80%％），低级浆液性（<5％），粘液性（3％），子宫内膜样癌（10％）和透明细胞癌（10％）。

数据组成和方法思路

1. 数据组成：GEO数据库中的两个scRNA-seq数据集（GSE154600和GES158937）、TCGA中OC的bulk RNA-seq。

2. 方法思路：首先使用两个scRNA-seq数据集（GSE154600 和 GES158937）来表征OC 异质性，共聚类出20个细胞群。然后根据免疫细胞标记筛选出免疫相关的细胞群，包括T细胞或NK 细胞、B细胞和髓样细胞等。接下来，便是对TCGA-OC的bulk RNA-seq数据的分析，包括有免疫浸润、WGCNA 分析、亚型分析、OS相关基因鉴定、预后风险模型构建。最后的实验验证部分包括qPCR和IHC 分析。

结果

1. ScRNA-seq聚类分析

研究人员首先对下载的scRNA-seq数据集去批次（利用Seurat包中的SCTransform() 函数），对其进行降维。再通过函数FindCluster ()对细胞进行聚类，获得20个聚类。然后再根据标记基因的表达情况注释免疫相关细胞类群，如PTPRC是免疫细胞标记，EPCAM是上皮细胞标记，COL1A2是成纤维细胞标记（）。

2. 免疫细胞分析

在这一部分，首先是对前面注释的免疫细胞群基于免疫细胞标记基因进行亚群聚类和注释。然后是分别对T细胞、B 细胞和髓样细胞进行GSVA富集分析，并描绘各类的细胞相应的特征情况（）。其中，对于肿瘤浸润的髓样细胞，分析了M1-like和M2-like髓样细胞的活性，结果发现M1和 M2相关基因在具有数据集GSE154600的P3和P4患者中显著上调。综上所述，研究人员通过分析两个OC scRNA-seq数据集以及基于免疫细胞亚型的肿瘤和髓样细胞之间的差异相互作用探索了卵巢癌的异质性。

3. CIBERSORT

CIBERSORT是一种常用的计算免疫细胞浸润的方法，它利用线性支持向量回归的原理对免疫细胞亚型的表达矩阵进行去卷积，来估计免疫细胞的丰度。在这里，为了进一步分析OC浸润免疫细胞的临床意义，研究人员通过CIBERSORT这种计算方法，根据TCGA中的bulk RNA-seq数据的counts预测22个免疫细胞的比例，同时也计算了这378名 TCGA-OC患者的bulk RNA-seq数据中M1-like TAMs（肿瘤相关巨噬细胞）丰度（A）。生存分析结果表明，M1-TAMS丰度高的患者生存率更高（B、C）。而不同丰度M2-like TAMS的患者之间的生存差异并不显著。因此研究人员后续又对M1-like TAMs进行了深入的分析。

4. WGCNA分析和免疫治疗预测

为了进一步探索M1-TAMS在OC中的潜在作用，研究人员对TCGA中的样本数据进行了WGCNA分析。最后总共生成了7个模块（A-D），其中棕色模块（内含3213个基因）与M1-like TAMs score的相关性最高（E）。为了筛选出这一模块中与OS相关的基因，他们通过计算基因与模块特征基因之间的相关性与基因表达与 OS 之间的相关性这两者之间的相关性，最后得到了45个关键基因（MM>0.5 和 GS>0.1）（F）。

根据IMvigor210CoreBiologies包中OC的免疫治疗队列数据，可将样本分为“免疫浸润型” 、 “免疫排斥型”和 “免疫沙漠型”三种表型。根据三种表型的预测AUC值，可以发现这些关键基因在用于区分免疫浸润型和免疫沙漠型样本时预测效果更好（图 5H）。这些结果表明鉴定出的这45个基因是免疫浸润的潜在预测标记。

5. 基于M1-Related基因的分子分型

筛选出关键基因以后该做什么呢？一般套路主要是有两种，一种是对样本进行亚型分类，一种是构建回归模型。这里一共筛选出了有45个基因，构建回归模型的话显然基因数量过多了，所以研究人员接着基于前面所筛选出的基因表达情况将所有样本分成两个亚型，分别研究这两个亚型之间的差异，以进一步识别生存相关标记基因。首先是提取这45个基因的表达量矩阵，然后利用R包NMF基于M1相关基因将TCGA中的OC样本分成两个亚型（cluster1和cluster2），C中就展示了这45个基因在所有样本中的表达情况。预后分析显示，cluster1的预后效果要更差（D）。两种OC分子亚型中cluster2中M1-like髓样细胞的比例更高（E）。最后还分别对两个亚型进行了KEGG富集分析（F-G）。

6. 构建预后风险模型

研究人员将TCGA中的OC样本按照1：1比例随机分成训练集和测试集，从658个DEGs挑选出101个蛋白编码基因，再通过cox回归和Lasso回归构建了预后风险模型（）。首先使用Cox识别出4个生存相关基因（CXCL13、PLA2G2D、IL26、CARD17），然后再结合Lasso回归，发现CXCL13和 IL26的低表达是一个风险因素。最后确定将这两个基因用于构建风险模型，结果发现CXCL13和IL26基因表达较低的患者的风险评分和死亡率越高，高风险组的预后也更差。

7. 预后风险模型的验证

构建完风险模型之后自然是对该模型的验证了。根据测试集数据和TCGA中所有OC样本数据分析，结果发现具有高 RiskScore的样本的风险评分明显小于具有低RiskScore的样本。同时，KM曲线显示低风险组和高风险组之间存在显著的生存差异（）。

8. 实验验证

这篇文章并不完全的纯生信分析，研究人员最后还是做了一些实验验证的。他们通过qPCR（A、B）和IHC（C、D）分析了 OC 患者临床样本中特征基因（CXCL13和 IL26）的表达情况，结果显示OC组织中CXCL13和 IL26的表达确实较低。

小结

这篇文章所用的分析数据都源于公共数据，只是在分析出关键基因后做了一点简单的实验验证。我们来回顾一下，其实文章后面基于bulk RNA-seq的分析都是常规套路，不过他们想到了和肿瘤的异质性联系在一起，所以加入了scRNA-seq数据用来表征异质性，这也正是文章的一大亮点。

参考文献

Liang L, Yu J, Li J, Li N, Liu J, Xiu L, Zeng J, Wang T, Wu L. Integration of scRNA-Seq and Bulk RNA-Seq to Analyse the Heterogeneity of Ovarian Cancer Immune Cells and Establish a Molecular Risk Model. Front Oncol. 2021 Sep 21;11:711020.