【学习笔记】统计学习方法——聚类方法

摘要： 聚类是针对给定的样本，依据它们特征的相似度或距离，将其归并到若干个“类”或“簇”的数据分析问题。一个类是样本的一个子集。直观上，相似的样本聚集在相同的类，不相似的样本分散在不同的类。这里，样本之间的相似度或距离起着重要作用。

聚类的基本概念

首先定义一个矩阵X用来表示n个样本的m个属性。
$$X=[x_{ij}{]}_{m\times n}=\left[\begin{array}{cccc}{x}_{11}& x_{12}& \cdots & x_{1n}\\ x_{21}& x_{22}& \cdots & x_{2n}\\ ⋮& ⋮& & ⋮\\ x_{m1}& x_{m2}& \cdots & x_{mn}\end{array}\right]$$
这里矩阵X的定义是按照书上来的，1列代表一个样本，当然还有另外一种常见的写法是1行代表一个样本，怎么写都行。

相似度或距离

1）闵可夫斯基距离（距离越大相似度越小，距离越小相似度越大）

给定样本集合X，X是m维实数向量空间$R^m$中点的集合，其中$x_i,x_j\in X，x_i=(x_{1i},x_{2i},\cdots ,x_{mi}{)}^T$，$x_j=(x_{1j},x_{2j},\cdots ,x_{mj}{)}^T$，样本$x_i$与样本$x_j$的闵可夫斯基距离定义为
$$d_{ij}={\left(\sum _{k=1}^m{\mid x_{ki}-x_{kj}\mid }^p\right)}^{\frac{1}{p}}$$

这里$p\ge 1$。当$p=2$时称为欧式距离，即
$$d_{ij}={\left(\sum _{k=1}^m{\mid x_{ki}-x_{kj}\mid }^2\right)}^{\frac{1}{2}}$$

当$p=1$时称为曼哈顿距离，即
$$d_{ij}=\sum _{k=1}^m|x_{ki}-x_{kj}|$$
当$p=\infty$时称为切比雪夫距离，取各个坐标数值差的绝对值的最大值，即
$$d_{ij}={\max }_k|x_{ki}-x_{kj}|$$
接下来举个简单的例子把，假设$x_1=(1,2)，x_2=(5,6)$

那么闵可夫斯基距离为：$\sqrt[p]{(5-1{)}^p+(6-2{)}^p}$

欧式距离为：$\sqrt[2]{(5-1{)}^p+(6-2{)}^2}=\sqrt[2]{32}$

曼哈顿距离为：$|5-1|+|6-2|=8$

切比雪夫距离为：$max(|5-1|,|6-2|)=4$

2）马哈拉诺比斯距离（距离越大相似度越小，距离越小相似度越大）

给定一个样本集合X，$X={(x_{ij})}_{m\times n}$，其协方差矩阵记作S。样本$x_i$与样本$x_j$之间的马哈拉诺比斯距离$d_{ij}$定义为
$$d_{ij}={\left[(x_i-x_j{)}^T{S}^{-1}(x_i-x_j)\right]}^{\frac{1}{2}}$$
其中
$$x_i=(x_{1i},x_{2i},\cdots ,x_{mi}{)}^T{x}_j=(x_{1j},x_{2j},\cdots ,x_{mj}{)}^T$$

当S为单位矩阵时，即样本数据的各个分量相互独立且各个分量的方差为1时，上式提到的马氏距离$d_{ij}$也就是欧式距离了，所以马氏距离是欧式距离的推广。

证明如下：

当S为单位矩阵时，那么$S^{-1}$也是单位矩阵，即$S=S^{-1}$。假设$x_i=(1,2,3),x_j=(5,4,6)$则：

$$\begin{gathered} d_{ij}={\left[{\left[\begin{array}{c}1-5\\ 2-4\\ 3-6\end{array}\right]}^T\left[\begin{array}{ccc}1& 0& 0\\ 0& 1& 0\\ 0& 0& 1\end{array}\right]\left[\begin{array}{c}1-5\\ 2-4\\ 3-6\end{array}\right]\right]}^{\frac{1}{2}} \\ =\sqrt[2]{(-4{)}^2+(-2{)}^2+(-3{)}^2} \end{gathered}$$
上式最后得到的结果也就是欧式距离的结果了。

现在来解释下什么叫做协方差矩阵：

方差和协方差的定义：

在统计学中，方差是用来度量单个随机变量的离散程度，而协方差则一般用来刻画两个随机变量的相似程度，其中，方差的计算公式为

