第14章 聚类方法

习题14.1

  试写出分裂聚类算法，自上而下地对数据进行聚类，并给出其算法复杂度。

解答：

解答思路：

1. 给出一般聚类方法概述
2. 给出分裂聚类的定义
3. 给出分裂聚类预先确定的要素
4. 写出分裂聚类算法，并计算复杂度
5. 自编程实现分裂聚类算法

解答步骤：

第1步：一般聚类方法概述

  根据书中第14章的聚类方法介绍：

  聚类是针对给定的样本，依据它们的特征的相似度或距离，将其归并到若干个“类”或“簇”的数据分析问题。一个类是给定样本集合的一个子集。直观上，相似的样本聚集在相同的类，不相似的样本分散在不同的类。常用的聚类算法包括层次聚类和$k$均值聚类。层次聚类又分为聚合（自下而上）和分裂（自上而下）两种方法。

第2步：分裂聚类的定义

  根据书中第14章的分裂聚类介绍：

  分裂法开始将所有样本分到一个类，之后将已有类中相距最远的样本分到两个新的类，重复此操作直到满足停止条件，得到层次化的类别。

第3步：分裂聚类预先确定的要素

  分裂聚类的具体过程如下：对于给定的样本集合，开始将所有样本分到一个类，然后按照一定分裂方法，例如类中距离最大，将最满足规则条件的样本分到两个新的类；如此反复进行，每次减少剩余未分类的样本，直到满足停止条件，如所有样本都分到了对应的类中。

  由此可知，分裂聚类需要预先确定下面三个要素：
  （1）选择要分裂的类方法；
  （2）分裂类的方法；
  （3）停止条件。

  根据这些要素的不同组合，就可以构成不同的分裂聚类方法。选择要分裂的类的方法可以是选择最大直径的类、样本最多的类、最大均值的类、最大SSE的类。分裂类的方法可以是使用贪心算法选择新类中心、K-means算法、Ward方法。停止条件可以是每个类包含一个样本、类的个数达到阈值。

第4步：写出分裂聚类算法，并计算复杂度

分裂聚类算法
输入：$n$个样本组成的样本集合及设定的样本类别数$k$；
输出：满足设定的样本类别数的样本集合的一个层次化聚类。
（1）构造1个类，该类包含全部样本。
（2）计算$n$个样本两两之间的欧氏距离$\{d_{ij}\}$。
（3）分裂类中距离最大的两个样本，将其分到两个类，并设置为各自的类中心。
（4）计算未分裂的样本与目前各个类中心的距离，分裂类中距离最大的样本作为新的类中心，构造一个新类。
（5）如果类的个数满足设定的样本类别数$k$，则继续执行第（6）步，否则回到第（4）步继续分裂样本。
（6）根据当前的类中心，按照最近邻的原则对整个数据集进行分类。

计算复杂度：

1. 时间复杂度为$O(k n^2m)$：整个算法计算过程中需要计算两两样本之间的复杂度，设$m$是样本的维度，总体样本为$n$个，则时间复杂度为$O(n^2 m)$；设定的样本类别数为$k$，则需要重复$k$次，故时间复杂度为$O(k n^2 m)$。
2. 空间复杂度为$O(n^2)$：需要一个$n \times n$的矩阵保存两两样本间的距离。

第4步：自编程实现分裂聚类算法

import numpy as np

