上周有个做风控的朋友拿着一份用户行为明细来找我,说 KMeans 跑了七八遍,每次聚出来的簇都不一样,问我 K 到底该填几。这个场景我太熟了——KMeans 的算法本体十分钟就能跑通,真正吃掉整个项目的,全在 K 的选取上。这篇就聊透 KMeans 的 K 如何选取这件事:为什么没有唯一正确答案,主流的五六种选 K 方法各自适合什么数据,怎么用一套可复现的流程把 K 定下来,以及那些文档里不写、只有踩过才知道的坑。不管你是刚接触聚类的学生,还是已经在做用户分群、图像压缩、异常检测的一线开发者,都能从这里拿走能直接跑的东西。
先说一个容易被忽略的前提:选 K 本质上不是单纯的数学问题,而是数学指标和业务可解释性之间的谈判。你算出来 K=7 轮廓系数最高,但业务方只能接受把用户分成 3 档去做运营,那最终答案大概率是 3。把这句话记在心里,后面所有的指标都只是你的谈判筹码。
1. 先搞清楚 K 到底在跟谁较劲
1.1 KMeans 的目标函数决定了它永远不会“告诉你”K
KMeans 优化的目标函数是簇内平方误差和,也就是常说的 SSE 或者 inertia:把每个样本到它所属簇中心的距离平方加起来,让这个和尽量小。这个目标函数有个非常讨厌的性质——它对 K 是单调的。K 从 2 增加到 10,SSE 一定下降,哪怕你的数据根本没有簇结构,只要 K 等于样本数,SSE 直接归零。
很多人第一次跑 KMeans 会困惑:为什么 K 越大 inertia 越小,那不是应该选最大的 K 吗?这就是把“拟合程度”和“结构合理性”混为一谈了。KMeans 只负责在你的 K 给定之后,找到一个局部最优的划分,它从不负责告诉你这个 K 合理不合理。换句话说,K 是输入,不是输出,算法本身没有能力自我评价。
所以选 K 的所有方法,本质上都是在做同一件事:给“SSE 随 K 下降”这条单调曲线加一个惩罚项,或者找一个拐点,或者换一个不以 SSE 为尺度的评价标准。理解了这个统一逻辑,你再看肘部法、轮廓系数、Gap 统计量,就不会觉得它们是零散的技巧,而是同一个问题的不同解法。
1.2 一个自问自答:什么样的 K 才算“对”
我在实际项目里会把“对的 K”拆成三个层次的判断,从弱到强:
| 层次 | 判断标准 | 是否必须满足 |
|---|---|---|
| 数学层 | 内部指标(轮廓系数、CH、DB)表现良好,且结果可复现 | 必须 |
| 结构层 | 簇的大小、密度大致均衡,没有明显被切碎的簇 | 强烈建议 |
| 业务层 | 每个簇能被一句话描述清楚,且对应可执行的策略 | 视场景而定 |
只满足数学层的 K,往往在汇报时会被问“所以第七类人是什么人”——答不上来就白做了。只满足业务层的 K,可能数学上是硬凑的,模型后续做预测或者做画像时会露馅。三个层次能同时满足的 K,才是最终答案。
还有一点必须提前说清楚:如果数据本身没有簇结构,那最正确的答案是“不要用 KMeans”,而不是找出一个 K。判断有没有簇结构,可以用 Hopkins 统计量先做一次聚类趋势检验,也可以用后面要讲的稳定性评估来间接判断。这个前置判断能帮你省掉大量无效调参。
2. 六种选 K 方法,各自的地盘在哪
2.1 肘部法:最常用,也最容易被误用
肘部法的思路极其朴素:把 K 从 2 到 10 依次跑一遍,画出 SSE 对 K 的折线,找那条曲线上“从陡变缓”的拐点。它的直觉是——在真正的簇数量之前,每增加一个 K,模型能消除大量误差;超过真实簇数之后,再增加 K,收益急剧衰减,曲线就拐了。
优点不用多说:计算成本极低,可视化一眼就能看懂,汇报时不会有人听不懂。缺点同样明显:拐点经常不明显,曲线是条平滑的凸曲线,你说是 K=3 拐也行,K=4 拐也行,全靠肉眼主观判断。数据稍微脏一点或者簇有重叠,肘部直接消失。
实操里我通常这么用:不画原始的 SSE,而是画 SSE 的下降幅度,或者画归一化后的曲线。下降幅度 K 从 2 到 3 掉了 800,从 3 到 4 掉了 300,从 4 到 5 掉了 260,那 4 就是拐点候选。这种写法比看原始曲线敏感得多。
import numpy as np from sklearn.cluster import KMeans from sklearn.preprocessing import StandardScaler X_scaled = StandardScaler().fit_transform(X) ks = range(2, 11) inertias = [] for k in ks: km = KMeans(n_clusters=k, n_init=20, random_state=42) km.fit(X_scaled) inertias.append(km.inertia_) # 相邻下降量,比看原始曲线更容易定位拐点 drops = [-np.diff(inertias)[i] for i in range(len(inertias) - 1)] for i, d in enumerate(drops): print(f"K={i+2} -> K={i+3} 的 SSE 下降量: {d:.1f}")必须提醒一句:肘部法对特征尺度极度敏感。如果你不做标准化,量纲大的那几维会主导距离计算,SSE 曲线会变成一条几乎看不出结构的线。标准化是肘部法的前置动作,不是可选项。
2.2 轮廓系数:我更愿意把它当作默认首选
轮廓系数衡量的是“簇内紧密度”和“簇间分离度”的比值。对每个样本 i,算它到同簇其他样本的平均距离 a(i),再算它到最近邻簇所有样本的平均距离 b(i),然后 s(i) = (b(i) - a(i)) / max(a(i), b(i))。所有样本的 s 取平均,就是整体轮廓系数,取值范围 [-1, 1],越接近 1 说明簇内越紧、簇间越远。
为什么我更喜欢它:它不像 SSE 那样单调,它是一个“有峰值”的指标,最大值所在的位置天然就是候选 K,不需要人肉找拐点。而且它对簇的形状有一定要求,能间接帮你发现“某个 K 把一个本该完整的簇切成两半”的情况——切碎之后簇间距离变小,轮廓系数会掉下来。
但也有雷区。轮廓系数的计算复杂度是 O(n²),一万条样本还能忍,十万条就基本跑不动了。大数据量下常用的做法是随机采样 5000 到 10000 条算轮廓系数,或者换用 Calinski-Harabasz 指数(CH 指数计算是 O(n) 的,代价低得多)。
from sklearn.metrics import silhouette_score, calinski_harabasz_score, davies_bouldin_score for k in range(2, 11): km = KMeans(n_clusters=k, n_init=20, random_state=42).fit(X_scaled) labels = km.labels_ sil = silhouette_score(X_scaled, labels, sample_size=8000, random_state=42) ch = calinski_harabasz_score(X_scaled, labels) db = davies_bouldin_score(X_scaled, labels) print(f"K={k:2d} silhouette={sil:.4f} CH={ch:.1f} DB={db:.4f}")注意 CH 和 DB 的方向:CH 越大越好,DB 越小越好,轮廓系数越大越好。把三个指标的极值位置列在一起看,如果它们指向同一个 K,这个 K 就很稳;如果指向三个不同的 K,说明数据的簇结构本身模糊,这时候要回到稳定性检验去交叉判断。
2.3 Gap 统计量:给肘部法补上统计学依据
Gap 统计量的思路很聪明:与其主观判断曲线在哪拐,不如拿数据跟“完全没有簇结构的均匀随机数据”比一比。具体做法是在数据的最小覆盖矩形内生成若干组均匀分布的参考数据,对每组参考数据也跑 KMeans,算出 log(SSE),取平均得到 E*[log W_k]。然后
Gap(K) = E*[log W_k] - log W_k
Gap 越大,说明真实数据的聚类结构比随机数据紧凑得越多。选 K 的规则是:取满足 Gap(K) >= Gap(K+1) - s_{K+1} 的最小 K,其中 s 是参考分布 log SSE 的标准差修正项。
这个方法的优势在于它有统计学的假设检验味道,不只是看图说话,理论上是目前最严谨的选 K 方法之一。代价是计算量大——每个 K 都要额外跑几十次参考数据的 KMeans,K 取 2 到 10、B 取 50 的话就是四五百次拟合,中等规模数据也要跑上几分钟。
我的使用经验是:当肘部法和轮廓系数给出了矛盾的答案,或者数据几乎没有明显簇结构时,用 Gap 统计量来做仲裁。日常调参不必每次都用它。
2.4 信息准则类:换个模型来帮你定 K
BIC 和 AIC 本来是给高斯混合模型(GMM)这类概率模型用的,它们衡量的是模型拟合优度与参数数量的权衡。GMM 与 KMeans 的关系很近——KMeans 可以看成是协方差矩阵固定为球形、且各簇权重相等时 GMM 的特例。所以很多人干脆用 GMM 跑一遍 BIC 曲线,取 BIC 最小的 K,再把这个 K 拿去跑 KMeans。
这个做法的合理性在于:BIC 自带参数惩罚项,K 越大参数越多,惩罚越重,天然抑制了“K 越大越好”的倾向。而且 GMM 允许簇有不同的形状和大小,比 KMeans 的球形假设更宽容,在某些椭球状分布的数据上,GMM 给出的 K 反而更接近真实结构。
from sklearn.mixture import GaussianMixture bic_list = [] for k in range(2, 11): gmm = GaussianMixture(n_components=k, covariance_type="full", n_init=5, random_state=42) gmm.fit(X_scaled) bic_list.append(gmm.bic(X_scaled)) print(f"K={k} BIC={bic_list[-1]:.1f}") best_k = 2 + int(np.argmin(bic_list)) print("BIC 推荐的 K:", best_k)要提醒的是,BIC 对样本量敏感,样本量小时容易偏向更小的 K,甚至选出 K=1。另外 GMM 的初始化比 KMeans 更依赖随机种子,n_init 不要设太小,否则 BIC 曲线会抖得很厉害。
2.5 X-means 与二分 KMeans:让算法自己决定要不要分裂
X-means 是 KMeans 的扩展版,从 K=2 开始,对每个簇用 BIC 判断“这个簇该不该再劈成两个”,能劈就劈,直到没有簇值得分裂为止。它的好处是完全自动化,不需要你先给一个 K 的范围,特别适合那种簇数量事先完全没谱的探索性场景。
二分 KMeans 是另一个思路:一开始所有样本在一个大簇里,每次挑一个“分裂后 SSE 下降最多”的簇一分为二,重复到你想要的簇数量。它最大的价值不是自动定 K,而是缓解 KMeans 对初始化的敏感——因为每次只做二分,不会出现全局初始化跑偏的情况。实际项目里我经常用它作为对照组:如果二分 KMeans 在 K=5 时的 SSE 和普通 KMeans 在 K=5 时的差距很大,说明普通 KMeans 陷在局部最优里,你需要加大 n_init。
2.6 稳定性评估:把“随机种子一换就变”变成量化指标
这可能是最被低估的一招。做法是:对原始数据做多次 bootstrap 重采样,每次重采样都跑一遍 KMeans,然后把不同重采样之间的聚类结果做一致性比对,用调整兰德指数(ARI)或者调整互信息(AMI)来衡量。ARI 越接近 1,说明这个 K 下的划分越稳定。
它的逻辑很直观:如果一个 K 是“真”的,那你在数据里抽掉 20% 的样本重新聚,剩下的样本应该还是分在同一堆里;如果 K 是硬凑的,那少几个样本,簇的边界就大幅漂移,ARI 会明显掉下来。
from sklearn.utils import resample from sklearn.metrics import adjusted_rand_score def stability_score(X, k, n_boot=20, frac=0.8, seed=42): rng = np.random.RandomState(seed) base = KMeans(n_clusters=k, n_init=20, random_state=seed).fit_predict(X) scores = [] n = X.shape[0] for _ in range(n_boot): idx = rng.choice(n, int(n * frac), replace=False) sub_labels = KMeans(n_clusters=k, n_init=20, random_state=seed).fit_predict(X[idx]) scores.append(adjusted_rand_score(base[idx], sub_labels)) return float(np.mean(scores)) for k in range(2, 11): print(f"K={k} 稳定性 ARI={stability_score(X_scaled, k):.4f}")这个方法计算量中等,但信息量极大,尤其在业务数据这种簇结构模糊的场景里,它常常是唯一能给出明确结论的指标。我的经验阈值是:ARI 平均在 0.85 以上的 K 才算稳定,0.7 到 0.85 之间属于勉强,低于 0.7 说明这个 K 的划分基本是噪声。
3. 一套能直接抄的选 K 流程
3.1 前置处理:这几步不做,后面全是白算
选 K 之前有三件事必须做完,顺序不能反。
第一是缺失值和异常点。KMeans 用欧氏距离,一个极端离群点足以单独占掉一个簇,导致你以为需要更大的 K。我的做法是先算每个样本到全局中位数的距离,把超过 99.5 分位的样本单独拎出来看,确认是录入错误就删,是真实极端值就单独建一个“异常簇”,不参与主聚类。
第二是标准化。所有涉及距离的算法都对量纲敏感,这是常识,但真正踩坑的地方在于:类别型特征 one-hot 之后,多个 0/1 列加起来的“距离贡献”可能压过连续特征。比如你有 20 个 one-hot 列,每个样本在这 20 维上的欧氏距离平方期望值是 0.5 左右,而一个标准化后的连续特征贡献只有 1 左右。这时候聚类结果几乎完全由类别特征决定,连续特征形同虚设。解决办法是给 one-hot 块整体降权,或者改用 Gower 距离配合层次聚类。
第三是降维的可视化确认。把数据用 PCA 降到两维画散点图,你大概能看出有几个团。这一步纯粹是给你自己一个先验判断,后面指标算出来跟你的目测差太远时,你就知道要去查指标本身的问题了。
3.2 把六种指标塞进一个循环里跑
实际项目里我不会只用一种方法,而是把能算的指标一次性算完,做成一张表横向对比。下面这个脚本可以直接改成你自己的版本:
import pandas as pd from sklearn.cluster import KMeans from sklearn.metrics import (silhouette_score, calinski_harabasz_score, davies_bouldin_score) def evaluate_k(X, k_range=range(2, 11), seed=42): rows = [] for k in k_range: km = KMeans(n_clusters=k, n_init=20, random_state=seed).fit(X) lb = km.labels_ rows.append({ "k": k, "inertia": km.inertia_, "silhouette": silhouette_score(X, lb, sample_size=8000, random_state=seed), "ch": calinski_harabasz_score(X, lb), "db": davies_bouldin_score(X, lb), "min_cluster_ratio": pd.Series(lb).value_counts(normalize=True).min(), }) df = pd.DataFrame(rows) df["inertia_drop"] = -df["inertia"].diff() return df report = evaluate_k(X_scaled) print(report.to_string(index=False))这张表里我特意加了min_cluster_ratio这一列,它是“最小簇占总样本的比例”。这一列极其关键但经常被忽略——如果 K=6 时最小簇只占 0.2%,那这个簇很可能只有几十个样本,它在业务上八成是个噪声集合,不是真正的分群。我一般会设一个下限:最小簇不低于 3% 到 5%,低于这个值就把这个 K 排掉。
3.3 怎么读这张表:三步交叉裁决
拿到表之后不要直接取 silhouette 的最大值,按下面的顺序走:
第一步,先看稳定性。如果 ARI 那一列(需要另外算)在某个 K 之后断崖式下跌,那个断崖点之前的 K 才是可信区间。这一步先砍掉一半的候选。
第二步,在可信区间里比轮廓系数和 CH。两个指标都指向同一个 K,基本就定了。如果错开一位,看肘部的下降量曲线更支持哪一个。
第三步,看最小簇比例和业务可解释性。跑出各簇的中心点,用原始尺度还原,逐个描述每个簇的特征。描述不出来的簇,说明它没有业务含义,要么合并,要么索性少一个 K。
这套流程我在几个不同类型的数据上走过——用户行为、设备传感器、文本嵌入向量,结论都比较一致:真正稳定的 K 通常比轮廓系数的最大值所在位置略小一到两个。轮廓系数偏爱把簇切得更细,因为它只看紧密度和分离度的比值,而稳定性看的是结构的一致性,后者更接近“真实结构”。
3.4 关于 n_init 和 random_state 的两个坑
n_init参数在 sklearn 的旧版本默认是 10,新版本改成了"auto"。如果你用的是旧版本且没注意,跑出来的 K 曲线会在不同次运行之间跳动,你以为是数据本身的问题,其实是初始化次数太少。
我现在的习惯是统一写成n_init=20,并且固定random_state。固定种子不是为了“作弊”,而是为了让不同 K 之间的对比建立在同一初始条件附近,减少随机噪声对指标比较的干扰。真要评估稳定性,恰恰要反过来——换不同的 seed 跑十次,看聚类结果 ARI 的平均值和方差,这才是稳定性的正确测法。
第二个坑是max_iter。默认 300 次迭代在大数据上经常不够,尤其 K 比较大、簇的形状比较扁的时候,你会看到 warning 说未收敛。未收敛的 inertia 是不可比的,因为不同 K 的迭代停止位置不一样。看到这个 warning,第一件事是把max_iter提到 1000 再看指标,不要直接拿结果去解读。
4. 排查实录:那些让你怀疑人生的场景
4.1 肘部法完全没有拐点
这是最常见的求助场景。曲线是一条平滑的下凸线,斜率缓慢衰减,看不出任何转折。可能的原因有三个:
一是数据本身没有簇结构,是均匀分布或者单一方向的连续变化。这种情况下任何 K 都是硬切,解决办法是换模型,比如用密度聚类或者干脆不用聚类,改用分位数分箱做分层。
二是簇的密度差异极大,一个大而松的簇包着几个小而紧的簇。KMeans 的 SSE 被大簇主导,小簇的贡献完全被淹没。可以试着对样本做加权,或者先对数据进行一次对数变换压缩量纲差异。
三是特征维度太高,距离集中在很窄的范围内(这就是维度灾难里的“距离集中”现象),所有点看起来都差不多远。降维到 10 到 20 维再试,或者改用余弦距离配合球形 KMeans。
4.2 轮廓系数最大值出现在 K=2,但业务要求至少五档
这个问题我遇到过不止一次。K=2 时轮廓系数最高,因为把数据一分为二永远是最容易获得高分离度的划分,这是轮廓系数的固有偏好——它在 K 较小的时候天然占优。
处理办法不是强行否定指标,而是换一个评价口径。具体做法是:只在 K 大于等于业务下限的范围里比较指标。比如业务要求至少 5 档,那就只在 K 从 5 到 12 的区间里找轮廓系数和 CH 的极值,把 K=2 到 4 视为“业务不可行区”直接排除。这在逻辑上是站得住的——选 K 本来就是在可行域内寻优,不是在无约束空间里寻优。
同时也要反思一下:如果 K=2 的轮廓系数远高于 K=5,而你又必须做 5 档,那说明这 5 档里至少有 3 档是人为切出来的。这时候可以考虑分层聚类——先用 K=2 做粗分,再在每个大簇内部单独做一次小 K 的聚类,形成两级结构。这种分层方案在用户分层运营里非常常见,比硬凑五档效果好得多。
4.3 高维稀疏数据:所有指标都失灵
文本向量、one-hot 编码的交易特征,这类数据的维度动辄几百上千,KMeans 的很多假设都不成立了。现象是:不管 K 取几,轮廓系数都在 0.1 上下,CH 指数跳动剧烈,不同 seed 跑出来的 ARI 只有 0.4 左右。
这种数据要先降维再聚类。常用的路径是先做 TruncatedSVD(稀疏矩阵友好)降到 50 到 100 维,或者用自编码器压到低维嵌入,再跑 KMeans。降维之后再算指标,通常会看到一个清晰的峰值。另外这类数据上,球形 KMeans 配合余弦归一化往往比标准 KMeans 效果好,因为方向信息比长度信息更有意义。
4.4 各簇样本量极度不均衡
KMeans 有个隐含假设:各簇大小差不多、形状接近球形。当真实结构是一个大簇加几个很小的簇时,KMeans 会把小簇硬并进大簇,或者用多个簇去拼凑大簇的形状,导致指标显示需要比真实更多的 K。
识别这个问题的信号很明确:inertia 曲线在某个 K 之后仍然大幅下降,但轮廓系数不升反降。这说明新增的簇不是在切分真实结构,而是在给大簇“打补丁”。
5. 一张速查表收尾
5.1 问题与对策对照
| 现象 | 最可能的原因 | 优先尝试的解法 |
|---|---|---|
| 肘部曲线无拐点 | 无簇结构 / 密度差异大 / 维度灾难 | 先做 Hopkins 检验,再降维或换模型 |
| 轮廓系数峰值在 K=2 | 该指标对 K 小的天然偏好 | 限定业务可行区间后再比较 |
| 换 seed 结果大变 | n_init 太小 / 簇结构模糊 | n_init 提到 20 以上,看稳定性 ARI |
| 最小簇占比过低 | 异常点吸走了一个簇 | 先剔除离群点再跑 |
| 指标之间互相矛盾 | 数据簇结构本身弱 | 以稳定性指标为准,回归业务判断 |
| 高维数据指标普遍偏低 | 距离集中效应 | 先降维到 50 维以内再聚类 |
5.2 我个人的默认起手式
如果不想一上来就折腾六种方法,我给一个最省事的默认路径:标准化 → 算 K 从 2 到 10 的 inertia 和 silhouette → 算稳定性 ARI → 三个结果的交集就是答案。这个流程跑一遍大概三到五分钟,能覆盖八成以上的常见场景。剩下的两成,通常是有特殊结构或者特殊业务约束,才需要动用到 Gap 统计量、GMM 的 BIC 或者分层聚类。
最后分享一个我踩过很多次才悟出来的点:选 K 的过程中,最贵的成本不是算力,是反复推翻重来。所以我现在会先花十分钟把候选 K 范围内每个簇的中心点用原始尺度打出来,逐个读一遍,看能不能用一句话概括。能概括的 K,哪怕指标不是最优,我也会优先选它——因为后面所有基于这个分群的分析、策略、报表,都要建立在这些簇能被解释的基础上。指标决定候选,可解释性决定最终答案。