${\sigma }_x^2=\frac{1}{n-1}{\sum }_{i=1}^n(x_i-\bar{x}{)}^2$

其中，$n$表示样本量，符号$\bar{x}$表示观测样本的均值。

在此基础上，协方差的计算公式被定义为

$\sigma (x,y)=\frac{1}{n-1}{\sum }_{i=1}^n\left(x_i-\bar{x}\right)\left(y_i-\bar{y}\right)$

在公式中，符号 $\bar{x},\bar{y}$ 分别表示两个随机变量所对应的观测样本均值，据此, 我们发现：方差 ${\sigma }_x^2$ 可视作随机变量 $x$ 关于其自身的协方差 $\sigma (x,x)$。

从方差/协方差到协方差矩阵：

根据方差的定义, 给定 $d$ 个随机变量 $x_k,k=1,2,\dots ,d,$ 则这些随机变量的方差为
$\sigma \left(x_k,x_k\right)=\frac{1}{n-1}{\sum }_{i=1}^n{\left(x_{ki}-{\bar{x}}_k\right)}^2,k=1,2,\dots ,d$

其中, 为方便书写, $x_{ki}$ 表示随机变量 $x_k$ 中的第 $i$ 个观测样本, $n$ 表示样本量，每个随机变量所对应的观测样本数量均为 $n$ 。

对于这些随机变量，我们还可以根据协方差的定义，求出两两之间的协方差, 即
$$\sigma \left(x_m,x_k\right)=\frac{1}{n-1}\sum _{i=1}^n\left(x_{mi}-{\bar{x}}_m\right)\left(x_{ki}-{\bar{x}}_k\right)$$
因此, 协方差矩阵为
$$\Sigma =\left[\begin{array}{ccc}\sigma \left(x_1,x_1\right)& \cdots & \sigma \left(x_1,x_d\right)\\ ⋮& \ddots & ⋮\\ \sigma \left(x_d,x_1\right)& \cdots & \sigma \left(x_d,x_d\right)\end{array}\right]\in {\mathbb{R}}^{d\times d}$$
其中，对角线上的元素为各个随机变量的方差, 非良角线上的元素为两两随机变量之间的协方差, 根据协方差的定义，我们可以认定：矩阵 $\Sigma$ 为对称矩阵(symmetric matrix)，其大小为 $d\times d$。

3）相关系数（相关系数的绝对值越接近于1，表示样本越相似；越接近于0，表示样本越不相似）

样本$x_i$与样本$x_j$之间的相关系数定义为
$$r_{ij}=\frac{{\sum }_{k=1}^m\left(x_{ki}-{\bar{x}}_i\right)\left(x_{kj}-{\bar{x}}_j\right)}{{\left[{\sum }_{k=1}^m{\left(x_{ki}-{\bar{x}}_i\right)}^2{\sum }_{k=1}^m{\left(x_{kj}-{\bar{x}}_j\right)}^2\right]}^{\frac{1}{2}}}$$
其中
$${\bar{x}}_i=\frac{1}{m}\sum _{k=1}^m{x}_{ki},{\bar{x}}_j=\frac{1}{m}\sum _{k=1}^m{x}_{kj}$$

4）夹角余弦（夹角余弦越接近于1，表示样本越相似；越接近于0，表示样本越不相似）

样本$x_i$与样本$x_j$之间的夹角余弦定义为
$$s_{ij}=\frac{{\sum }_{k=1}^m{x}_{ki}{x}_{kj}}{{\left[{\sum }_{k=1}^m{x}_{ki}^2{\sum }_{k=1}^m{x}_{kj}^2\right]}^{\frac{1}{2}}}$$

类或簇

类或簇的四个定义：

1）设 $T$ 为给定的正牧, 若集合 G 中任意两个样本 $x_i,x_j,$ 有
$$d_{ij}⩽T$$
则称G为一个类或簇。

2）设 T 为给定的正牧，若对集合 G 的任意样本 $x_i,$ 一定存在 G 中的另 个样本 $x_j,$ 使得
$$d_{ij}⩽T$$
则称G为一个类或簇。

3）设 $T$ 为给定的正数，若对集合 $G$ 中任意一个样本 $x_i,G$ 中的另一个 样本 $x_j$ 满足
$$\frac{1}{n_G-1}\sum _{x_j\in G}{d}_{ij}⩽T$$
$\text{ 其中 }{n}_G\text{ 为 }G\text{ 中样本的个数, 则称 }G\text{ 为一个类或籍。 }$