class DivisiveClustering:
    def __init__(self, num_class):
        # 聚类类别个数
        self.num_class = num_class
        # 聚类数据集
        self.cluster_data = []
        if num_class > 1:
            self.cluster_data = [[] for _ in range(num_class)]

    def fit(self, data):
        """
        :param data: 数据集
        """
        num_sample = data.shape[0]

        if self.num_class == 1:
            # 如果只设定了一类，将所有数据放入到该类中
            for d in data:
                self.cluster_data.append(d)
        else:
            # (1) 构造1个类，该类包含全部样本
            # 初始化类中心
            class_center = []

            # (2) 计算n个样本两两之间的欧氏距离
            distance = np.zeros((num_sample, num_sample))
            for i in range(num_sample):
                for j in range(i + 1, num_sample):
                    distance[j, i] = distance[i, j] = np.linalg.norm(data[i, :] - data[j, :], 
                                                                     ord=2)

            # (3) 分裂距离最大的两个样本，并设置为各自的类中心
            index = np.where(np.max(distance) == distance)
            class_1 = index[1][0]
            class_2 = index[1][1]
            # 记录已经分裂完成的样本
            finished_data = [class_1, class_2]
            class_center.append(data[class_1, :])
            class_center.append(data[class_2, :])

            num_class_temp = 2
            # (5) 判断类的个数是否满足设定的样本类别数
            while num_class_temp != self.num_class:
                # (4.1) 计算未分裂的样本与目前各个类中心的距离
                data2class_distance = np.zeros((num_sample, 1))
                for i in range(num_sample):
                    # 计算样本到各类中心的距离总和
                    data2class_sum = 0
                    for j, center_data in enumerate(class_center):
                        if i not in finished_data:
                            data2class_sum += np.linalg.norm(data[i, :] - center_data)
                    data2class_distance[i] = data2class_sum

                # (4.2) 分裂类间距离最大的样本作为新的类中心，构造一个新类
                class_new_index = np.argmax(data2class_distance)
                num_class_temp += 1
                finished_data.append(class_new_index)
                # 添加到类中心集合中
                class_center.append(data[class_new_index, :])

            # 根据当前的类中心，按照最近邻的原则对整个数据集进行分类
            for i in range(num_sample):
                data2class_distance = []
                for j, center_data in enumerate(class_center):
                    # 计算每个样本到类中心的距离
                    data2class_distance.append(np.linalg.norm(data[i, :] - center_data))

                # 将样本划分到最近的中心的类中
                label = np.argmin(data2class_distance)
                self.cluster_data[label].append(data[i, :])

# 使用书中例14.2的样本数据集
dataset = np.array([[0, 2],
                    [0, 0],
                    [1, 0],
                    [5, 0],
                    [5, 2]])

num_class = 2
divi_cluster = DivisiveClustering(num_class=num_class)
divi_cluster.fit(dataset)
print("分类数:", num_class)
for i in range(num_class):
    print(f"class_{i}:", divi_cluster.cluster_data[i])

分类数: 2
class_0: [array([0, 0]), array([1, 0]), array([5, 0])]
class_1: [array([0, 2]), array([5, 2])]

习题14.2

  证明类或簇的四个定义中，第一个定义可推出其他三个定义。

解答：

解答思路：

1. 给出定义14.5
2. 给出定义14.6，并由定义14.5推导定义14.6
3. 给出定义14.7，并由定义14.5推导定义14.7
4. 给出定义14.8，并由定义14.5推导定义14.8

解答步骤：

第1步：定义14.5

  根据书中第14.1.2节的定义14.5：

  用$G$表示类或簇，用$x_i,x_j$表示类中的样本，$N_G$表示$G$中样本的个数，用$d_{ij}$表示样本$x_i$与样本$x_j$之间的距离。

定义14.5 设$T$为给定的正数，若集合$G$中任意两个样本$x_i, x_j$，有

$$$$

则称$G$为一个类或簇。

  由定义14.5可知，$T$大于等于集合$G$中两两样本间的距离最大值：

  则集合$G$中任取两个样本$x_i, x_j$，满足

  定义14.6得证。

第3步：由定义14.5推导定义14.7

证明：

  根据书中第14.1.2节的定义14.7：

定义14.7 设$T$为给定的正数，若对集合$G$中任意一个样本$x_i$，$G$中的另一个样本$x_j$满足

$$$$

\sum_{x_j \in G} d_{ij} \leqslant (n_G - 1) T

\frac{1}{n_G - 1} \sum_{x_j \in G} d_{ij} \leqslant T

则称$G$为一个类或簇。

  根据定义14.5可知：$d_{ij} \leqslant T$，则不妨令$V=T$，对于$G$中任意两个样本$x_i, x_j$，则有$d_{ij} \leqslant V, d_{ij} \leqslant T$，结合定义14.7，则有

  所以

  定义14.8得证。

习题14.3

  证明式(14.21)成立，即$k$均值的可能解的个数是指数级的。

解答：

解答思路：

1. 给出书中公式（14.21）
2. 给出可区分球到不可区分盒子的分配策略，得到指数的生成函数
3. 利用指数生成函数进行二项式展开，推导并证明公式

