实战--利用HierarchicalClustering 进行基因表达聚类分析

利用建立分级树对酵母基因表达数据进行聚类分析

一、原理

根据基因表达数据，得出距离矩阵

　　↓

• 最初，每个点都是一个集合
• 每次选取距离最小的两个集合，将他们合并，然后更新这个新集合与其它点的距离
  • 新集合与别的集合距离的计算方法
  • ①两个集合之间的最短距离
  • ②两个集合所有点之间求距离求平均   →
• 把这个新集合加入距离矩阵中，原来的两个小集合就被替换掉
• 如此循环，直到剩下一个集合，那就建立了一棵树
• 在树的某一处横断，就可以得到6类

230个酵母基因表达数据

http://bioinformaticsalgorithms.com/data/realdatasets/Clustering/230genes_log_expression.txt

二、python下利用Sklearn包实现

Sklearn包的安装

参照 https://scikit-learn.org/stable/install.html

代码（python 3.7环境）

from sklearn.cluster import AgglomerativeClustering
import numpy as np
from os.path import dirname
import numpy as np
import math
import random
import matplotlib.pyplot as plt

def InputData(dataset):
    dataset = [line.split() for line in dataset]
    name = [item[1] for item in dataset[1:]]
    points = []
    for line in dataset[1:]:
        if(len(line)==10):
            points.append(list(map(float,line[3:])))
        elif(len(line)==9):
            points.append(list(map(float,line[2:])))
    return points

if __name__ == '__main__':

    f = open("output","w")

    INF = 999999
    dataset = open(dirname(__file__)+'230genes_log_expression.txt').read().strip().split('\n')

    X = np.array(InputData(dataset))

    clustering = AgglomerativeClustering(n_clusters=6).fit(X)

    #print(clustering.labels_)

    lb = clustering.labels_

    clusters = [[] for i in range(6)]

    for i in range(len(X)):
        clusters[lb[i]].append(i)

    fig=plt.figure()

    x = [i for i in range(1,8)]
    for c in range(6):
        plt.subplot(231+c)
        for i in clusters[c]:
            plt.plot(x,X[i],linewidth=0.5,linestyle='-',marker='')
    plt.show()

结果

！注意

初次使用sklearn的时候python3.7解释器会报错

/Library/Frameworks/Python.framework/Versions/3.7/lib/python3.7/site-packages/sklearn/externals/joblib/externals/cloudpickle/cloudpickle.py:47: DeprecationWarning: the imp module is deprecated in favour of importlib; see the module's documentation for alternative uses import imp

是因为从python3.4开始就不再支持 import imp ， 需要按照错误信息找到相关的文件，将import imp 改为 import importlib，问题即解决

