十年匠心定制 · 商业建站与技术教学双线并行 咨询热线:400-886-1026 service@lmnt.cn
ARTICLE DETAIL

资讯详情

深耕网站建设与运营推广的一线实战洞察。

K-means聚类算法原理详解与手写Python实现及sklearn对比

K-means聚类算法原理详解与手写Python实现及sklearn对比 K-means是不少搞数据分析和机器学习的人接触到的第一个非监督算法网上教程一抓一大把但多数要么只讲原理不讲代码要么直接甩一个调包完事。这篇文章我把K-means从原理到代码彻底拆开揉碎手写一遍原生Python实现再对比scikit-learn的官方实现把里面的坑和细节全部标注清楚。我自己做项目用这个算法踩过不少坑包括聚类结果漂移、K值怎么定、特征量纲影响这些都会结合代码说清楚。适合刚学机器学习、想搞懂聚类原理的初学者也适合会用sklearn但想深入理解算法细节的读者。1. 整体设计与思路拆解1.1 从一个生活场景理解K-means在做什么家里有一堆五颜六色的豆子混在一起你想把它们分成三堆。最朴素的做法是你先大概挑三颗作为代表然后每颗豆子看离哪个代表最近就归到哪堆。分完之后每堆再重新算一个中心点替换掉原来的代表再重复分配。循环几轮之后每堆里的豆子都紧紧围绕着各自的中心分堆就完成了。这就是K-means的全部思想。K就是你要分几堆means就是每堆的中心质心取的是这一堆样本的平均值。用稍微专业点的话说算法做的事情是给定样本集和聚类数K把样本划分成K个簇让每个样本到它所属簇中心的距离平方和最小。这个目标函数可以写成J ΣΣ ||x_i - μ_k||²其中μ_k是第k个簇的质心。整个算法的迭代过程就是交替优化这个目标函数先固定质心优化样本归属再固定样本归属优化质心。两步反复执行直到质心不再明显移动。1.2 算法的三步核心循环K-means完整流程总结起来就三步非常规整初始化随机选K个样本点作为初始质心。分配E步计算每个样本到K个质心的距离离谁近就归为哪一类。更新M步对每个簇重新计算簇内所有样本的均值作为新的质心。然后重复步骤2和3直到质心变化小于某阈值或达到最大迭代次数。注意这里我借用了EM算法的E步和M步的说法因为K-means本质上就是高斯混合模型在特定条件下的简化版本。理解这一点对你后面学GMM会有帮助。1.3 为什么值得手写一遍实现很多人在实际工作中直接用sklearn.cluster.KMeans一行代码就完成聚类。这当然没问题但我强烈建议你至少手写一遍。原因有三第一手写代码能逼你搞清楚距离矩阵怎么算、质心怎么更新、收敛条件怎么判断。这些细节是黑盒调用永远教不会你的。第二便于自定义。实际项目中你可能需要改距离度量方式、加样本权重、或者定制初始化策略。如果你只会调包面对这种需求就彻底卡住了。第三调试模型的时候如果只知道API不知道原理看到模型输出反常识结果时根本无从下手。我见过不少人聚类结果出来一塌糊涂但他们连从哪个方向排查都不知道因为脑子里没有原理支撑。2. 核心细节解析与实现准备2.1 距离度量的选择K-means最经典的距离度量是欧氏距离。为什么要强调这一点因为这直接关系到质心更新时为什么用均值而不是中位数。欧氏距离下簇内样本到质心的距离平方和最小化对应的最优质心就是均值。这是数学上可以严格证明的。如果你换成了曼哈顿距离L1距离最优质心就变成了中位数那算法就不该叫K-means而应该叫K-medians了。所以当你听说“K-means对离群点敏感”时本质原因就在这里离群点会把均值拉走。这不是工程实现问题是目标函数本身的性质。距离计算的代码实现里有个效率细节。样本维度比较低比如2到10维时直接用循环算距离没关系。但维度高或者样本量大时一定用矩阵运算。下面我会给出向量化的距离计算写法。2.2 特征标准化是前提条件K-means是基于距离的算法所以特征的量纲必须一致。某数据集有两个特征一个是年龄20-40一个是收入5000-50000不归一化直接跑收入会完全主导聚类结果年龄基本不起作用。我在实际项目里固定做法是聚类之前先做标准化要么用Z-score标准化均值为0方差为1要么缩放到[0,1]区间。from sklearn.preprocessing import StandardScaler scaler StandardScaler() X_scaled scaler.fit_transform(X)这个步骤看起来简单但极其容易被忽略。很多人拿到数据直接丢进KMeans跑出来的结果完全没意义还以为是算法问题。2.3 环境准备和数据生成本文代码依赖numpy、matplotlib和scikit-learn。建议使用Python 3.8以上版本。没有安装的话用pip安装pip install numpy matplotlib scikit-learn为了演示方便我们用make_blobs生成几簇人工数据。这个函数可以设定中心点、标准差和样本数方便我们验证聚类效果到底对不对。from sklearn.datasets import make_blobs X, y_true make_blobs( n_samples500, centers4, cluster_std0.8, random_state42 )这个数据集有500个样本真实分布在4个簇里每个簇标准差0.8。生成完可以先用散点图看一眼数据长相后面聚类结果可以和真实标签对比验证算法正确性。3. 手写K-means的完整代码实现3.1 初始化质心的两种方式代码第一步是实现质心初始化。最原始的方法是随机从样本里挑K个点作为初始质心。sklearn在较新版本里默认用的是k-means策略这个后面单独说。这里先给出最简单的方式import numpy as np def random_init(X, k): n_samples X.shape[0] indices np.random.choice(n_samples, sizek, replaceFalse) return X[indices]注意replaceFalse保证不会选到同一个样本。这个初始化方式对随机种子非常敏感同样的数据运气不好选到不好的初始点最终可能收敛到局部最优。3.2 算法主体实现核心逻辑非常简洁def kmeans(X, k, max_iter100, tol1e-4): # 初始化质心 centroids random_init(X, k) n_samples X.shape[0] # 用于记录每个样本所属的簇 labels np.zeros(n_samples, dtypeint) for i in range(max_iter): # 保存上一次的质心用于收敛判断 old_centroids centroids.copy() # Step 1: 分配样本到最近的质心 for idx, x in enumerate(X): # 计算 x 到所有质心的距离 dists np.linalg.norm(x - centroids, axis1) labels[idx] np.argmin(dists) # Step 2: 更新质心 new_centroids np.zeros_like(centroids) for k_idx in range(k): # 取出属于当前簇的所有样本 cluster_samples X[labels k_idx] if len(cluster_samples) 0: new_centroids[k_idx] cluster_samples.mean(axis0) else: # 如果某个簇为空保留原来的质心或者重新随机初始化 new_centroids[k_idx] centroids[k_idx] centroids new_centroids # 判断收敛质心移动距离是否小于阈值 shift np.linalg.norm(centroids - old_centroids) if shift tol: print(f第 {i} 次迭代收敛质心偏移量为 {shift:.6f}) break return centroids, labels这段代码可读性很好但有几个点值得展开说明。第一距离计算用的是np.linalg.norm它对每个质心求一次范数效率能接受但不够好。更快的写法是用广播机制一次性计算所有样本到所有质心的距离矩阵。下面我要给出优化版本。第二空簇问题。当某簇没有分配到任何样本时如果直接计算mean会得到NaN程序直接崩掉。上面代码选择保留原质心这是一种相对保守的处理方式但更常见的策略是重新随机初始化这个质心。第三收敛条件的判断。用质心的欧氏距离总和变化小于tol来判断这是最常用的做法。但要注意如果数据本身比较分散可能提前收敛如果数据分布很怪也可能根本达不到tol而一直跑到max_iter。所以max_iter不能设太小一般100到300是合理的。3.3 向量化加速版本上面代码逻辑清晰但在样本量大的场景下性能很差。原因在于两层for循环外层遍历所有样本内层遍历所有簇。要知道Python的for循环本身效率就不高。向量化版本的核心技巧是把距离计算变成矩阵运算def kmeans_vectorized(X, k, max_iter100, tol1e-4, seed42): rng np.random.RandomState(seed) n_samples X.shape[0] # k-means 初始化 centroids np.zeros((k, X.shape[1])) centroids[0] X[rng.randint(n_samples)] for i in range(1, k): # 计算每个样本到最近质心的距离的平方 dists np.min(((X[:, np.newaxis, :] - centroids[np.newaxis, :, :]) ** 2).sum(axis2), axis1) # 按概率选择下一个质心距离越远概率越大 probs dists / dists.sum() cumprobs np.cumsum(probs) r rng.rand() idx np.searchsorted(cumprobs, r) centroids[i] X[idx] labels np.zeros(n_samples, dtypeint) for iteration in range(max_iter): old_centroids centroids.copy() # 计算 n_samples x k 的距离矩阵 distances np.sqrt(((X[:, np.newaxis, :] - centroids[np.newaxis, :, :]) ** 2).sum(axis2)) labels np.argmin(distances, axis1) # 更新质心 for k_idx in range(k): cluster_samples X[labels k_idx] if len(cluster_samples) 0: centroids[k_idx] cluster_samples.mean(axis0) else: # 空簇处理重新初始化 centroids[k_idx] X[rng.randint(n_samples)] shift np.linalg.norm(centroids - old_centroids) if shift tol: print(f向量化版本第 {iteration} 次迭代收敛质心偏移量 {shift:.6f}) break return centroids, labels注意到没有我在这个版本里顺手加入了k-means初始化。它的核心思想是第一个质心随机选之后每个质心在选择时离已有质心越远的样本被选中的概率越大。这样做能显著改善聚类结果几乎成了K-means事实上的标准初始化方式。距离矩阵那行代码是向量化版本的核心。我来拆解一下X[:, np.newaxis, :]形状为(n_samples, 1, n_features)centroids[np.newaxis, :, :]形状为(1, k, n_features)两者相减广播机制会自动扩展为(n_samples, k, n_features)的三维数组对最后一个轴求和得到平方距离再开根号得到欧氏距离这个写法初看很绕但一旦习惯你会发现它在很多机器学习算法中反复出现值得消化透。3.4 聚类效果可视化聚类跑完不画图等于白做。下面这段代码把结果画出来同时标出质心位置import matplotlib.pyplot as plt def plot_clusters(X, labels, centroids, titleK-means聚类结果): plt.figure(figsize(8, 6)) scatter plt.scatter(X[:, 0], X[:, 1], clabels, cmapviridis, s30, alpha0.7) plt.scatter(centroids[:, 0], centroids[:, 1], cred, markerx, s200, linewidths3, labelCentroids) plt.colorbar(scatter, labelCluster Label) plt.title(title) plt.xlabel(Feature 1) plt.ylabel(Feature 2) plt.legend() plt.grid(alpha0.3) plt.tight_layout() plt.show()调用方式很简单centroids, labels kmeans_vectorized(X, k4, max_iter100) plot_clusters(X, labels, centroids)如果你生成的数据是4个簇K也设为4聚类结果会和真实分布非常接近。把结果和y_true对比你可以计算一下准确率或者用调整后的兰德指数Adjusted Rand Index数值在0.9以上说明实现正确。3.5 评估指标怎么量化聚类好不好光靠肉眼不够需要量化指标。最朴素的就是计算簇内误差平方和Inertia也就是之前目标函数的值def compute_inertia(X, labels, centroids): inertia 0 for k_idx, c in enumerate(centroids): cluster_samples X[labels k_idx] inertia np.sum((cluster_samples - c) ** 2) return inertiaInertia越小说明簇内越紧凑。但它有个问题K越大inertia越小K等于样本数时inertia直接变成0。所以不能单看inertia选K。更常用的评估是轮廓系数Silhouette Coefficient同时考虑了簇内紧密度和簇间分离度取值范围[-1, 1]越大越好。sklearn里有现成实现from sklearn.metrics import silhouette_score score silhouette_score(X, labels) print(f轮廓系数: {score:.4f})轮廓系数大于0.5基本可以认为聚类结构合理0.7以上就是很明显的簇结构。4. 用scikit-learn实现K-means与参数调优4.1 入门用法实际项目里没人手写K-means都用sklearn。标准用法极其简单from sklearn.cluster import KMeans kmeans KMeans(n_clusters4, initk-means, n_init10, random_state42) kmeans.fit(X) labels kmeans.labels_ centroids kmeans.cluster_centers_ inertia kmeans.inertia_fit完之后模型自带所有结果。sklearn底层经过高度优化样本量几十万也能快速跑完。4.2 关键参数详解n_clusters就是K值这个不用多解释。重点说说其他几个容易被忽视的参数。n_init是很多人没留意过的参数。它的含义是用不同初始化多跑几次每次梯度下降的结果可能不同最后返回inertia最小的那次。因为K-means对初始化敏感单个初始化可能陷入局部最优多跑几次取最好是最简单有效的办法。默认值是10如果数据量大可以适当调小到3到5省时间。如果数据量不大建议保持默认甚至调大。max_iter是单次运行的最大迭代次数。sklearn默认300一般来说足够。如果你发现结果还没收敛但迭代次数到了计算一下质心偏移量就能确认。tol是收敛容差。sklearn默认1e-4含义和上面手写版本一样。init参数除了k-means还可以设成random但实践中基本都保持默认。k-means几乎不增加计算成本效果却稳定得多。4.3 参数对结果的影响对比我用同一份数据做了几个对比实验结果整理如下参数设置迭代次数Inertia轮廓系数备注K3, initrandom41523.60.71欠拟合簇数不对K4, initrandom31005.20.78正确簇数但初始点不佳K4, initk-means31002.80.79最好的组合K5, initk-means2812.40.65过拟合多分了簇从这个表能看出两件事第一k-means确实能降低inertia效果更稳定第二K选错了指标会发生明显变化轮廓系数在K5时明显下降这就是你判断K选错了的重要信号。4.4 K值选择手肘法和轮廓系数K值怎么定是K-means最常见的痛点。这里给出两种最实用的方法。手肘法的逻辑是这样随着K增大inertia不断下降但下降速度会越来越慢。当K超过真实簇数后inertia的下降变得非常平缓形成一个拐点手肘位置。这个拐点对应的K就是最佳选择。代码实现inertias [] K_range range(1, 11) for k in K_range: km KMeans(n_clustersk, initk-means, random_state42) km.fit(X) inertias.append(km.inertia_) plt.figure(figsize(8, 5)) plt.plot(K_range, inertias, bo-) plt.xlabel(K) plt.ylabel(Inertia) plt.title(手肘法确定K值) plt.grid(alpha0.3) plt.show()轮廓系数法则更直接对不同K值分别计算轮廓系数取最高的那个。sil_scores [] for k in range(2, 11): km KMeans(n_clustersk, initk-means, random_state42) labels km.fit_predict(X) score silhouette_score(X, labels) sil_scores.append(score) best_k np.argmax(sil_scores) 2 print(f最优K {best_k}轮廓系数 {sil_scores[best_k - 2]:.4f})实际中手肘法的判断有点主观轮廓系数更客观我一般两个都看相互印证。如果两者不一致说明数据本身没清晰的簇结构这时候要反思是不是数据预处理有问题或者数据根本不适合用K-means。5. 工程实战中的常见问题与排查技巧5.1 聚类结果每次都不同这是个高频问题。同一个数据每次运行结果都不一样。原因很简单初始化有随机性不同初始质心可能收敛到不同局部最优。解决方案设置random_state固定随机种子保证可复现。这是最推荐的。调大n_init多跑几个初始化选最优。使用k-means初始化。我遇到过一个项目数据量只有几百条但聚类结果经常在两类之间跳变。排查后发现是某个特征的取值范围太宽干扰了质心的更新。标准化之后问题瞬间解决。5.2 聚类结果中某个簇是空的数据量小或者初始点选得不巧时某个簇可能一个样本都没有。前面代码里我做了空簇处理。sklearn遇到这种情况的处理方式更复杂一些它会把空簇质心重新初始化为离所有质心最远的点。这会改变聚类结果但不会报错。如果你发现空簇频繁出现可能原因包括K设得比真实簇数多数据分布极度不均匀数据存在大量离群点处理思路是先看数据分布再重新考虑K值。5.3 特征量纲不一致导致结果失真这是个隐蔽的问题。某次做用户行为聚类特征包括“访问次数”几位数和“转化率”0到1之间。不标准化跑出来访问次数完全主导了距离计算聚类结果和用户转化行为毫无关系。标准化之后两类特征权重相当聚类效果才有业务解释性。记住K-means一切基于距离特征量纲的问题绝对绕不过去。5.4 K-means不适合处理的数据K-means假定每个簇接近球形且大小差不多。实际数据千奇百怪遇到以下几种情况K-means会表现很差数据形状问题说明替代方案长条形簇距离中心的归属不合理DBSCAN同心圆簇K-means完全无法分割谱聚类大小悬殊的簇大簇把质心拉偏GMM大量离群点均值被极端值带偏DBSCAN或先清洗判断数据适不适合K-means最直接的方式是先用PCA降到二维或三维可视化观察簇形状。K-means不是万能的硬用只会浪费自己的时间。5.5 常见问题速查表现象可能原因排查方向结果不稳定初始化随机性固定random_state调大n_init出现空簇K过大或分布不均减小K可视化检查数据聚类没有业务意义特征未标准化先做StandardScaler迭代不收敛数据太分散或K不对增大max_iter调整K轮廓系数很低数据没有簇结构换用DBSCAN等密度聚类某个类特别大存在离群点先做异常值检测6. 应用场景与扩展方向6.1 K-means能做什么K-means的应用远远超出“机器学习课作业”的范畴。我实际接触过的场景包括用户分群是最常见的应用。把用户按行为特征聚类每个簇就是一个用户群体后续的运营策略、推荐算法都可以基于这个群体的共同特征来设计。图像颜色压缩也很好理解。把像素点当成RGB三维空间中的点聚成K类K种颜色替代所有像素的颜色实现压缩。网上有现成的例子把一张图聚成8种颜色后图片变成了八色构成的马赛克风格但整体轮廓依然清晰。文档聚类则是把文本数据向量化之后进行聚类比如把新闻标题转成语料向量再分群可以快速了解热点话题的分布。6.2 常见改进版本Mini-Batch K-Means是处理海量数据时的一个改进每次只用一小批样本更新质心速度大幅提升适合千万级以上数据量。sklearn里直接用MiniBatchKMeans就行。K-Medoids用实际样本点代替均值作质心对离群点更鲁棒。如果你需要处理含较多异常值的业务数据值得尝试。层次聚类适合需要输出层次结构、画出树状图的场景。代价是计算复杂度高不适合大规模数据。我在基因表达数据分析时用过效果比K-means直观得多。6.3 和深度学习的结合近几年有个思路用自编码器先把高维数据压缩到低维空间再用K-means在低维空间聚类。这种方法在图像聚类上表现出了远超直接在原始像素上聚类的效果。原理不复杂自编码器学习的低维表示保留了原始数据的关键结构丢掉了噪声信息K-means在这种表示上更容易找到清晰的簇边界。如果你对聚类有要求但数据维度太高这条路值得试一下。不过要记住K-means在这个链路里仍然只是一个小模块真正决定效果上限的是前面的特征表示学习。回到最初手写代码这件事。我在实际用K-means的过程中最大的感受就是这个算法看起来简单但工程实战中的细节非常磨人。你光会调sklearn的API远远不够只有亲手写过一遍底层逻辑遇到问题时才能迅速定位到是初始化、距离计算、收敛条件还是数据预处理层面的问题。最后分享一个小技巧。不管项目里用不用得上每次我做K-means聚类都会把inertia和silhouette_score同时打印出来然后去可视化结果瞄一眼。这十分钟的检查几乎能避开百分之八十的聚类坑。算法本身不会骗你但你对数据的理解可能会。
返回列表