kmeans 算法原理与调整兰德系数(ARI)

452 阅读2分钟
1.基本原理

是一个比较简单的算法

随机抽取k个点作为最初的质心,进入循环

1.将每个样本点分配到离它最近的质心,生成k个簇

2.对于每个簇,计算所有被分配到该点的样本的新的质心,新质心取簇的平均坐标位置

3.当质心的位置变化小于某个数值的时候,迭代停止,跳出循环

from sklearn.datasets import make_blobs
from sklearn.cluster import KMeans
import numpy as np


class Kmeans:
    def __init__(self, k, max_iter, eps):
        self.k = int(k)
        self.max_iter = int(max_iter)
        self.eps = float(eps)
        self.features = None
        self.centers = None
        self.labels = None
        self.m, self.n = 0, 0

    def fit(self, features):
        self.features = features
        self.m, self.n = features.shape
        self.centers = np.array(self.features[np.random.choice(self.m, self.k, replace=False)])

        for _ in range(self.max_iter):
            distances = np.sqrt((self.features - self.centers[:, np.newaxis]) ** 2).sum(axis=2)
            self.labels = np.argmin(distances, axis=0)
            new_centers = np.array([self.features[self.labels == k].mean(axis=0) for k in range(self.k)])

            if np.all((abs(new_centers - self.centers) < self.eps)):
                break
            self.centers = new_centers

    def get_args(self):
        return self.centers, self.labels


# making data
x, y = make_blobs(
    n_samples=500,
    n_features=2,
    centers=4,
    random_state=1
)

mk = Kmeans(4, 1000, 1e-5)
mk.fit(x)
center, label = mk.get_args()

print(center)

# compare with sklearn.cluster.KMeans
kmeans = KMeans(n_clusters=4, n_init=10)
kmeans.fit(x)
centers = kmeans.cluster_centers_
print(centers)

机器学习中的常见距离

欧几里得距离:

d(x,c)=i=1m(xici)2d(x,c) = \sqrt{\sum\limits_{i=1}^{m}(x_i - c_i)^2}

曼哈顿距离:

d(x,c)=i=1mxicid(x,c) = \sum\limits_{i=1}^{m}|x_i - c_i|

切比雪夫距离:

d(x,c)=maxi=1mxicid(x,c) = \max\limits_{i=1}^{m}|x_i - c_i|

额大概都用的是欧式距离吧...

2. 常用参数
1.n_clusters

表示需要得到k个簇,找到目标的位置

2.n_init

n_init 参数指定了算法运行的次数,具体来说,是用不同的随机质心初始化运行算法的次数

这个参数控制算法用不同的随机质心初始化并运行的次数。由于K-means算法对初始质心的选择非常敏感,不同的初始化可能会导致不同的聚类结果。n_init参数允许算法多次运行(每次都使用不同的随机初始化质心),然后选择最佳输出(即具有最小inertia_属性的输出)

3.init

选择初始化聚类中心的方式

k-means++, random 等等,k-means++random 更加均匀

4. max_iter

最大迭代次数,如果到达某个循环次数还未到达指定误差就终止

5.tol

表示的是可容忍的误差范围,有的地方也会写成eps

3.对无监督学习的结果进行评价

调整兰德指数(Adjusted Rand Index)可以表示

比如

from sklearn.metrics import adjusted_rand_score

true_labels = [0, 0, 1, 1, 2, 2]
clustered_labels = [1, 1, 0, 0, 2, 2]

# 计算调整后的调整兰德指数
ari = adjusted_rand_score(true_labels, clustered_labels)
print("Adjusted Rand Index:", ari)

在之前写的代码中可以用以下方式表示

from sklearn.metrics import adjusted_rand_score

ari = adjusted_rand_score(label, kmeans.labels_)
print("Adjusted Rand Index:", ari)

兰德指数的数学公式,同样应该是很好代码实现的

  1. 真正相同(True Positive, TP):A和B在真实聚类中属于同一类,在算法产生的聚类中也属于同一类。
  2. 真正不同(True Negative, TN):A和B在真实聚类中属于不同类,在算法产生的聚类中也属于不同类。
  3. 假正相同(False Positive, FP):A和B在真实聚类中属于不同类,但在算法产生的聚类中属于同一类。
  4. 假正不同(False Negative, FN):A和B在真实聚类中属于同一类,在算法产生的聚类中属于不同类。

对于非二分类问题也可以使用混淆矩阵实现目标

RI=TP+TNTP+TN+FP+FNRI = \frac{TP + TN}{TP + TN + FP + FN}

from itertools import combinations

# 定义聚类结果
cluster_a = [1, 1, 0, 0, 2, 2]
cluster_b = [0, 0, 1, 1, 2, 2]

# 初始化计数器
TP = FP = FN = TN = 0

# 遍历所有可能的数据点对
for (i, j) in combinations(range(len(cluster_a)), 2):
    same_cluster_a = cluster_a[i] == cluster_a[j]
    same_cluster_b = cluster_b[i] == cluster_b[j]
    if same_cluster_a and same_cluster_b:
        TP += 1
    elif same_cluster_a and not same_cluster_b:
        FP += 1
    elif not same_cluster_a and same_cluster_b:
        FN += 1
    elif not same_cluster_a and not same_cluster_b:
        TN += 1