4）设 T 和 V 为给定的两个正牧，如果集合 G 中任意两个样本 $x_i,x_j$ 的 距离 $d_{ij}$ 满足
$$\begin{array}{l}\frac{1}{n_G\left(n_G-1\right)}{\sum }_{x_i\in G}{\sum }_{x_j\in G}{d}_{ij}⩽T\\ d_{ij}⩽V\end{array}$$
则称 G 为一个类或籍。

类的三个特征：

1）类的均值 ${\bar{x}}_G$, 又称为类的中心
$${\bar{x}}_G=\frac{1}{n_G}\sum _{i=1}^{n_G}{x}_i$$
式中 $n_G$ 是类 $G$ 的样本个数。

2）类的直径 (diameter) $D_G$
类的直径 D_ 是类中任意两个样本之间的最大距离, 即
$$D_G=\underset{x_i,x_j\in G}{\max }{d}_{ij}$$
3）类的样本散布矩阵$A_G$与样本协方差矩阵$S_G$

类的样本散布矩阵$A_G$为
$$A_G=\sum _{i=1}^{n_G}\left(x_i-{\bar{x}}_G\right){\left(x_i-{\bar{x}}_G\right)}^{\mathrm{T}}$$
样本协方差矩阵 SG 为
$$\begin{array}{rl}{S}_G& =\frac{1}{m-1}{A}_G\\ & =\frac{1}{m-1}\sum _{i=1}^{n_C}\left(x_i-{\bar{x}}_G\right){\left(x_i-{\bar{x}}_G\right)}^{\mathrm{T}}\end{array}$$
其中 m 为样本的维数（样本属性的个数）。

类与类之间的距离

设类 $G_p$ 包含 $n_p$ 个样本, $G_q$ 包含 $n_q$ 个样本, 分别用 ${\bar{x}}_p$ 和 ${\bar{x}}_q$ 表示 $G_p$ 和 $G_q$ 的 均值，即类的中心。

1）最短距离或单连接

定义类 $G_p$ 的样本与 $G_q$ 的样本之间的最短距离为两类之间的距离
D p q = min ⁡ { d i j ∣ x i ∈ G p , x j ∈ G q } D_{p q}=\min \left\{d_{i j} \mid x_{i} \in G_{p}, x_{j} \in G_{q}\right\} Dpq​=min{dij​∣xi​∈Gp​,xj​∈Gq​}
2）最长距离或完全连接

定义类 $G_p$ 的样本与 $G_q$ 的样本之间的最长距离为两类之间的距离
D p q = max ⁡ { d i j ∣ x i ∈ G p , x j ∈ G q } D_{p q}=\max \left\{d_{i j} \mid x_{i} \in G_{p}, x_{j} \in G_{q}\right\} Dpq​=max{dij​∣xi​∈Gp​,xj​∈Gq​}
3）中心距离

定义类 $G_p$ 与类 $G_q$ 的中心 ${\bar{x}}_p$ 与 ${\bar{x}}_q$ 之间的距离为两类之同的距离
$$D_{pq}=d_{x_p{x}_q}$$
（4）半均距隨
定义类 $G_p$ 与类 $G_q$ 任套两个样本之但距离的平均值为两类之间的距离
$$D_{pq}=\frac{1}{n_p{n}_q}\sum _{x_i\in G_p}\sum _{x_j\in G_q}{d}_{ij}$$

聚类算法

层次聚类

输入: $n$ 个样本组成的样本集合及样本之间的距离： 输出: 对样本集合的一个层次化聚类。

（1）计算 n 个样本两两之间的欧氏距离 { d i j } , \left\{d_{i j}\right\}, {dij​}, 记作矩阵 $D={\left[d_{ij}\right]}_{n\times n}$
（2）构造 n 个类，每个类只包含一个样本。
（3）合并类间距离最小的两个类，其中最短距离为类间距离，构建一个新类。
（4）计算新类与当前各类的距离。若类的个数为 1，终止计算, 否则回到步 (3)。

可以看出聚合层次政类算法的复杂度是 $O\left(n^{3}m\right),$ 其中 $m$ 是样本的维数，n 是样本个数。

k均值聚类

输入: $n$ 个样本的集合 X:

输出: 样本集合的款类 C”。