解答步骤：

第1步：书中公式（14.21）

  根据书中第14.3.2节的描述：

  $k$均值聚类就是求解最优化问题：

$$$$

C^* &= \arg \min \limits_{C} W(C) \ &= \arg \min \limits_{C} \sum_{l=1}^k \sum_{C(i)=l} | x_i - \bar{x}_l |^2 \end{aligned}

这个数字是指数级的。

第2步：给出可区分球到不可区分盒子的分配策略

  以下示例来自于《应用组合数学（第2版）》第216页5.5.3节：

  $S(n, k)$的定义是把$n$个可区分球分配到$k$个可区分盒子里，其每一个盒子不为空时的方法数量。 首先，考虑把$n$个可区分球放入标有标签$1, 2, \cdots, k$的可区分盒子里，且每一个盒子都不为空的方法数量$T(n,k)$的问题，注意

$$$$

其中$T(n, k)$是$H(x)$的展开式中$\displaystyle \frac{x^n}{n!}$的系数。

  根据《应用组合数学（第2版）》第213页指数生成函数定义：

  记$P(n, k)$是$n$集合的$k$排序的数量，固定$n$，则$P(n, k)$的普通生成函数由下式给出：

$$$$

G(x) &= P(n, 0) x^0 + P(n, 1) x^1 + P(n, 2) x^2 + \cdots + P(n, n) x^n \ &= \sum_{k=0}^{\infty} \frac{(n-k)!}{n!} x^k \end{aligned}

通过二项式展开可以化简为$(1 + x)^n$，则有：

$$$$

其中$P(n, k)$是$(1 + x)^n$的展开式中$\displaystyle \frac{x^k}{k!}$的系数。
  如果$(a_k)$是任意序列，则一序列的指数生成函数是

$$$$

T(n, k) = k!S(n, k)

H(x) = \left( x + \frac{x^2}{2!} + \frac{x^3}{3!} + \cdots \right)^k = (e^x - 1)^k

用$(k-i) x$取代$x$，得到

$$$$

H(x) &= \sum_{i=0}^k \left( \begin{array}{c} k \ i \end{array} \right) (-1)^i \sum_{n=0}^{\infty} \frac{1}{n!} (k-i)^n x^n \ &= \sum_{n=0}^{\infty} \frac{x^n}{n!} \sum_{i=0}^k (-1)^i \left( \begin{array}{c} k \ i \end{array} \right) (k-i)^n \end{aligned}

结合$T(n, k) = k!S(n, k)$，可得

$$$$

S(n,k) = \frac{1}{k!} \sum_{l=0}^k (-1)^{k-l} \left( \begin{array}{c} k \ l \end{array} \right) l^n

就是说在类中心确定的情况下，将每个样本分到一个类中，使样本和其所属类的中心之间的距离总和最小。求解结果，将每个样本指派到与其最近的中心$m_l$的类$G_l$中。
  然后，对给定的划分$C$，再求各个类的中心$(m_1, m_2, \cdots , m_k)$，使得目标函数极小化：

$$$$

  重复以上两个步骤，直到划分不再改变，得到聚类结果。

第2步：高斯混合模型参数估计的EM算法

  根据书中第9章的定义9.2的高斯混合模型：

定义9.2（高斯混合模型） 高斯混合模型是指具有如下形式的概率分布模型：

$$$$

  根据书中第9.3.2节的高斯混合模型参数估计的EM算法：

假设观测数据$y_1, y_2, \cdots, y_N$由高斯混合模型生成，

$$$$

$\gamma_{jk}$是0-1随机变量。
完全数据的对数似然函数为

$$$$

3.EM算法的M步：计算新一轮迭代的模型参数

$$$$

\hat{\sigma}^2 = \frac{\displaystyle \sum_{j=1}^N \hat{\gamma}{jk} (y_j - u_k)^2}{\displaystyle \sum{j=1}^N \hat{\gamma}{jk}}， \quad k = 1,2,\cdots, K \ \hat{\alpha}k = \frac{\displaystyle \sum{j=1}^N \hat{\gamma}{jk} }{N}, \quad k = 1,2,\cdots, K