# 计算兰德指数
rand_index = (TP + TN) / (TP + FP + FN + TN)
print(TP, FP, FN, TN, rand_index)

调整兰德指数的公式

ARI=IE(I)MaxIE(I)ARI = \frac{I - E(I)}{Max I - E(I)}

其中的E(I)E(I) 表示RR期望

排列组合符号标记C(n,r)=(nr)=n!r!(nr)!C(n, r) = \binom{n}{r} = \frac{n!}{r!(n-r)!}

对于混淆矩阵CCCijC_{ij} 表示ii行,jj列的元素,也就是类别ii被分配到类别jj 的数量

I=ij(Cij2)I = \sum\limits_{ij}\binom{C_{ij}}{2}

E(I)=i(ai2)j(bj2)(N2)E(I) = \frac{\sum_{i}\binom{a_i}{2} \cdot\sum_{j}\binom{b_j}{2}}{\binom{N}{2}}

MaxI=12(i(ai2)+j(bj2))Max I = \frac{1}{2} (\sum_i\binom{a_i}{2} + \sum_j\binom{b_j}{2})

from scipy.special import comb
from sklearn.metrics import adjusted_rand_score
import numpy as np


def compute_confusion_matrix(true_labels, predicted_labels, num_classes):
    confusion_matrix = [[0 for _ in range(num_classes)] for _ in range(num_classes)]

    for true, pred in zip(true_labels, predicted_labels):
        confusion_matrix[true][pred] += 1

    return np.array(confusion_matrix)


def calculate_ari(confusion_matrix):
    # 求和项
    sum_comb_Cij = np.sum([comb(Cij, 2) for Cij in confusion_matrix.flatten()])
    ai = np.sum(confusion_matrix, axis=1)
    bj = np.sum(confusion_matrix, axis=0)
    sum_comb_ai = np.sum(comb(ai, 2))
    sum_comb_bj = np.sum(comb(bj, 2))
    N = np.sum(confusion_matrix)

    # 索引、期望索引和最大索引
    Index = sum_comb_Cij
    Expected_Index = sum_comb_ai * sum_comb_bj / comb(N, 2)
    Max_Index = 0.5 * (sum_comb_ai + sum_comb_bj)

    # ARI计算
    ARI = (Index - Expected_Index) / (Max_Index - Expected_Index)

    return ARI


# 示例使用
labels_true = [1, 1, 1, 2, 0, 0]
labels_pred = [0, 0, 1, 1, 2, 2]
ari = calculate_ari(compute_confusion_matrix(labels_true, labels_pred, 3))
print(f"sklearn ari:{adjusted_rand_score(labels_true, labels_pred)}")
print(f"test ari: {ari}")

4.如何选择合适的K值

inertia属性kmeans.inertia_是表示簇内方差的属性

1.Elbow方法

当k值到某个特定的值之后inertia属性减小的速度会突然减小,图形类似人的手肘,也名为肘部法

2. 轮廓分析法

样本点ii的轮廓系数定义为

a(i)a(i) 为点ii与其同一聚类中所有其他点之间的平均距离,它衡量了点ii与其聚类内部的紧密度

b(i)b(i)为点ii与下一个最近聚类中所有点之间的平均最小距离,它衡量了点ii与最近的另一个聚类的分离度

s(i)=b(i)a(i)max(a(i)b(i))s(i) = \frac{b(i) - a(i)}{max{(a(i) - b(i))}}

a(i)<b(i)a(i)<b(i) 的时候s(i)s(i)是正的,意味着样本更靠近于自己的聚类中心而不是其他的聚类的点

a(i)>b(i)a(i)>b(i) 的时候s(i)s(i) 是负的,意味着样本远离自己的聚类中心

a(i)=b(i)a(i)=b(i) 的时候s(i)=0s(i) = 0 ,意味着,样本在分界线上

当k值到某个特定的值之后轮廓系数会有一个局部最小值, 这个最小值对应的k系数就是最适合的k值

这玩意的范围在[1,1][-1,1] 之间

可以用以下代码表示当前数据的

from sklearn.cluster import KMeans
from sklearn.datasets import make_blobs
from sklearn.metrics import silhouette_score
import matplotlib.pyplot as plt

X, y = make_blobs(
    n_samples=500,
    n_features=2,
    centers=4,
    random_state=1
)

# 肘部法则
sse = []
for k in range(1, 11):
    kmeans = KMeans(n_clusters=k, n_init=10).fit(X)
    sse.append(kmeans.inertia_)
plt.plot(range(1, 11), sse, marker='o')
plt.title('Elbow Method')
plt.xlabel('Number of clusters')
plt.ylabel('SSE')
plt.show()

# 轮廓分析
silhouette_scores = []
for k in range(2, 11):  # 轮廓系数需要至少两个聚类
    kmeans = KMeans(n_clusters=k, random_state=42, n_init=10)
    kmeans.fit(X)
    score = silhouette_score(X, kmeans.labels_)
    silhouette_scores.append(score)
plt.plot(range(2, 11), silhouette_scores, marker='o')
plt.title('Silhouette Score')
plt.xlabel('Number of clusters')
plt.ylabel('Silhouette Score')
plt.show()