跟着小鱼头学单细胞测序-如何整合公共单细胞RNA数据？

导语

GUIDE ╲

随着单细胞研究的火热，公共的单细胞RNA-seq数据也迅速增加。虽然随着技术的进步，单细胞实验逐渐普遍，但其实验成本还是相对偏高。那么如何利用公共的单细胞RNA数据进行整合分析来挖掘更多的信息也是生信研究的一个重要方向。我们在这个专题会跟大家分享一些此方向的文章，希望能给大家带来一些新的研究思路。

文章简介

今天给大家带来的文章是发表在Bioinformatics上的名为“VirtualCytometry”的webserver，通过创建一个网络计算平台，利用单细胞基因表达差异来评估免疫细胞的分化。这篇文章所用到的数据均为公共scRNA-seq数据。

背景介绍

免疫系统，众所周知，包含多种细胞类型。这些细胞通过各种信号通路和特定条件下的转录调控被激活或者分化为其他亚型，因此分析这些信号分子和转录因子对免疫细胞功能产生的作用是一个重要的研究方向。常见的分析方法有流式细胞技术与质谱技术 (flow and mass cytometry)，但是此技术最大的局限是通量较低，单次实验只能检测少量的蛋白质。而单细胞RNA测序的发展实现了对单个细胞全基因组表达的定量分析，从而为研究免疫细胞分化提供的一种新的分析思路。并且由于近年来，公共数据库中人类和小鼠的免疫细胞scRNA-seq数据量大大增加，因此可以通过分析和整合来研究免疫细胞分化状态。

处于不同发育状态的细胞往往对应着不同的转录群体，因此可以通过在细胞群体间比较与分化相关的基因表达来鉴定细胞的分化状态。这种常规的方法已经有所实现。但是局限在于，处于分化状态的细胞往往表现出一种连续的状态（a continuous spectrum of cellular states），而不是具有高度差异的细胞亚群。因此我们需要别的方法来找出与分化相关的基因。由于标记基因（marker genes）能够用来来鉴别不同的细胞类型和状态，本文的研究假设是，根据标记基因表达水平来区分出的细胞亚群是跟细胞分化相关的，能够帮助我们找到与分化相关的差异表达基因。

作者开发了一款名为“VirtualCytometry”的计算平台，其主要的计算思路是根据用户选择的标记基因表达水平将所有细胞分为两个子集，找出两个子集的差异表达基因。然后通过correlation analysis，找出与标记基因高度相关的差异表达基因，则被认为是跟免疫细胞类型和功能有关的基因。

数据介绍

这篇文章用到的数据来自于JingleBells 和GEO。JingleBells 是一个专门收录跟免疫相关的scRNA-seq数据库。GEO相信大家都很熟悉了，他们在GEO中通过搜索“single cell”和“Expression profiling by high throughput sequencing”来获取数据。值得注意的是，这个平台是被用于检测已知细胞类型中的分化情况，因此筛选数据的时候要求数据集带有类型标签（cell type label），并且要求提供quantitative values/count values（一般是cell by gene matrix）。他们搜索的范围包括人类和小鼠的免疫细胞，不限制组织来源和疾病状态。

数据处理

数据处理流程如上图所示，包括了数据标准化，特征选择以及可视化。数据处理是基于R package: Seurat。他们的主要步骤包括：

1. 标准化： 对于提供raw count table的数据集，作者以count matrix 为起点，筛除低质量细胞并进行数据标准化。对于提供标准化数值的数据集（例如CPM, TPM, RPKM or FPKM）, 则进行log2(X+1)的变换。经过预处理的数据根据疾病状态，组织来源等进行分类。

2. 特征选择：作者选择了top 2000个高度可变基因。

3. 可视化：通过PCA来降维并用UMAP和t-SNE两种方法使数据可视化，PC数目的选择范围在10到50 之间。

4. 细胞聚类：通过Louvain聚类算法对细胞进行聚类。

5. 差异表达分析：通过差异表达分析找出每个细胞集的标记基因，满足条件的基因需要在至少10%的细胞中表达，并且log FC>=0.25并且 FDR<0.05.

识别与免疫细胞分化相关的基因

该平台通过模块“Discovery module”来识别与免疫细胞分化相关的基因。具体步骤包括

1: Selection of a cell population for analysis: 该平台共包含226个免疫相关的scRNA数据集，有人类和小鼠的。每一个数据集都经过了重新处理（如parr3描述），用户需要先选择一个数据集。

2. Interactive visual analysis of marker genes for subset segregation: 用户根据先验经验来选择一个标记基因，如果用户没有提供标记基因的话，平台会通过对数据集原先定义的细胞类型做差异表达分析，提供候选基因。确定标记基因之后，用户可以通过交互界面选择一个表达阈值，使得细胞能被“很好的”分为两个子集。理想情况下，表达值大于阈值的细胞在图中应该是聚在一起的，而不是分散的。如下图中黑色的点分布。

3. Identification of DEGs between the cell subsets: 细胞分为两个子集之后，通过差异表达分析得到差异表达基因。选择的阈值在前面的处理流程中提到了。他们还通过从数据库中收集转录因子和细胞表面受体信息来对差异表达基因进行注释。

4. Output plots to display the relationship between the marker gene and DEGs in the given cell population: 通过图形来显示标注基因和差异表达基因的相关性，如下图所示。

应用实例

实例演示中，VirtualCytometry被用于识别与人类肿瘤中衰竭CD8+ T 细胞（exhausted CD8+ T cells，一种t细胞功能异常的状态）相关的转录因子和信号分子。CD8+ T细胞是免疫系统攻击癌症细胞的主要利器。但是在肿瘤微环境里，CD8+ T细胞经常会进入机能缺陷或衰竭状态，从而不能有效地阻止癌症的进展。找到与CD8+ T细胞衰竭状态相关的基因能帮助我们了解细胞的分化程序，进而找到癌症相关的免疫疗法。

作者提出一个假设，与衰竭CD8+ T细胞相关的基因在衰竭早期和晚期的表达是有显著差异的，为了验证这一想法，他们选择了一个来源于肝癌的细胞毒性T细胞（cytotoxic T cells）的数据集。他们选择的其中一个标记基因为TCF7 (TCF1) 并且表达阈值为1。TCF7已知是CD8+ T 细胞（PD-1 低表达）的标记基因。从下图中可以看到，该基因表达值高于1 的细胞主要分布在t-SNE图的上半部分。通过对两个子集的差异分析，得到的其中一个显著差异基因SATB1, 被报道能够负向调节PDCD1的表达，从而防止T细胞的衰竭。而且从相关性图中，可以看出标记基因TCF7与SATB1高度相关。

小编总结

这篇文章所用到的数据都是公共数据，并且使用的是Seurat的标准处理流程，并没有用到任何花哨的分析方法。作者非常巧妙的运用了标记基因，开发了这个平台。其实这篇文章的思路和方法都非常常规，这也告诉我们，整合分析单细胞RNA数据其实并不一定意味着需要非常复杂的方法。随着单细胞RNA数据的增多，这种挖掘整合也可以应用到多种细胞类型中。有关挖掘公共单细胞RNA数据的项目咨询，欢迎大家与我们联系！

参考文献

Kyungsoo Kim, Sunmo Yang, Sang-Jun Ha, Insuk Lee, VirtualCytometry: a webserver for evaluating immune cell differentiation using single-cell RNA sequencing data, Bioinformatics, Volume 36, Issue 2, 15 January 2020, Pages 546–551, https://doi.org/10.1093/bioinformatics/btz610