单细胞中的流形（一）： 理解 tSNE中的perplexity

单细胞聚类的一般步骤有：数据质控，基因筛选，PCA降维，以及基于PCA预降维空间的聚类和进一步的降维可视化。之前的文章介绍了单细胞数据的PCA降维和聚类算法。本次文章介绍基于PCA空间的降维可视化算法tSNE，将单细胞数据降至2维或者3维可视化。

 tSNE(t-distributed stochastic neighbor embedding)是一种非线性降维的算法，可以归入流形学习的范畴。流形学习的经典方法有：Isomap，LLE，LE和diffusion maps等[1]。流行学习中的“流形”，指的是一种局部可以近似为欧式空间的拓扑空间。通俗的讲，流形可以理解为高维空间中的低维映射结构，如下图中的瑞士卷，用三维表示是比较冗余的。非线性降维的方法可以用更低维度的坐标表示高维空间中的流形结构[2]。

tSNE可视化的利与弊：

优点：

1. 相对于PCA有更高的区分度，不同细胞类型分得开。

2. 保留局部结构，相类似的细胞聚集在一起，同时不会过于重叠扎堆。

缺点：

1. 相对于UMAP，不能很好的展示全局结构。

2. 经典的tSNE算法不适用于较大的数据，耗时较长O(N²).。

tSNE算法原理

        在单细胞研究中，细胞群体的二维展示图可以说是她的衣服，可视化效果的好坏直接影响整个研究的观感。tSNE的可视化效果受参数设置的影响，如困惑度（perplexity）等。所以有必要了解一下tSNE的原理，选择合适的参数，提高单细胞可视化的效果。

tSNE降维过程分成两大部分[3]：

1. 高纬度数据的欧式距离转成概率分布（高斯分布）；

2. 用t分布初始化tSNE的二维点，KL散度表示高斯分布和t分布的差异，进行梯度优化，使t分布逼近高斯分布，最终确定tSNE二维点。

欧氏距离转概率分布

        我们以四个细胞x50个维度的数据为例。首先以细胞a为对象，用正态分布的公式转化细胞a与细胞b/c/d的距离，然后归一化，把距离变成概率。

        相同的，以b/c/d为对象，它们也有自己的距离概率分布。在做距离转换时，每个细胞的概率分布的标准差σ是未知的，需要用我们设定的tSNE参数，perplexity，来估算σ。

参数perplexity的作用和原理

        困惑度（perplexity）可以表示细胞的邻近个数，在tSNE图上的直观反映是细胞点的分布是否紧凑。perplexity设置越大，细胞分布越紧凑。

        当然perplexity有合适的取值范围，之前提到perplexity的值设置越大，细胞分布越紧凑。但perplexity有一个上限，R包“Rtsne”中规定perlexity* 3 < 细胞数 – 1，如果按照默认参数perlexity=30，来降维我们的4个细胞数据，会报错：“Perplexity too large”。有文章总结perplexity设置为：细胞数 / 50 比较合适[6]。根据个人经验，对于细胞数比较少的数据（如胚胎细胞SmartSeq2数据），用接近上限的perlexity值可视化效果比较好，如31个细胞，用(31-1)/3 – 1 = 9 设置perplexity。

        细胞a距离分布的perplexity值由2^香农熵（Shannon Entropy）计算（下图），perplexity和香农熵正相关。香农熵用细胞a距离分布相关的公式表示，此时细胞a的距离是带有未知数σ的，根据我们给定的perplexity值（默认30），用二分搜索的方法，迭代一定次数，就可以确定σ。

香农熵的原理

        香农熵H可以表示细胞a距离分布参差不齐的程度。香农熵越大，perplexity越大，距离分布标准差σ越小，细胞a的邻近细胞个数越多，越紧凑。

        香农熵在生物信息领域也有其他的应用。如衡量免疫组库TCR/BCR种类多样性的指标。华大基因张伟开发的免疫组库分析流程中，就是用香农熵来表示CDR3序列种类的多样性[4,5]。CDR3多样性越高，香农熵越大。其中，CDR3的数据表示用reads count标准化后的概率值，这跟tSNE距离转换后的数据形式相似。

T分布的梯度优化

        确定高维数据的概率分布之后，用t分布初始化tSNE的二维点，KL散度表示高斯分布和t分布的差异，进行梯度优化，缩小cost函数的值，使t分布逼近高斯分布，最终确定tSNE二维点。由于这一步在单细胞数据tSNE展示的参数设置作用不大，算法的实现过程就不细讲了。一般默认梯度优化的迭代次数是1000。

T分布的作用

        t分布使得相近的细胞更加紧凑，较远的细胞更加疏远，同时保持了细胞之间的差异，不会过于重叠。不过用t分布做梯度优化的过程中有一个特点，随着cost函数的变小，会倾向于用t分布中较大的概率值（距离较小）表示正态分布中较小的概率值（距离较大）。换句话说，tSNE倾向于保留数据中的局部结构。距离较远的两个细胞在tSNE图上的表示可能失真。

Barnes-Hut-SNE加快tSNE的运行速度

        经典的tSNE算法不适用于较大的数据，耗时较长O(N²)。后来tSNE的作者在2013年推出新方法Barnes-Hut-SNE，用四叉树（quadtree）的方法减少运行时间O(N*logN) [7]。目前R包“Rstne”正是基于这种方法。

总结：

1．用tSNE展示单细胞数据有比较高的区分度；

2．tSNE的参数设置：perplexity < (细胞数-1)/3，建议perplexity = 细胞数 / 50；

3．tSNE倾向于保留数据的局部结构。

      鉴于tSNE的缺点，近年来umap逐渐在单细胞数据的可视化方法中占据位置。umap相对tSNE来说，保留了数据的全局结构，而且大数据的运行时间较短。

下期预告：

流形学习专题2：umap

校对人：QQQ

作者才疏学浅，有错误之处敬请勘正。

参考文献：

[1].WIREs Comput Stat 2012, 4:439–446. doi: 10.1002/wics.1222

[2]. Tenenbaum,Joshua B., Vin De Silva, and John C. Langford. "A global geometricframework for nonlinear dimensionality reduction." science 290.5500(2000): 2319-2323.

[3]. Maaten,Laurens van der, and Geoffrey Hinton. "Visualizing data using t-SNE."Journal of machine learning research 9.Nov (2008): 2579-2605.

[4].Carlson, C. S., R. O. Emerson, A. M. Sherwood, C. Desmarais, M. W.Chung et al.,2013 Using synthetic templates to design an unbiased multiplex PCR assay. Nat.Commun. 4: 2680.

[5]. Zhang,Wei, et al. "IMonitor: a robust pipeline for TCR and BCR repertoireanalysis." Genetics 201.2 (2015): 459-472.

[6]. Kobak,Dmitry, and Philipp Berens. "The art of using t-SNE for single-celltranscriptomics." Nature communications 10.1 (2019): 1-14.

[7]. VanDer Maaten, Laurens. "Barnes-hut-sne." arXiv preprint arXiv:1301.3342(2013).