三、手动实现建立 Hierarchical树

 '''
 coder Lokwongho 2018.11

 '''

 from os.path import dirname
 import numpy as np
 import math
 import random
 import matplotlib.pyplot as plt

 ######## Get the Distance Matrix ########
 def PearsonCorrelationDistance(Expres):
     n = len(Expres)
     l = len(Expres[0])
     Distance = np.zeros(shape=[n,n],dtype=float)

     u = [ sum(i)/l for i in Expres ]
     Sigma = []
     for i in range(n):
         sig=0;
         for x in range(l):
             sig += (Expres[i][x]-u[i])**2;
         Sigma.append(sig)

     for i in range(n):
         for j in range(i,n):
             sigU=0
             for x in range(l):
                 sigU += (Expres[i][x]-u[i])*(Expres[j][x]-u[j])
             Distance[i][j]=1-(sigU)/math.sqrt(Sigma[i]*Sigma[j])
             Distance[j][i]=Distance[i][j]                

     return Distance

 ######## Build the Hierarchical Tree ########
 def locateMin(size): # calculate the diatance between two cluster: minimum distance method
     minVal = 999999
     minId = [-1,-1]
     for x in range(size):
         for y in range(size):
             if x!=y and matrix[x][y]<minVal and delete[x]==False and delete[y]==False:
                 minVal = matrix[x][y]
                 minId = [x,y]
     #print(minId,size)
     return [minId,minVal]

 def buildTree():
     Tree = [[]for i in range(maxn)]
     global matrix
     global delete
     ptr = n
     delete = [False for i in range(maxn)]
     num = [0 for i in range(maxn)]
     for i in range(n):
         num[i] = 1
     weight = [0 for i in range(maxn)]
     while(ptr<maxn):

         [[minX,minY],minVal] = locateMin(ptr)
         Tree[ptr].append([minX,minVal/2-weight[minX]])
         Tree[ptr].append([minY,minVal/2-weight[minY]])
         Tree[minX].append([ptr,minVal/2-weight[minX]])
         Tree[minY].append([ptr,minVal/2-weight[minY]])
         weight[ptr]=minVal/2
         delete[minX]=True
         delete[minY]=True
         #print(minX,minY,minVal)
         #print(delete)
         for i in range(ptr+1):
             if delete[i]==False:
                 tmp=(matrix[minX][i]*num[minX]+matrix[minY][i]*num[minY])/(num[minX]+num[minY])
                 matrix[ptr][i] = tmp
                 matrix[i][ptr] = tmp
         num[ptr] = num[minX]+num[minY]
         ptr += 1
     return Tree

 def HierarchicalClustering(DMatrix):
     global maxn
     global n
     global Tree
     global matrix

     n = len(DMatrix)
     maxn = 2*n-1
     matrix = np.zeros(shape=(maxn,maxn))
     for i in range(n):
         matrix[i][:n] = DMatrix[i]
     Tree = buildTree()

 ######## Input and output Data ########
 def InputData(dataset):
     dataset = [line.split() for line in dataset]
     name = [item[1] for item in dataset[1:]]
     #print(m)
     # print(m,k)
     points = []
     for line in dataset[1:]:
         if(len(line)==10):
             points.append(list(map(float,line[3:])))
         elif(len(line)==9):
             points.append(list(map(float,line[2:])))
     return [points,name]

 '''
 Reference to 'whoami_T' in csdn for print the tree
 https://blog.csdn.net/weixin_39722498/article/details/81534247
 '''
 def tree(lst):
     # 树状图输出列表
     l = len(lst)
     if l == 0:
         print('-' * 3)
     else:
         for i, j in enumerate(lst):
             if i != 0:
                 f.write(tabs[0])
                 print(tabs[0], end='')
             if l == 1:
                 s = '=' * 3
             elif i == 0:
                 s = '┬' + '-' * 2
             elif i + 1 == l:
                 s = '└' + '─' * 2
             else:
                 s = '├' + '─' * 2
             f.write(s)
             print(s, end='')
             if isinstance(j, list) or isinstance(j, tuple):
                 if i + 1 == l:
                     tabs[0] += blank[0] * 3
                 else:
                     tabs[0] += '│' + blank[0] * 2
                 tree(j)
             else:
                 print(name[j])
                 f.write(name[j] + "\n")
     tabs[0] = tabs[0][:-3]

 def traversalTree(v):

     global visited
     if v < n:
         return [v]
     visited[v] = 1
     items = []
     for i in Tree[v]:
         if visited[i[0]] == 0:
             items += traversalTree(i[0])

     return [items]    

 def outPut():
     global visited
     global tabs
     global blank

     blank = [
             chr(183)]  ##此处为空格格式;Windows控制台下可改为chr(12288) ;linux系统中可改为chr(32)【chr(32)==' ' ;chr(183)=='·' ;chr(12288)=='　'】
     tabs = ['']
     visited = [0 for i in range(maxn)]
     TreeList = traversalTree(maxn-1)
     #print(Tree)
     tree(TreeList)
 ###########################

 if __name__ == '__main__':

     f = open("output","w")

     INF = 999999
     dataset = open(dirname(__file__)+'230genes_log_expression.txt').read().strip().split('\n')

     [points,name] = InputData(dataset)

     DMatrix = PearsonCorrelationDistance(points)

     HierarchicalClustering(DMatrix)

     outPut()

