• 编写层次聚类算法
    • 从文件中读取数据
    • 初始化优先队列
      • 距离相等的问题以及为何要使用元组
      • 距离相等的另一个问题
    • 重复下述步骤，直到仅剩一个分类

    编写层次聚类算法

    编写层次聚类算法 - 图1

    我们可以使用优先队列来实现这个聚类算法。

    什么是优先队列呢？

    普通的队列有“先进先出”的规则，比如向队列先后添加Moa、Suzuka、Yui，取出时得到的也是Moa、Suzuka、Yui：

    编写层次聚类算法 - 图2

    而对于优先队列，每个元素都可以附加一个优先级，从队列中取出时会得到优先级最高的元素。比如说，我们定义年龄越小优先级越高，以下是插入过程：

    编写层次聚类算法 - 图3

    取出的第一个元素是Yui，因为她的年龄最小：

    编写层次聚类算法 - 图4

    我们看看Python中如何使用优先队列：

    1. >>> from Queue import PriorityQueue           # 加载优先队列类
    2. >>> singersQueue = PriorityQueue()            # 创建对象
    3. >>> singersQueue.put((16, 'Suzuka Nakamoto')) # 插入元素
    4. >>> singersQueue.put((15, 'Moa Kikuchi'))
    5. >>> singersQueue.put((14, 'Yui Mizuno'))
    6. >>> singersQueue.put((17, 'Ayaka Sasaki'))
    7. >>> singersQueue.get() # 获取第一个元素，即最年轻的歌手Yui。
    8. (14, 'Yui Mizuno')
    9. >>> singersQueue.get()
    10. (15, 'Moa Kikuchi')
    11. >>> singersQueue.get()
    12. (16, 'Suzuka Nakamoto')
    13. >>> singersQueue.get()
    14. (17, 'Ayaka Sasaki')

    在进行聚类时，我们将分类、离它最近的分类、以及距离插入到优先队列中，距离作为优先级。比如上面的犬种示例，Border Collie最近的分类是Portuguese WD，距离是0.232：

    编写层次聚类算法 - 图5

    我们将优先队列中距离最小的两个分类取出来，合并成一个分类，并重新插入到优先队列中。比如下图是将Border Collie和Portuguese WD合并后的结果：

    编写层次聚类算法 - 图6

    重复这个过程，直到队列中只有一个元素为止。当然，我们插入的数据会复杂一些，请看下面的讲解。

    从文件中读取数据

    数据文件是CSV格式的（以逗号分隔），第一行是列名，第一列是犬种，第二列之后是特征值：

    编写层次聚类算法 - 图7

    我们用Python的列表结构来存储这些数据，data[0]用来存放所有记录的分类，如data[0][0]是Border Collie，data[0][1]是Boston Terrier。data[1]则是所有记录的高度，data[2]是重量。

    特征列的数据都会转换成浮点类型，如data[1][0]是20.0，data[2][0]是45.0等。在读取数据时就需要对其进行标准化。此外，我们接下来会使用“下标”这个术语，如第一条记录Border Collie的下标是0，第二条记录Boston Terrier下标是1等。

    初始化优先队列

    以Border Collie为例，我们需要计算它和其它犬种的距离，保存在Python字典里：

    1. {1: ((0, 1), 1.0244),  # Border Collie（下标为0）和Boston Terrier（下标为1）之间的距离为1.0244
    2.  2: ((0, 2), 0.463),   # Border Collie和Brittany Spaniel（下标为2）之间的距离为0.463
    3.  ...
    4.  10: ((0, 10), 2.756)} # Border Collie和Yorkshire Terrier的距离为2.756

    此外，我们会记录Border Collie最近的分类及距离：这对犬种是(0, 8)，即下标为0的Border Collie和下标为8的Portuguese WD，距离是0.232。

    距离相等的问题以及为何要使用元组

    你也许注意到了，Portuguese WD和Standard Poodle的距离是0.566，Boston Terrier和Brittany Spaniel的距离也是0.566，

    如果我们通过最短距离来取，很可能会取出Standard Poodle和Boston Terrier进行组合，这显然是错误的，所以我们才会使用元组来存放这对犬种的下标，以作判断。比如说，Portuguese WD的记录是：

    1. ['Portuguese Water Dog', 0.566, (8, 9)]

    它的近邻Standard Poodle的记录是：

    1. ['Standard Poodle', 0.566, (8, 9)]

    我们可以通过这个元组来判断这两条记录是否是一对。

    距离相等的另一个问题

    在介绍优先队列时，我用了歌手的年龄举例，如果他们的年龄相等，取出的顺序又是怎样的呢？

    编写层次聚类算法 - 图8

    可以看到，如果年龄相等，优先队列会根据记录中的第二个元素进行判断，即歌手的姓名，并按字母顺序返回，如Avaka会比Moa优先返回。

    在犬种示例中，我们让距离成为第一优先级，下标成为第二优先级。因此，我们插入到优先队列的一条完整记录是这样的：

    编写层次聚类算法 - 图9

    重复下述步骤，直到仅剩一个分类

    我们从优先队列中取出两个元素，对它们进行合并。如合并Border Collie和Portuguese WD后，会形成一个新的分类：

    1. ['Border Collie', 'Portuguese WD']

    然后我们需要计算新的分类和其它分类之间的距离，方法是对取出的两个分类的距离字典进行合并。如第一个分类的距离字段是distanceDict1，第二个分类的是distanceDict2，新的距离字段是newDistanceDict：

    1. 初始化newDistanceDict
    2. 对于distanceDict1的每一个键值对:
    3.     如果这个键在distanceDict2中存在:
    4.         如果这个键在distanceDict1中的距离要比在distanceDict2中的距离小:
    5.             将distanceDict1中的距离存入newDistanceDict
    6.         否则:
    7.             将distanceDict2中的距离存入newDistanceDict

    编写层次聚类算法 - 图10

    经过计算后，插入到优先队列中的新分类的完整记录是：

    编写层次聚类算法 - 图11

    代码实践

    你能将上面的算法用Python实现吗？你可以从hierarchicalClustererTemplate.py这个文件开始，完成以下步骤：

    1. 编写init方法，对于每条记录：
      1. 计算该分类和其它分类之间的欧几里得距离；
      2. 找出该分类的近邻；
      3. 将这些信息放到优先队列的中。
    2. 编写cluster方法，重复以下步骤，直至剩下一个分类：
      1. 从优先队列中获取两个元素；
      2. 合并；
      3. 将合并后的分类放回优先队列中。

    编写层次聚类算法 - 图12

    解答

    注意，我的实现并不一定是最好的，你可以写出更好的！

    1. from queue import PriorityQueue
    2. import math
    4. """
    5. 层次聚类示例代码
    6. """
    8. def getMedian(alist):
    9.     """计算中位数"""
    10.     tmp = list(alist)
    11.     tmp.sort()
    12.     alen = len(tmp)
    13.     if (alen % 2) == 1:
    14.         return tmp[alen // 2]
    15.     else:
    16.         return (tmp[alen // 2] + tmp[(alen // 2) - 1]) / 2
    18. def normalizeColumn(column):
    19.     """计算修正的标准分"""
    20.     median = getMedian(column)
    21.     asd = sum([abs(x - median) for x in column]) / len(column)
    22.     result = [(x - median) / asd for x in column]
    23.     return result
    25. class hClusterer:
    26.     """该聚类器默认数据的第一列是标签，其它列是数值型的特征。"""
    28.     def __init__(self, filename):
    29.         file = open(filename)
    30.         self.data = {}
    31.         self.counter = 0
    32.         self.queue = PriorityQueue()
    33.         lines = file.readlines()
    34.         file.close()
    35.         header = lines[0].split(',')
    36.         self.cols = len(header)
    37.         self.data = [[] for i in range(len(header))]
    38.         for line in lines[1:]:
    39.             cells = line.split(',')
    40.             toggle = 0
    41.             for cell in range(self.cols):
    42.                 if toggle == 0:
    43.                    self.data[cell].append(cells[cell])
    44.                    toggle = 1
    45.                 else:
    46.                     self.data[cell].append(float(cells[cell]))
    47.         # 标准化特征列（即跳过第一列）
    48.         for i in range(1, self.cols):
    49.                 self.data[i] = normalizeColumn(self.data[i])
    51.         ###
    52.         ###  数据已经读入内存并做了标准化，对于每一条记录，将执行以下步骤：
    53.         ###     1. 计算该分类和其他分类的距离，如当前分类的下标是1，
    54.         ###        它和下标为2及下标为3的分类之间的距离用以下形式表示：
    55.         ###        {2: ((1, 2), 1.23),  3: ((1, 3), 2.3)... }
    56.         ###     2. 找出距离最近的分类；
    57.         ###     3. 将该分类插入到优先队列中。
    58.         ###
    60.         # 插入队列
    61.         rows = len(self.data[0])              
    63.         for i in range(rows):
    64.             minDistance = 99999
    65.             nearestNeighbor = 0
    66.             neighbors = {}
    67.             for j in range(rows):
    68.                 if i != j:
    69.                     dist = self.distance(i, j)
    70.                     if i < j:
    71.                         pair = (i,j)
    72.                     else:
    73.                         pair = (j,i)
    74.                     neighbors[j] = (pair, dist)
    75.                     if dist < minDistance:
    76.                         minDistance = dist
    77.                         nearestNeighbor = j
    78.                         nearestNum = j
    79.             # 记录这两个分类的配对信息
    80.             if i < nearestNeighbor:
    81.                 nearestPair = (i, nearestNeighbor)
    82.             else:
    83.                 nearestPair = (nearestNeighbor, i)
    85.             # 插入优先队列
    86.             self.queue.put((minDistance, self.counter,
    87.                             [[self.data[0][i]], nearestPair, neighbors]))
    88.             self.counter += 1
    90.     def distance(self, i, j):
    91.         sumSquares = 0
    92.         for k in range(1, self.cols):
    93.             sumSquares += (self.data[k][i] - self.data[k][j])**2
    94.         return math.sqrt(sumSquares)
    96.     def cluster(self):
    97.          done = False
    98.          while not done:
    99.              topOne = self.queue.get()
    100.              nearestPair = topOne[2][1]
    101.              if not self.queue.empty():
    102.                  nextOne = self.queue.get()
    103.                  nearPair = nextOne[2][1]
    104.                  tmp = []
    106.                  ##  我从队列中取出了两个元素：topOne和nextOne，
    107.                  ##  检查这两个分类是否是一对，如果不是就继续从优先队列中取出元素，
    108.                  ##  直至找到topOne的配对分类为止。
    109.                  while nearPair != nearestPair:
    110.                      tmp.append((nextOne[0], self.counter, nextOne[2]))
    111.                      self.counter += 1
    112.                      nextOne = self.queue.get()
    113.                      nearPair = nextOne[2][1]
    115.                  ## 将不处理的元素退回给优先队列
    116.                  for item in tmp:
    117.                      self.queue.put(item)
    119.                  if len(topOne[2][0]) == 1:
    120.                      item1 = topOne[2][0][0]
    121.                  else:
    122.                      item1 = topOne[2][0]
    123.                  if len(nextOne[2][0]) == 1:
    124.                      item2 = nextOne[2][0][0]
    125.                  else:
    126.                      item2 = nextOne[2][0]
    127.                  ##  curCluster即合并后的分类
    128.                  curCluster = (item1, item2)
    130.                  ## 对于这个新的分类需要做两件事情：首先找到离它最近的分类，然后合并距离字典。
    131.                  ## 如果item1和元素23的距离是2，item2和元素23的距离是4，我们取较小的那个距离，即单链聚类。
    132.                  minDistance = 99999
    133.                  nearestPair = ()
    134.                  nearestNeighbor = ''
    135.                  merged = {}
    136.                  nNeighbors = nextOne[2][2]
    137.                  for (key, value) in topOne[2][2].items():
    138.                     if key in nNeighbors:
    139.                         if nNeighbors[key][1] < value[1]:
    140.                              dist =  nNeighbors[key]
    141.                         else:
    142.                             dist = value
    143.                         if dist[1] < minDistance:
    144.                              minDistance =  dist[1]
    145.                              nearestPair = dist[0]
    146.                              nearestNeighbor = key
    147.                         merged[key] = dist
    149.                  if merged == {}:
    150.                     return curCluster
    151.                  else:
    152.                     self.queue.put( (minDistance, self.counter,
    153.                                      [curCluster, nearestPair, merged]))
    154.                     self.counter += 1
    156. def printDendrogram(T, sep=3):
    157.     """打印二叉树状图。树的每个节点是一个二元组。这个方法摘自：
    158.     http://code.activestate.com/recipes/139422-dendrogram-drawing/"""
    160.     def isPair(T):
    161.         return type(T) == tuple and len(T) == 2
    163.     def maxHeight(T):
    164.         if isPair(T):
    165.             h = max(maxHeight(T[0]), maxHeight(T[1]))
    166.         else:
    167.             h = len(str(T))
    168.         return h + sep
    170.     activeLevels = {}
    172.     def traverse(T, h, isFirst):
    173.         if isPair(T):
    174.             traverse(T[0], h-sep, 1)
    175.             s = [' ']*(h-sep)
    176.             s.append('|')
    177.         else:
    178.             s = list(str(T))
    179.             s.append(' ')
    181.         while len(s) < h:
    182.             s.append('-')
    184.         if (isFirst >= 0):
    185.             s.append('+')
    186.             if isFirst:
    187.                 activeLevels[h] = 1
    188.             else:
    189.                 del activeLevels[h]
    191.         A = list(activeLevels)
    192.         A.sort()
    193.         for L in A:
    194.             if len(s) < L:
    195.                 while len(s) < L:
    196.                     s.append(' ')
    197.                 s.append('|')
    199.         print (''.join(s))    
    201.         if isPair(T):
    202.             traverse(T[1], h-sep, 0)
    204.     traverse(T, maxHeight(T), -1)
    207. filename = '/Users/raz/Dropbox/guide/data/dogs.csv'
    209. hg = hClusterer(filename)
    210. cluster = hg.cluster()
    211. printDendrogram(cluster)

    运行结果和我们手算的一致：

    编写层次聚类算法 - 图13

    动手实践

    这里提供了77种早餐麦片的营养信息，包括以下几项：

    • 麦片名称
    • 热量
    • 蛋白质
    • 脂肪
    • 纳
    • 纤维
    • 碳水化合物
    • 糖
    • 钾
    • 维生素

    编写层次聚类算法 - 图14

    请对这个数据集进行层次聚类：

    • 哪种麦片和Trix最相近？
    • 与Muesli Raisins & Almonds最相近的是？

    数据集来自：http://lib.stat.cmu.edu/DASL/Datafiles/Cereals.html

    结果

    我们只需将代码中的文件名替换掉就可以了，结果如下：

    编写层次聚类算法 - 图15

    因此Trix和Fruity Pebbles最相似（你可以去买这两种麦片尝尝）。Muesli Raisins & Almonds和Muesli Peaches & Pecans最相似。

    编写层次聚类算法 - 图16

    好了，这就是层次聚类算法，很简单吧！