生物_医药_科研 / 待分类 / 单细胞聚类方法大比拼(III)

分享

   

单细胞聚类方法大比拼(III)

2020-06-07  生物_医药...
嘿,大家好我又来了!回忆一下,在前两期的文章中我已经给大家介绍了常见的聚类如层次聚类、k-means等,以及为单细胞专门开发的SC3聚类。今天带来本系列的最后一个聚类方法——谱聚类(Spectral Clustering)的介绍。听这个名字大家可能会觉得很陌生,但是如果你有过单细胞的分析经历,可能会在无意中使用过它。没错,它就是大名鼎鼎的RSeurat自带的默认聚类方法。

单细胞聚类方法大比拼(III

简单来说,Seurat中的谱聚类基于共享最近邻图(SNN)和模块化优化的聚类算法识别细胞簇首先计算k最近邻(k-nearest neighbors)并构造SNN 然后优化模块化功能以确定具体类群听起来好像有点难懂,没办法,高大上的方法背后一般都是深奥的理论作为支撑,并且这个算法最早居然在2013年发表在物理学的期刊The European Physical Journal B上。我会尽量简化其中所涉及到的数学相关知识,更多地传达方法背后生物学意义。

第一个难点是谱聚类是一种基于图论方法,如果你傻傻地分不清图论的“图”和图像的“图”,建议先去学习一下基本的离散数学。谱聚类是从图论中演化出来的算法,它的主要思想是把所有的细胞看为高维空间中的点,这些点之间可以用边连接起来。距离较远的两个点之间的边权重值较低;距离较近的两个点之间的边权重值较高。然后通过对所有数据点组成的图进行“切图”,让切图后不同的子图间边权重和尽可能的低,而子图内的边权重和尽可能的高,从而达到聚类的目的。

谱聚类中的图采用的是无向有权图,也就是说,每个细胞(图中的点)之间的连接是没有方向性的。但是每条边上是有权重信息的。类比一下地图,两地之间的路是没有方向(双向的),但是地点与地点之间路的远近距离是不同的。对于一个无向有权图G,一般用点的集合V和边的集合E来描述即为G(V,E)。其中V即为我们数据集里面所有的点(v1,v2,...vn)的集合。对于V中的任意两个点,定义权重wij为点vi和点vj之间的权重由于是无向图,所以wij=wji如果两点之间没有边连接,那么wij=0OK,如果你看到这都能懂,那么接下来也没啥问题了。

第二个难点如何把我们的单细胞的数据构建成为无向有权图,需要明确的是,我们手中最初只有细胞的表达量矩阵,这个每个细胞中各个基因的表达情况,属于“节点信息”。那么哪些细胞之间需要有边相连,边的权重又该如何计算呢?这里引入了邻接矩阵(W)的概念,我们首先定义一种计算细胞之间距离度量的方法,来计算每条边之间的权重wij,又称相似性sij,然后获得整个的邻接矩阵W。构建邻接矩阵W常用方法有三类ϵ-邻近法,K邻近法和全连接法Seurat中选用的是k-近邻法,那我们也就把k-近邻法作为例子讲解。首先利用KNN算法遍历所有的样本点,取每个细胞最近的k细胞作为近邻,然后设置只有与此细胞距离最近的k个点之间的wij>0,其他全设置为0。但是这种方法会造成重构之后的邻接矩阵W非对称,我们后面的算法需要对称邻接矩阵。为了解决这种问题,一般采取下面两种方法之一(1)只要第一个点在另一个点的k近邻中即可保留sij2)必须两个点互为k近邻才可保留sij。由此我们就能用计算出的权重矩阵构建出如下所示的图。

 

然后是度矩阵(D)的概念:对于图中的任意一个点vi,它的度di定义为和它相连的所有边的权重之和,即

 从每个点度的定义,我们可以得到一个n x n的度矩阵D,它是一个对角矩阵,只有主对角线有值,对应第i行的第i个点的度的值

好了,最后一个重要的矩阵:拉普拉斯矩阵(L)登场。拉普拉斯矩阵定义很简单:L = D -W , W上面所说的邻接矩阵D是度矩阵。但是拉普拉斯矩阵有着很重要的性质。

终于到了“切图的环节,实际就是将图切成相互没有连接的k个子图。目标就是让子图内的点权重和高,子图间的点权重和低。由于这里的求解过程十分复杂全是公式,具体的过程就不讲了,实际就是一个最优化问题。可以理解为通过DL这两个矩阵所给的信息,在图的部分边进行切分,让完整的图变成许多子图,每一个子图就构成了一个类群,如下图所示

 
当然,如果你只想在Seurat中简单的用一下,直接用FindClusters()这个函数就好了。就是简单的一句话’pbmc <- FindClusters(pbmc, resolution = int)’。最常调整参数有resolutionpc.usepc.use决定哪些维度用来进行聚类,如果之前做过PCA,会选一般前10-20维度。resolution参数决定类群的“粒度”,官网建议将此参数设置在0.4-1.2之间通常可为包含大约3K个cell的数据集返回良好的聚类结果。对于较大的数据集,最佳resolution通常会增加。值得一提的是,FindClusters()函数是不需要输入希望划分类群数量的,由算法自动划分。调整resolution会得到不同的类群的结果,越高细胞类群数量越多,想获得更多类群,可以尝试高于1.0的值

最后加上一句umap或者t-SNE的绘图命令就能可视化聚类结果了。下图是我用Seurat自带的pbmc数据的聚类结果,可以发现效果还是非常不错的。 


好啦,以上就是谱聚类的所有内容了。如果你不太懂算法原理,可以不用过于纠结,能在R语言中灵活运用才是解决生物信息问题的关键。

    本站是提供个人知识管理的网络存储空间,所有内容均由用户发布,不代表本站观点。请注意甄别内容中的联系方式、诱导购买等信息,谨防诈骗。如发现有害或侵权内容,请点击一键举报。

    0条评论

    发表

    请遵守用户 评论公约

    类似文章 更多
    喜欢该文的人也喜欢 更多

    ×
    ×

    ¥.00

    微信或支付宝扫码支付:

    开通即同意《个图VIP服务协议》

    全部>>