===┬--┬--┬--┬--YJL109C
···│··│··│··└──YNL174W
···│··│··└──┬--YLR129w
···│··│·····└──┬--YGR264C
···│··│········└──┬--YLR449W
···│··│···········└──YNL002C
···│··└──┬--┬--YOL039W
···│·····│··└──┬--YNL067W
···│·····│·····└──┬--YDR382W
···│·····│········└──YGL031C
···│·····└──┬--YBR238C
···│········└──┬--YNL065W
···│···········└──┬--YNL303W
···│··············└──┬--┬--┬--YLR249W
···│·················│··│··└──YPL131W
···│·················│··└──┬--┬--YBR032W
···│·················│·····│··└──┬--YNL069C
···│·················│·····│·····└──┬--YBL027W
···│·················│·····│········└──YNL119W
···│·················│·····└──┬--┬--┬--┬--┬--YML063W
···│·················│········│··│··│··│··└──┬--YHL015W
···│·················│········│··│··│··│·····└──YDL136w
···│·················│········│··│··│··└──┬--YLR062C
···│·················│········│··│··│·····└──┬--YLL045c
···│·················│········│··│··│········└──YBR181C
···│·················│········│··│··└──┬--YPL220W
···│·················│········│··│·····└──┬--┬--YDR418W
···│·················│········│··│········│··└──YIL018W
···│·················│········│··│········└──┬--YER131w
···│·················│········│··│···········└──YLR198C
···│·················│········│··└──┬--┬--YGL102C
···│·················│········│·····│··└──YJL190C
···│·················│········│·····└──┬--┬--YBR189W
···│·················│········│········│··└──┬--YJL136C
···│·················│········│········│·····└──YBR191W
···│·················│········│········└──┬--YMR121C
···│·················│········│···········└──┬--YLR076C
···│·················│········│··············└──YLR325C
···│·················│········└──┬--┬--┬--YJR145C
···│·················│···········│··│··└──YLR344W
···│·················│···········│··└──┬--YHL033C
···│·················│···········│·····└──YOL040C
···│·················│···········└──┬--┬--┬--YDR064W
···│·················│··············│··│··└──┬--YHR089C
···│·················│··············│··│·····└──YDL208W
···│·················│··············│··└──┬--YIL069C
···│·················│··············│·····└──┬--YNL096C
···│·················│··············│········└──┬--YOR234C
···│·················│··············│···········└──┬--YDR417C
···│·················│··············│··············└──YGR148C
···│·················│··············└──┬--┬--YKL009W
···│·················│·················│··└──┬--YNL301C
···│·················│·················│·····└──YOL120C
···│·················│·················└──┬--┬--YLR340W
···│·················│····················│··└──┬--YJR123W
···│·················│····················│·····└──YLR048w
···│·················│····················└──┬--┬--YGR214W
···│·················│·······················│··└──YOR309C
···│·················│·······················└──┬--┬--YOR310C
···│·················│··························│··└──YOR312C
···│·················│··························└──┬--YEL054c
···│·················│·····························└──YKR059W
···│·················└──┬--┬--┬--┬--┬--YGL076C
···│····················│··│··│··│··└──YNL175C
···│····················│··│··│··└──┬--YPR137W
···│····················│··│··│·····└──YDL083C
···│····················│··│··└──┬--┬--YJL148W
···│····················│··│·····│··└──┬--YGL078C
···│····················│··│·····│·····└──YAL012W
···│····················│··│·····└──┬--YMR217W
···│····················│··│········└──┬--YGR103W
···│····················│··│···········└──YLR196W
···│····················│··└──┬--YLR186W
···│····················│·····└──┬--YDR025W
···│····················│········└──YBR247C
···│····················└──┬--┬--┬--┬--┬--YHR128W
···│·······················│··│··│··│··└──YKL081W
···│·······················│··│··│··└──┬--┬--YLR180W
···│·······················│··│··│·····│··└──YNL141W
···│·······················│··│··│·····└──┬--YDR060w
···│·······················│··│··│········└──YLR413W
···│·······················│··│··└──┬--┬--┬--YMR290C
···│·······················│··│·····│··│··└──YPL012W
···│·······················│··│·····│··└──┬--YDR144C
···│·······················│··│·····│·····└──YIL053W
···│·······················│··│·····└──┬--┬--YJL177W
···│·······················│··│········│··└──YBR249C
···│·······················│··│········└──┬--YLL044W
···│·······················│··│···········└──YLR448W
···│·······················│··└──┬--┬--YBR048W
···│·······················│·····│··└──YMR131C
···│·······················│·····└──┬--┬--YLR339C
···│·······················│········│··└──YNR053C
···│·······················│········└──┬--YPL226W
···│·······················│···········└──┬--YDR398W
···│·······················│··············└──YAL003W
···│·······················└──┬--YLR355C
···│··························└──┬--YGR160W
···│·····························└──YMR093W
···└──┬--YFR053C
······└──┬--┬--┬--YDL085w
·········│··│··└──┬--YBR117C
·········│··│·····└──┬--┬--YBR116C
·········│··│········│··└──┬--YKL187C
·········│··│········│·····└──YLR267W
·········│··│········└──┬--YDL199c
·········│··│···········└──┬--YBL043W
·········│··│··············└──YBL049W
·········│··└──┬--┬--YBL108W
·········│·····│··└──YCL025C
·········│·····└──┬--┬--YBR051W
·········│········│··└──┬--┬--┬--YNL117W
·········│········│·····│··│··└──YCR010C
·········│········│·····│··└──┬--YER024w
·········│········│·····│·····└──YGR067C
·········│········│·····└──┬--YDL215C
·········│········│········└──┬--YLR377C
·········│········│···········└──┬--YER065c
·········│········│··············└──YKR097W
·········│········└──┬--┬--YLR142w
·········│···········│··└──┬--YIL125W
·········│···········│·····└──┬--YNL195C
·········│···········│········└──┬--YJL089W
·········│···········│···········└──YAL054C
·········│···········└──┬--┬--YJR095W
·········│··············│··└──┬--YOL084W
·········│··············│·····└──┬--YLR174W
·········│··············│········└──┬--YHR096C
·········│··············│···········└──YJL045W
·········│··············└──┬--YEL012w
·········│·················└──┬--YBL045C
·········│····················└──┬--YGL259W
·········│·······················└──YLR312C
·········└──┬--┬--YCR021c
············│··└──┬--YDR343C
············│·····└──┬--YER053c
············│········└──YMR250W
············└──┬--┬--┬--┬--┬--YIL136W
···············│··│··│··│··└──┬--YMR090W
···············│··│··│··│·····└──YNL305C
···············│··│··│··└──┬--┬--YHR087W
···············│··│··│·····│··└──YKL142W
···············│··│··│·····└──┬--YDR031w
···············│··│··│········└──YJR096W
···············│··│··└──┬--YDR516C
···············│··│·····└──YIL111W
···············│··└──┬--┬--YDR342C
···············│·····│··└──YBR183W
···············│·····└──┬--┬--┬--┬--YBR072W
···············│········│··│··│··└──┬--YKL085W
···············│········│··│··│·····└──YLR327C
···············│········│··│··└──┬--┬--YGR008C
···············│········│··│·····│··└──┬--YNL200C
···············│········│··│·····│·····└──┬--YNL173C
···············│········│··│·····│········└──YOL053C
···············│········│··│·····└──┬--YGR248W
···············│········│··│········└──YNL015W
···············│········│··└──┬--YLR258W
···············│········│·····└──┬--YKL103C
···············│········│········└──┬--YGL037C
···············│········│···········└──YLR178C
···············│········└──┬--┬--┬--YHL021C
···············│···········│··│··└──YML128C
···············│···········│··└──┬--YFR015C
···············│···········│·····└──YMR105C
···············│···········└──┬--┬--YFL014W
···············│··············│··└──YOR374W
···············│··············└──┬--YGR244C
···············│·················└──┬--YLL041c
···············│····················└──┬--YER150w
···············│·······················└──YNL160W
···············└──┬--┬--┬--┬--┬--YLL026w
··················│··│··│··│··└──YDR258C
··················│··│··│··└──┬--YDR258C
··················│··│··│·····└──YLR149C
··················│··│··└──┬--┬--YBR139W
··················│··│·····│··└──YLR304C
··················│··│·····└──┬--┬--YEL024w
··················│··│········│··└──┬--YDR529C
··················│··│········│·····└──YPR149W
··················│··│········└──┬--YDR178W
··················│··│···········└──YEL011w
··················│··└──┬--┬--┬--YHR051W
··················│·····│··│··└──┬--YNL052W
··················│·····│··│·····└──YNR001C
··················│·····│··└──┬--YKL141W
··················│·····│·····└──YOR065W
··················│·····└──┬--YBL100C
··················│········└──YPR184W
··················└──┬--┬--YIL162W
·····················│··└──┬--┬--┬--YKL193C
·····················│·····│··│··└──┬--YIL113W
·····················│·····│··│·····└──YMR170C
·····················│·····│··└──┬--┬--YOL032W
·····················│·····│·····│··└──┬--YKL151C
·····················│·····│·····│·····└──┬--YDR272W
·····················│·····│·····│········└──YCL035C
·····················│·····│·····└──┬--YFL054C
·····················│·····│········└──┬--YNL274C
·····················│·····│···········└──┬--YLR270W
·····················│·····│··············└──┬--YDR272W
·····················│·····│·················└──YHR104W
·····················│·····└──┬--YLR356W
·····················│········└──YBR241C
·····················└──┬--YDL204w
························└──┬--┬--┬--YBL015W
···························│··│··└──┬--YKL217W
···························│··│·····└──YMR107W
···························│··└──┬--YMR191W
···························│·····└──┬--YKL109W
···························│········└──YML054C
···························└──┬--YGL191W
······························└──┬--┬--┬--YGR243W
·································│··│··└──┬--YBL048W
·································│··│·····└──YGR236C
·································│··└──┬--YDR070c
·································│·····└──┬--YBR147W
·································│········└──┬--YNL134C
·································│···········└──YNL194C
·································└──┬--┬--┬--YBL064C
····································│··│··└──YGR043C
····································│··└──┬--┬--YDR171W
····································│·····│··└──YBL078C
····································│·····└──┬--YKL026C
····································│········└──┬--YOR215C
····································│···········└──┬--YFR033C
····································│··············└──YGR088W
····································└──┬--YER067w
·······································└──┬--YDR533C
··········································└──YOR178C

Hierarchical Tree