FindNeighbors是KNN+SNN聚类 KNN计算最近邻,SNN计算共享最近邻-均是计算的过程,可以认为是将细胞进行连线的过程
FindClusters是Louvain算法分群 将计算好的细胞划分成不同的clusters,可以认为这是一个切割的过程,即将高维空间中连好线的细胞群切割成群
FindNeighbors耗内存,时间慢,需要先经过PCA降维,去除差异小的基因
目前每个细胞中的基因表达量都可以用我们熟知的单细胞转录组测序(scRNA-Seq)方法获得,可用于分析物种组织的细胞异质性,鉴定稀有细胞亚群或分化发育过程中的中间态细胞,反映细胞的发育阶段和分化状态等多种分析点。scRNA-Seq数据分析过程中最关键的步骤就是根据基因表达模式将属于相同细胞类型的细胞分组,使得组内相似、组间不同。然而单细胞转录组测序的表达量数据因具有以下特点使得细胞聚类变得困难:
(1)与传统转录组测序相比,可以获得精度更高、维度更广的每个细胞中的每个基因的表达量用于更多下游的个性化分析;
(2)由于对每个细胞中的每个基因表达量进行测序,然而现实生物体中不是所有细胞中的所有基因都会转录表达,因此scRNA-Seq的表达量数据是稀疏矩阵,平均稀疏度可以高达50%(比如一个细胞的很多基因表达量都是0,甚至0表达量的基因达到了50%);
(3)由于一个项目可以得到几万到几十万的细胞样本,基因数量通常也是几万个之间,这样高维度的数据集是难以分析细胞类型间基因表达模式的差异。因此,scRNA-Seq数据细胞聚类将面临以下计算挑战:稀疏性、高噪音、高维度。
|常见方法
细胞聚类的目标是根据细胞中各个基因表达模式的相似性(或距离)将一组细胞分组变成大类,使得这些大类成为有数学意义的亚群,这是scRNA-Seq数据挖掘的一个重要步骤和目标。目前有几种不同类型的scRNA-Seq细胞聚类算法,其中常用的方法是K-means,它迭代发现k-cluster的中心,并通过欧氏距离最小化将每个单元分配到最近的单元。然而,它是一种贪心算法,可能只会导致局部最优。它对异常值敏感,倾向于寻找大小相等的细胞簇,可能无法检测到罕见的细胞类型。集群的数量k必须是预先定义的,k的不同值将导致不同的集群。
scRNA-Seq聚类的常见方法是基于图(graph-based clustering)的方法。 顾名思义,这种方法有两个步骤组成,第一步是画图,第二步是识别图。画图通常由k-最近邻(KNN,k-nearest neighbor)和共享最近邻(SNN,shared nearest neighbor)两步来组成。 这张图可以看作细胞之间相似性计算的结果。识别图通常由社区检测算法——Louvain算法来实现,在图中挑选相似度最高的一群细胞作为一个细胞亚群。KNN和Louvain检测算法的结合已经产生了几个软件工具,包括PhenoGraph,Seurat和Scanpy。Louvain算法的主要优点是它的速度和可扩展性,集群的数量不需要预先定义。
下面就对KNN、SNN、Louvain算法的具体原理进行介绍,因为他们最常用也是seurat包中应用的方法,帮助您了解seurat的细胞聚类方法。
>1.KNN(k-nearest neighbor)
邻近算法,也叫k最近邻分类算法,它是数据挖掘分类技术中最简单的方法之一。所谓k最近邻,就是k个最近邻居的意思,说的是每个样本都可以用它最接近的k个邻居来代表。用通俗易懂的话说“人以类聚,物以群分”、“近朱者赤,近墨者黑”。那么最近邻居怎么理解呢?其实,KNN就是当预测一个未知的新细胞时,根据新细胞距离它最近的K个细胞是什么类别来判断这个新细胞属于什么类别。下面用一个示意图来轻松理解吧!
图1 当K=4时,判定新细胞属于红色圈亚群
如图1所示,图中黑色四边形就是要预测的新细胞,假设k=4 时,KNN算法就会找到与它距离最近的4个细胞(图中虚线圆圈),判断该圈中哪种类型的细胞多,就判定新细胞为细胞多的那个类型,例如图中所示,圆圈中红色类型3个,绿色1个,所以判定新细胞属于红色细胞类型。
图2 当k=7时,判定新细胞属于绿色圈亚群
但是,当k=7,预测的新细胞结果就和k=4的结果不一样了,变成了绿色亚群。从这里我们就可以看出K的取值太重要了。
还有一个决定结果变化的重要因素就是点的距离,下面对K值的选取和点的距离计算做一个详细说明。
|距离计算
要度量空间中点距离的话,有好几种度量方式,比如常见的曼哈顿距离计算,欧式距离计算等等。不过通常KNN算法中使用的是欧式距离(全称:欧几里得度量),先拿二维平面为例,二维空间两个点的欧式距离计算公式如下:
其实它和我们高中数学中在坐标系中计算两点之间的距离一样,比较简单。如果拓展到多维空间,那公式如下:
应用到单细胞转录组细胞聚类中,欧式距离的公式都代表什么呢?其实d(x,y)表示两个细胞之间的距离,根号下的(x1-y1)表示两个细胞中某一个基因的表达相减再平方,如果有10,000个基因的表达量,那这里相当于10,000维空间的欧式距离计算。
这样我们就明白了如何计算距离,KNN算法最简单粗暴的就是将预测点与所有点距离进行计算,然后保存并排序,选出前面k个值看看哪些类别比较多。---搜嘎!
|K值选取
通过上面的原理介绍我们知道k的取值比较重要,那么该如何确定k取多少比较好呢?这里可以通过交叉验证(将样本数据按照一定比例,拆分出训练用的数据和验证用的数据,比如6:4拆分出部分训练数据和验证数据),从选取一个较小的k值开始,不断增加k值,然后计算验证集合的方差,最终找到一个比较合适的k值。
通过交叉验证计算方差后你会得到下面这样的图:
图3 k值与错误率曲线图
如图3所示,当k值增大时,一般错误率会先降低,因为有周围更多的样本可以借鉴了,分类效果会变好。但注意,和K-means不一样,当k值更大的时候,错误率会更高。这也很好理解,比如说你一共就35个样本,当你K增大到30的时候,KNN基本上就没意义了。
所以选择k点的时候可以选择一个较大的临界k点,当它继续增大或减小的时候,错误率都会上升,比如图3中的k=18时,分类效果最佳。
|KNN算法优点
简单易用,相比其他算法,KNN算是比较简洁明了的算法,即使没有很好的数学基础也能理解它的原理;
模型训练时间快;
预测效果好;
对异常值不敏感。
|KNN算法缺点
对内存要求较高;预测阶段可能很慢;对不相关的功能和数据规模敏感。
也正是因为KNN具有的固有缺陷,我们在图论聚类的过程中,KNN并不承担分类器的职能,而是仅用于寻找每个细胞距离最近的k个细胞。此外,为了提高运算速度,降低背景噪音,主成分分析(PCA)会优先于KNN进行,如此,欧式距离的计算就是基于主成分特征值完成,降低了计算的维度。
2.SNN (shared nearest neighbor)
SNN是一种基于共享最近邻的聚类算法,它通过使用数据点间共享最近邻的个数作为相似度来处理密度不同的聚类问题,从而可以在含有噪音并且高维的数据集中发现各不相同的空间聚类。
那SNN是怎么计算的呢?它是在KNN的基础上,通过计算数据对象之间共享最近邻相似度得到SNN图。计算过程如图4 所示。
图4 SNN算法计算过程
首先在细胞对之间使用欧式距离来计算相似性矩阵(例如,点代表细胞,并且使用细胞基因表达水平向量来计算两个细胞之间的距离)。接下来,对于每个数据点Xi,我们使用相似性矩阵列出k个最近邻居(KNN),其中Xi本身作为列表中第一个条目。然后,为了构造SNN图,对于一个细胞对Xi和Xj,仅当 Xi和Xj具有至少一个共享节点时,我们才分配边e(Xi,Xj)(图5)。
图5 e(Xi,Xj)概念
边e(Xi,Xj)的权重定义为k与共有的KNN的最大平均排名之间的差,具体公式如下:
其中k是最近邻居列表的大小,在Seurat包中k默认为20,就是说利用KNN算法找到某一个细胞的20个邻居,并按照欧式距离从小到大排序。
其中rank(v,Xi)表示共享邻居节点v在Xi所有邻居中的排位。注意,相近邻居v 的排位越高,而rank(v,Xi)越低,例如rank(Xi,Xi)=1,是因为Xi最先排列在Xi的最近邻居列表中第一个位置;
该SNN图捕获了两个节点之间在其连通邻域中的相似性,换句话说,与之前相似性不同,在我们的计算过程中,两个节点之间的相似性需要通过它们与其他节点(常见的最近邻居)的接近度来确认。
SNN背后的基本原理是节点的排名通常在高维空间中仍然有意义,尽管之前的距离计算可能不是这样。并且期待同一群组中两个节点的共享邻居的排名很高,从而使得边上的连接权重较高。相反,我们期待来自不同群组的两个节点的共享邻居的排名较低,使得边上连接权重较低。同时共享最近邻图通常是稀疏的,因此允许扩展到大数据集,适合单细胞数据集分析。
3.社区检验算法—Louvain
前两步KNN和SNN主要计算了细胞之间的联通关系。Louvain算法在单细胞聚类中的作用主要是对前两步的计算结果进行亚群划分,也就是起分类器的作用,该算法引入网络去噪,能够有效提升细胞聚类性能。
Louvain算法是基于模块度(Modularity)的社区发现算法,该算法在效率和效果上都表现比较好,并且能够发现层次性的社区结构。
Louvain算法的总体思路是:
(1)刚开始时,每个细胞或者由20个细胞组成的小簇;
(2)通过模块度函数值优化,将每个节点归到聚类效果“最好”的簇中。该步骤直到所有的节点所属的簇不再变化时,程序才停止;(
3)将网络图进行分解和简化,把一个簇内的所有节点抽象的看作是一个节点,看抽象后的网络图是否还有优化的可能性。如果有,对抽象之后的点继续进行第(2)步。
这里我们可以看到,理解Louvain算法的核心就是要理解模块度的计算过程。
模块度的定义就是描述社区内紧密程度的值,或者说是评估一个社区网络划分好坏的度量方法,记作Q,取值范围为 [ -1/2 , 1 ),公式如下:
看到这个公式是不是懵了呢?经过对这个公式的简化理解,我们对中括号中间关键的分类部分做一个说明(公式的其他部分对分群的影响较小,可以忽略)。直接看下面这张实例图。
图6 随机权重网络图
如图6所示,数字是两个节点的权重值,该图是初始的graph,A-F一共6个节点。初始阶段,节点之间的合并以节点A为例:
以A-B为例,
5 为A和B之间的权重;
30 是整个graph的权重之和;
10 是A和所有和它连接的边的权重之和;
7 是B和所有和它连接的边的权重之和。
这样可以对公式的中括号部分有了一点理解。Aij是实际的节点A指向社区B的权重,后面那一项表示的是平均意义之下节点A指向社区B的期望权重,如果真实权重大于期望权重,则代表了节点A和社区B的连接关系是存在显著的意义的;
同理,可以得到所有节点的最大模块度增益的计算结果:
根据公式说明,小于0则不进行社区划分,大于0则取最大的模块度增益对应节点进行合并,因此AB合并、CD合并、EF合并,合并之后我们得到(图7所示):
图7 合并后的结果
图7为初始阶段的合并结果,然后算法继续运行:
这里需要注意的地方!graph的权重发生改变了,m=10而不是30,因为此时划分出来的社区已经拿走了一部分边的权重了,所以实际上整个graph的权重应该去除社区内的边的权重重新计算!
根据上面的计算可以知道模块度增益都是小于0的,迭代停止,此时即为最终的社区划分的结果。
利用Louvain算法对6个节点分群,最终得到3个亚群的结果。