（1）初始化。令 $t=0,$ 随 机选择 $k$ 个样本点作为初始聚类中心 $m^{(0)}=$ $\left(m_1^{(0)},\cdots ,m_l^{(0)},\cdots ,m_k^{(0)}\right)$
（2）对样本进行麥类。对固定的类中心 $m^{(t)}=\left(m_1^{(t)},\cdots ,m_l^{(t)},\cdots ,m_k^{(t)}\right),$ 其中
$m_l^{(t)}$ 为类 $G_l$ 的中心，计算每个样本到类中心的距离，将每个样本指派到与其最近的 中心的类中，构成聚类结果 $C^{(t)}$
（3）计算新的类中心。对聚类结果 C ${}^{(t)},$ 计算当前各个类中的样本的均值，作为新
的类中心 $m^{(t+1)}=\left(m_1^{(t+1)},\cdots ,m_l^{(t+1)},\cdots ,m_k^{(t+1)}\right)$
（4）如果迭代收琉或符合停止条件，输出 $C^{\ast }=C^{(t)}$。否则， $\iff t=t+1,$ 返回步 (2)

k 均值聚类算法的复杂度是 $O(mnk),$ 其中 $m$ 是样本维数，n 是样本个数，k 是类别个数。

代码

# 定义聚类数的节点

class ClusterNode:
    def __init__(self, vec, left=None, right=None, distance=-1, id=None, count=1):
        """
        :param vec: 保存两个数据聚类后形成新的中心
        :param left: 左节点
        :param right:  右节点
        :param distance: 两个节点的距离
        :param id: 用来标记哪些节点是计算过的
        :param count: 这个节点的叶子节点个数
        """
        self.vec = vec
        self.left = left
        self.right = right
        self.distance = distance
        self.id = id
        self.count = count
        
 def euler_distance(point1: np.ndarray, point2: list) -> float:
    """
    计算两点之间的欧拉距离，支持多维
    """
    distance = 0.0
    for a, b in zip(point1, point2):
        distance += math.pow(a - b, 2)
    return math.sqrt(distance)

# 层次聚类（聚合法）

class Hierarchical:
    def __init__(self, k):
        self.k = k
        self.labels = None
        
    def fit(self, x):
        nodes = [ClusterNode(vec=v, id=i) for i, v in enumerate(x)]
        distances = {}
        point_num, feature_num = x.shape
        self.labels = [-1] * point_num
        currentclustid = -1
        while(len(nodes)) > self.k:
            min_dist = math.inf
            nodes_len = len(nodes)
            closest_part = None
            for i in range(nodes_len - 1):
                for j in range(i+1, nodes_len):
                    d_key = (nodes[i].id, nodes[j].id)
                    if d_key not in distances:
                        distances[d_key] = euler_distance(nodes[i].vec, nodes[j].vec)
                    d = distances[d_key]
                    if d < min_dist:
                        min_dist = d
                        closest_part = (i, j)
                        
            part1, part2 = closest_part
            node1, node2 = nodes[part1], nodes[part2]
            new_vec = [ (node1.vec[i] * node1.count + node2.vec[i] * node2.count ) / (node1.count + node2.count)
                        for i in range(feature_num)]
            new_node = ClusterNode(vec=new_vec,
                                   left=node1,
                                   right=node2,
                                   distance=min_dist,
                                   id=currentclustid,
                                   count=node1.count + node2.count)
            currentclustid -= 1
            del nodes[part2], nodes[part1]
            nodes.append(new_node)
            
        self.nodes = nodes
        self.calc_label()
        
    def calc_label(self):
        """
        调取聚类的结果
        """
        for i, node in enumerate(self.nodes):
            # 将节点的所有叶子节点都分类
            self.leaf_traversal(node, i)

    def leaf_traversal(self, node: ClusterNode, label):
        """
        递归遍历叶子节点
        """
        if node.left == None and node.right == None:
            self.labels[node.id] = label
        if node.left:
            self.leaf_traversal(node.left, label)
        if node.right:
            self.leaf_traversal(node.right, label)
            
# kmeans

class MyKmeans:
    def __init__(self, k, n=20):
        self.k = k
        self.n = n
        
    def fit(self, x, centers=None):
        # 第一步，随机选择 K 个点, 或者指定
        if centers is None:
            idx = np.random.randint(low=0, high=len(x), size=self.k)
            centers = x[idx]
        #print(centers)
        
        inters = 0
        while inters < self.n:
            #print(inters)
            #print(centers)
            points_set = {key: [] for key in range(self.k)}

            # 第二步，遍历所有点 P，将 P 放入最近的聚类中心的集合中
            for p in x:
                nearest_index = np.argmin(np.sum((centers - p) ** 2, axis=1) ** 0.5)
                points_set[nearest_index].append(p)

            # 第三步，遍历每一个点集，计算新的聚类中心
            for i_k in range(self.k):
                centers[i_k] = sum(points_set[i_k])/len(points_set[i_k])
                
            inters += 1

        
            
        return points_set, centers