python中基于基因表达矩阵的层次聚类

Question

我如何在Python中进行分层聚类(在本例中是针对基因表达数据)，以显示基因表达值矩阵和树状图？我的意思是像下面这样的例子：

http://www.mathworks.cn/access/helpdesk/help/toolbox/bioinfo/ug/a1060813239b1.html

在项目符号6(图1)之后显示，其中树状图绘制在基因表达矩阵的左侧，其中行已重新排序以反映聚类。

如何在Python中使用numpy/scipy或其他工具执行此操作？另外，用欧几里德距离作为度量，用大约11,000个基因的矩阵来做这件事，在计算上可行吗？

编辑:很多人建议使用聚类包，但我仍然不确定如何绘制上面在Python中链接的图像。例如，我如何使用Matplotlib将树状图覆盖在热图矩阵旁边？

谢谢。

Answer 1

包括scipy.cluster在内的许多聚类方法都是通过对所有成对距离进行排序开始的，在您的情况下约为6000万，不会太大。

完成以下操作需要多长时间？

import scipy.cluster.hierarchy as hier
import pylab as pl

def fcluster( pts, ncluster, method="average", criterion="maxclust" ):
    """ -> (pts, Y pdist, Z linkage, T fcluster, clusterlists)   
        ncluster = n1 + n2 + ... (including n1 singletons)
        av cluster size = len(pts) / ncluster
    """
    pts = np.asarray(pts)
    Y = scipy.spatial.distance.pdist( pts )  # ~ N^2 / 2
    Z = hier.linkage( Y, method )  # N-1                         
    T = hier.fcluster( Z, ncluster, criterion=criterion )
        # clusters = clusterlists(T)
    return (pts, Y, Z, T)

hier.dendrogram( Z )

3月份，here在《So》中问到了如何对矩阵进行置换和绘制，并给出了部分答案。

Answer 2

您可以使用scipy的cluster.hierarchy模块来实现这一点。这些命令实际上甚至非常相似。但是，您必须使用correlation而不是corr作为pdist的参数，并且函数scipy的集群模块的名称是fcluster而不是cluster。此外，对于树状图，函数在scipy中为dendrogram，而不是在Matlab中为clustergram。

您绝对可以使用欧几里得度量(认为它是pdist的默认度量)。我认为对11,000个基因这样做应该是可行的，因为这将是要计算的11000* (11000 -1)/2 = 60494500 (11000选择2)距离。这是一个很大的数字，但我认为肯定是可行的。

Answer 3

有几个人在使用scipy和matplotlib创建用于分层聚类和热图可视化的原型模块方面取得了一些进展：

How to get flat clustering corresponding to color clusters in the dendrogram created by scipy

我一直在修改这段代码，以制作一个完整的层次聚类模块，我可以将其集成到我的一个转录组分析包中。我对最终的产品非常满意，它将使用各种聚类度量和方法以及着色梯度来生成热图。代码和示例输出如下所示：

http://altanalyze.blogspot.com/2012/06/hierarchical-clustering-heatmaps-in.html