# 用Python实战K-Means聚类:从Iris数据集到可视化分析(附完整代码)
如果你刚开始接触机器学习,可能会觉得那些复杂的算法和数学公式有点让人望而生畏。但别担心,今天我们要聊的K-Means聚类算法,可以说是机器学习世界里最“接地气”的算法之一。它不需要你事先给数据打上任何标签,就能自动发现数据中隐藏的自然分组,这个过程本身就充满了探索的乐趣。想象一下,你手头有一堆客户数据,不知道他们有什么共同特征,K-Means能帮你自动把他们分成几个有意义的群体,比如“高价值客户”、“价格敏感型客户”等等,这比一个个手动分析要高效得多。
我刚开始做数据分析的时候,第一次用K-Means分析用户行为数据,那种看着散乱的数据点自动聚集成几个清晰簇群的感觉,至今记忆犹新。它不像深度学习那样需要庞大的算力和数据,几行Python代码,一个经典的数据集,就能让你直观地理解聚类的核心思想。今天,我们就以经典的鸢尾花(Iris)数据集为舞台,手把手带你走完K-Means聚类的完整流程:从数据加载、预处理,到模型训练、结果评估,最后用生动的可视化把抽象的结果呈现出来。我会分享一些实际编码中容易踩的“坑”和解决技巧,确保你读完就能在自己的项目里用起来。
## 1. 理解K-Means:核心思想与算法流程
K-Means算法的目标其实非常直观:把一堆数据点分成K个组,让同一个组内的数据点彼此尽量相似,不同组之间的数据点尽量不同。这里的“相似”通常用距离来衡量,最常用的就是欧几里得距离。你可以把它想象成在教室里给一群学生分组,目标是让每个小组内部的学生兴趣、能力相近,而不同小组之间则有明显的差异。
算法的名字“K-Means”就揭示了它的两个关键点:“K”代表你希望得到的簇的数量,“Means”代表每个簇的中心(即质心)是通过计算簇内所有点的平均值得到的。这个质心就像一个引力中心,不断吸引着周围相似的数据点。
> **注意**:K-Means是一种“硬聚类”算法,这意味着每个数据点有且只能属于一个簇。这与一些“软聚类”算法(如高斯混合模型)不同,后者允许数据点以一定的概率属于多个簇。
整个算法是一个迭代优化的过程,可以概括为以下几个步骤:
1. **初始化**:随机选择K个数据点作为初始的簇质心。
2. **分配阶段**:对于数据集中的每一个点,计算它与K个质心的距离,并将其分配给距离最近的那个质心所在的簇。
3. **更新阶段**:对于每一个新形成的簇,重新计算它的质心(即该簇所有点的平均值)。
4. **迭代**:重复步骤2和步骤3,直到满足停止条件。停止条件通常是质心的位置不再发生显著变化,或者达到了预设的最大迭代次数。
下面这个简单的伪代码可以帮助你理解这个循环:
```python
# 伪代码:K-Means核心流程
初始化 K 个质心
while 质心变化大于阈值或未达到最大迭代次数:
for 数据集中每个点:
计算该点到所有质心的距离
将该点分配给距离最近的质心对应的簇
for 每个簇:
将该簇的质心更新为簇内所有点的均值
```
这个迭代过程本质上是在优化一个目标函数,即**簇内误差平方和**。这个值越小,说明簇内点越紧密,聚类效果越好。算法通过不断调整质心位置来最小化这个值。
**K-Means的优缺点一览**
在实际使用前,了解它的长处和短板至关重要。我整理了一个简单的对比表格,方便你快速把握:
| 特性 | 优点 | 缺点与挑战 |
| :--- | :--- | :--- |
| **原理与实现** | 原理直观,实现简单,收敛速度通常较快。 | 需要预先指定簇数K,K值选择不当会严重影响结果。 |
| **结果质量** | 对于球形、大小密度相近的簇,效果很好,可解释性强。 | 对非球形簇、大小密度差异大的簇效果不佳;对异常值敏感。 |
| **计算性能** | 计算复杂度相对较低,可扩展到大样本量。 | 初始质心随机选择可能导致结果不稳定,容易陷入局部最优。 |
| **数据要求** | 适用于数值型数据,对特征尺度差异敏感(需标准化)。 | 难以处理分类变量,高维数据下距离度量可能失效(“维度诅咒”)。 |
了解了这些,你就能明白为什么我们常说“没有免费的午餐”。K-Means不是万能的,但在符合其假设的场景下,它是一个极其强大且高效的工具。接下来,我们就进入实战环节,看看如何用Python和`scikit-learn`库来驾驭它。
## 2. 环境搭建与Iris数据集初探
工欲善其事,必先利其器。我们首先需要准备好Python环境。我强烈建议使用Anaconda来管理你的Python环境和包,它能省去很多依赖冲突的麻烦。如果你还没有安装,可以去Anaconda官网下载安装。接下来,我们通过终端或Anaconda Prompt创建并激活一个专门用于本教程的虚拟环境。
```bash
# 创建一个名为ml_cluster的新环境,指定Python版本
conda create -n ml_cluster python=3.9
# 激活该环境
conda activate ml_cluster
```
环境激活后,我们来安装必要的库。核心是`scikit-learn`,它提供了K-Means的实现以及我们要用的数据集。`pandas`和`numpy`用于数据处理,`matplotlib`和`seaborn`用于可视化。
```bash
# 安装核心库
pip install scikit-learn pandas numpy matplotlib seaborn
```
安装完成后,就可以在Python脚本或Jupyter Notebook中导入它们了。我习惯在开头一次性导入所有需要的库:
```python
# 导入所需库
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import seaborn as sns
# 从sklearn中导入数据集、聚类算法和评估工具
from sklearn import datasets
from sklearn.cluster import KMeans
from sklearn.preprocessing import StandardScaler
from sklearn.metrics import silhouette_score, silhouette_samples
# 设置绘图风格,让图表更好看
sns.set_style("whitegrid")
plt.rcParams['font.sans-serif'] = ['SimHei'] # 用来正常显示中文标签
plt.rcParams['axes.unicode_minus'] = False # 用来正常显示负号
```
现在,让我们请出今天的主角——**鸢尾花(Iris)数据集**。这个数据集在机器学习领域就像“Hello World”一样经典。它包含了150朵鸢尾花的测量数据,每朵花有4个特征:花萼长度、花萼宽度、花瓣长度、花瓣宽度。这些花实际分属3个品种:山鸢尾、变色鸢尾和维吉尼亚鸢尾。在无监督学习中,我们暂时“忘记”这些品种标签,看看K-Means能否仅凭测量数据重新发现这些自然类别。
用`sklearn`加载这个数据集非常简单:
```python
# 加载Iris数据集
iris = datasets.load_iris()
# 将数据和特征名转换为DataFrame,方便查看
iris_df = pd.DataFrame(iris.data, columns=iris.feature_names)
# 添加真实品种标签列(仅用于后续对比验证,聚类本身不用)
iris_df['species'] = iris.target
# 查看前5行数据
print(iris_df.head())
```
运行这段代码,你会看到一个清晰的表格。是时候先直观地感受一下数据了。我们可以先看看特征的统计摘要,并用散点图矩阵观察两两特征之间的关系。
```python
# 查看数据的基本统计信息
print(iris_df.describe())
# 绘制特征间的散点图矩阵,用颜色区分真实品种
sns.pairplot(iris_df, hue='species', palette='Set2', diag_kind='kde')
plt.suptitle('鸢尾花数据集特征散点图矩阵(按真实品种着色)', y=1.02)
plt.show()
```
观察这些图表,你会发现花瓣长度和花瓣宽度这两个特征似乎能很好地区分三个品种。这给了我们一个初步的信心:数据本身存在可分组的结构。在开始聚类前,还有一个至关重要的步骤——**数据预处理**。K-Means基于距离计算,如果特征的量纲(单位)和尺度差异很大(比如一个特征范围是0-1,另一个是1000-10000),那么数值大的特征会主导距离计算,导致聚类结果失真。因此,我们通常需要对数据进行标准化,使每个特征均值为0,方差为1。
```python
# 分离特征和标签(聚类时我们只使用特征)
X = iris.data
y = iris.target # 真实标签,仅用于评估
# 数据标准化
scaler = StandardScaler()
X_scaled = scaler.fit_transform(X)
print("标准化后的前5行数据:")
print(X_scaled[:5])
```
准备工作全部就绪,数据已经清洗并标准化,接下来就是最激动人心的部分:构建和训练我们的第一个K-Means模型。
## 3. 实战:用scikit-learn实现K-Means聚类
现在,我们将正式使用`sklearn.cluster.KMeans`类来创建模型。这里我们面临第一个关键决策:**K值应该设为多少?** 对于Iris数据集,我们知道真实类别数是3,但在实际项目中,这个数字往往是未知的。我们先假设知道K=3,完成一次完整的聚类,然后再探讨如何寻找最佳K值。
```python
# 实例化KMeans模型,设置n_clusters=3,random_state确保结果可复现
kmeans = KMeans(n_clusters=3, random_state=42)
# 在标准化后的数据上拟合模型
kmeans.fit(X_scaled)
# 获取聚类标签(每个点被分配到的簇编号)
cluster_labels = kmeans.labels_
# 获取最终的3个质心坐标
centroids = kmeans.cluster_centers_
print("聚类标签(前10个):", cluster_labels[:10])
print("质心坐标:\n", centroids)
```
模型训练完成后,我们可以把聚类标签添加回原来的DataFrame,方便对比分析。
```python
# 将聚类结果添加到DataFrame中
iris_df['cluster'] = cluster_labels
# 查看添加了聚类结果的前几行
print(iris_df[['sepal length (cm)', 'sepal width (cm)', 'cluster', 'species']].head(10))
```
仅仅看数字还不够直观,让我们把聚类结果画出来。由于数据有4个维度,我们无法在二维平面完全展示,但可以选取两个特征(比如花瓣长度和宽度)来绘制散点图。
```python
# 使用花瓣长度和花瓣宽度两个特征进行可视化
plt.figure(figsize=(10, 6))
# 绘制数据点,按聚类标签着色
scatter = plt.scatter(X_scaled[:, 2], X_scaled[:, 3], c=cluster_labels, cmap='viridis', s=50, alpha=0.7, edgecolor='k')
# 绘制质心
plt.scatter(centroids[:, 2], centroids[:, 3], c='red', marker='X', s=200, label='质心', edgecolor='k', linewidth=1.5)
plt.xlabel('花瓣长度(标准化后)')
plt.ylabel('花瓣宽度(标准化后)')
plt.title('K-Means聚类结果(K=3)')
plt.legend()
plt.colorbar(scatter, label='簇标签')
plt.show()
```
这张图能让你清晰地看到三个簇的分布以及质心的位置。你可以尝试更换不同的特征组合(如花萼长度和宽度)来观察,可能会发现有些视角下的分离效果不如这个好。这就是多维数据的特性。
> **提示**:`KMeans`类中的`random_state`参数非常重要。因为初始质心是随机选择的,不同的随机种子可能导致不同的最终聚类结果(尤其是当数据本身聚类不明显时)。设置`random_state`可以确保每次运行得到相同的结果,这对于实验的可复现性至关重要。
现在,我们来回答之前留下的问题:**如果不知道K值怎么办?** 这是应用K-Means时最常遇到的挑战。有几种常用的方法可以帮助我们确定一个相对合理的K值:
- **肘部法则**:计算不同K值对应的簇内误差平方和,绘制曲线,寻找拐点(像手肘一样)。
- **轮廓系数法**:计算不同K值下所有样本的平均轮廓系数,选择系数最大的K。
让我们用代码来实现这两种方法,并对比结果。
```python
# 方法一:肘部法则
inertias = []
K_range = range(1, 11) # 测试K从1到10
for k in K_range:
kmeans_test = KMeans(n_clusters=k, random_state=42)
kmeans_test.fit(X_scaled)
inertias.append(kmeans_test.inertia_) # inertia_属性即簇内误差平方和
plt.figure(figsize=(12, 4))
plt.subplot(1, 2, 1)
plt.plot(K_range, inertias, 'bo-')
plt.xlabel('簇的数量 K')
plt.ylabel('簇内误差平方和 (Inertia)')
plt.title('肘部法则')
plt.grid(True)
# 方法二:轮廓系数法
silhouette_scores = []
for k in K_range:
if k == 1:
silhouette_scores.append(0) # 当K=1时,轮廓系数无定义,设为0
continue
kmeans_test = KMeans(n_clusters=k, random_state=42)
cluster_labels_test = kmeans_test.fit_predict(X_scaled)
score = silhouette_score(X_scaled, cluster_labels_test)
silhouette_scores.append(score)
plt.subplot(1, 2, 2)
plt.plot(K_range, silhouette_scores, 'ro-')
plt.xlabel('簇的数量 K')
plt.ylabel('平均轮廓系数')
plt.title('轮廓系数法')
plt.grid(True)
plt.tight_layout()
plt.show()
```
观察左边的肘部图,你会发现当K从1增加到3时,误差平方和下降得非常快;超过3之后,下降速度明显变缓。这个“拐点”出现在K=3附近,暗示3可能是一个合适的簇数。再看右边的轮廓系数图,K=3时轮廓系数达到了峰值。两种方法都指向了K=3,这与数据的真实情况吻合。在实际项目中,你需要结合业务理解和这些指标共同判断。
## 4. 评估聚类结果与高级可视化
模型训练好了,K值也确定了,但我们怎么知道聚类效果好不好呢?在无监督学习中,因为没有真实的标签作为标准答案,评估通常更具挑战性。不过,我们仍然有一些内部评估指标和方法。
**内部评估指标**
我们已经用到的**轮廓系数**就是一个非常好的内部评估指标。它结合了簇内凝聚度和簇间分离度。对于单个样本点i,其轮廓系数s(i)的计算公式为:
> s(i) = (b(i) - a(i)) / max{a(i), b(i)}
其中,a(i)是点i到同簇其他点的平均距离(凝聚度),b(i)是点i到最近其他簇中所有点的平均距离的最小值(分离度)。s(i)的取值范围在[-1, 1]之间,越接近1表示聚类越合理。
我们可以计算整个数据集的平均轮廓系数,也可以为每个样本点计算,并通过可视化来查看每个簇的轮廓系数分布。
```python
# 计算当K=3时的轮廓系数
kmeans_final = KMeans(n_clusters=3, random_state=42)
cluster_labels_final = kmeans_final.fit_predict(X_scaled)
silhouette_avg = silhouette_score(X_scaled, cluster_labels_final)
print(f"K=3时,平均轮廓系数为: {silhouette_avg:.3f}")
# 绘制每个样本的轮廓系数图(轮廓分析图)
from sklearn.metrics import silhouette_samples
sample_silhouette_values = silhouette_samples(X_scaled, cluster_labels_final)
plt.figure(figsize=(10, 7))
y_lower = 10
for i in range(3):
# 获取属于簇i的所有样本的轮廓系数,并排序
ith_cluster_silhouette_values = sample_silhouette_values[cluster_labels_final == i]
ith_cluster_silhouette_values.sort()
size_cluster_i = ith_cluster_silhouette_values.shape[0]
y_upper = y_lower + size_cluster_i
color = plt.cm.viridis(float(i) / 3)
plt.fill_betweenx(np.arange(y_lower, y_upper),
0, ith_cluster_silhouette_values,
facecolor=color, edgecolor=color, alpha=0.7)
# 在图中标注簇的编号
plt.text(-0.05, y_lower + 0.5 * size_cluster_i, str(i))
y_lower = y_upper + 10 # 为下一个簇留出10个单位的间隔
plt.axvline(x=silhouette_avg, color="red", linestyle="--", label=f'平均轮廓系数: {silhouette_avg:.3f}')
plt.xlabel("轮廓系数值")
plt.ylabel("簇标签")
plt.title("各簇的轮廓系数分布")
plt.legend()
plt.show()
```
这张图能直观地展示每个簇的“健康”程度。理想情况下,每个簇的轮廓系数条形图都应该较宽且平均值较高,并且没有太多负值(负值表示该点可能被分错了簇)。
**与真实标签对比(外部评估)**
由于Iris数据集有真实标签,我们可以进行“作弊”式的外部评估,看看聚类结果与真实品种的匹配程度。但这在真正的无监督学习场景中是无法做到的。我们可以用交叉表来观察对应关系。
```python
# 创建聚类标签与真实标签的交叉表
ct = pd.crosstab(iris_df['cluster'], iris_df['species'])
print("聚类标签 vs 真实品种 交叉表:")
print(ct)
```
你可能会发现,聚类标签(0,1,2)与真实品种标签(0,1,2)并不一定一一对应,因为K-Means不知道真实的品种名。但交叉表能清晰地显示聚类结果是否捕获了真实的结构。通常,我们会看到绝大多数样本都被正确地分组了,可能只有一个簇存在少量混淆。
**高级可视化:降维与3D视图**
为了在二维平面上更好地展示四维数据的聚类效果,我们可以使用**主成分分析**将数据降到2维或3维,然后再绘图。
```python
from sklearn.decomposition import PCA
# 使用PCA降维到2维
pca = PCA(n_components=2)
X_pca = pca.fit_transform(X_scaled)
plt.figure(figsize=(10, 6))
scatter = plt.scatter(X_pca[:, 0], X_pca[:, 1], c=cluster_labels_final, cmap='viridis', s=50, alpha=0.7, edgecolor='k')
plt.xlabel('第一主成分')
plt.ylabel('第二主成分')
plt.title('PCA降维后的K-Means聚类结果(K=3)')
plt.colorbar(scatter, label='簇标签')
# 在降维后的空间也计算并绘制质心(需要将质心也进行PCA变换)
centroids_pca = pca.transform(kmeans_final.cluster_centers_)
plt.scatter(centroids_pca[:, 0], centroids_pca[:, 1], c='red', marker='X', s=200, label='质心', edgecolor='k', linewidth=1.5)
plt.legend()
plt.show()
# 解释主成分
print(f"前两个主成分的方差解释率: {pca.explained_variance_ratio_}")
print(f"累计方差解释率: {sum(pca.explained_variance_ratio_):.3f}")
```
如果累计方差解释率较高(比如>0.8),说明这两个主成分保留了原始数据的大部分信息,那么这个二维视图就具有较好的代表性。你还可以尝试降维到3维,并用`mpl_toolkits.mplot3d`绘制3D散点图,获得更立体的视角。
## 5. 进阶技巧与常见问题排错
掌握了基础流程后,我们来看看如何提升聚类效果,以及遇到问题时该如何调试。这部分内容往往是在实际项目中积累的经验。
**处理初始质心敏感性问题**
K-Means对初始质心的选择很敏感,可能收敛到局部最优解。`sklearn`的`KMeans`类默认使用更智能的`k-means++`初始化方法(通过`init='k-means++'`参数),这比纯随机初始化效果更好、更稳定。你还可以通过`n_init`参数让算法用不同的初始质心运行多次,并选择效果最好(惯性最小)的一次作为最终结果。
```python
# 使用k-means++初始化,并运行10次选择最佳结果
kmeans_robust = KMeans(n_clusters=3, init='k-means++', n_init=10, random_state=42)
kmeans_robust.fit(X_scaled)
print(f"使用k-means++和多次运行后的惯性值: {kmeans_robust.inertia_:.2f}")
```
**处理异常值**
K-Means对异常值非常敏感,因为质心是簇内点的平均值,一个远离群体的点会把质心“拉”向它。在聚类前,识别并处理异常值是很好的实践。你可以使用简单的统计方法(如基于IQR)或专门的异常检测算法。
```python
# 示例:使用简单的Z-score方法检测异常值(假设数据已标准化)
from scipy import stats
z_scores = np.abs(stats.zscore(X_scaled))
# 找出Z-score绝对值大于3的样本(通常认为是异常值)
outliers = np.where(z_scores > 3)[0]
print(f"检测到潜在异常值索引: {outliers}")
# 在实际项目中,你可以选择删除、修正或单独处理这些点
```
**当K-Means效果不佳时**
如果你发现轮廓系数很低,或者可视化结果显示簇的形状很奇怪(比如拉得很长),这可能意味着数据不适合用K-Means。这时可以考虑其他聚类算法:
- **DBSCAN**:基于密度的算法,能发现任意形状的簇,且能识别噪声点。
- **层次聚类**:不需要预先指定K值,可以生成树状图。
- **高斯混合模型**:一种软聚类算法,假设每个簇服从高斯分布。
**一个完整的、可复用的代码模板**
最后,我把自己常用的一个K-Means分析流程封装成了一个函数,你可以直接拿去用在你的数据集上。
```python
def kmeans_analysis_pipeline(data, feature_names=None, k_range=range(2, 11)):
"""
一个完整的K-Means分析管道。
参数:
data: 二维numpy数组或DataFrame,仅包含特征。
feature_names: 特征名称列表,可选。
k_range: 要探索的K值范围。
返回:
最佳K值,对应的模型,以及评估结果字典。
"""
# 1. 数据标准化
scaler = StandardScaler()
data_scaled = scaler.fit_transform(data)
# 2. 寻找最佳K值(肘部法则+轮廓系数)
inertias = []
silhouette_scores_list = []
for k in k_range:
kmeans = KMeans(n_clusters=k, random_state=42, n_init=10)
labels = kmeans.fit_predict(data_scaled)
inertias.append(kmeans.inertia_)
if k > 1: # 轮廓系数要求至少2个簇
silhouette_scores_list.append(silhouette_score(data_scaled, labels))
else:
silhouette_scores_list.append(0)
# 可视化评估指标
fig, (ax1, ax2) = plt.subplots(1, 2, figsize=(12, 4))
ax1.plot(k_range, inertias, 'bo-')
ax1.set_xlabel('簇的数量 K')
ax1.set_ylabel('簇内误差平方和')
ax1.set_title('肘部法则')
ax1.grid(True)
ax2.plot(k_range, silhouette_scores_list, 'ro-')
ax2.set_xlabel('簇的数量 K')
ax2.set_ylabel('平均轮廓系数')
ax2.set_title('轮廓系数法')
ax2.grid(True)
plt.tight_layout()
plt.show()
# 选择轮廓系数最大的K作为最佳K(也可以结合肘部法)
best_k_idx = np.argmax(silhouette_scores_list[1:]) + 1 # 跳过K=1
best_k = k_range[best_k_idx]
print(f"建议的最佳K值为: {best_k}")
# 3. 用最佳K训练最终模型
best_kmeans = KMeans(n_clusters=best_k, random_state=42, n_init=10)
best_labels = best_kmeans.fit_predict(data_scaled)
final_silhouette = silhouette_score(data_scaled, best_labels)
print(f"最佳模型(K={best_k})的平均轮廓系数: {final_silhouette:.3f}")
# 4. 使用PCA降维进行可视化
pca = PCA(n_components=2)
data_pca = pca.fit_transform(data_scaled)
centroids_pca = pca.transform(best_kmeans.cluster_centers_)
plt.figure(figsize=(8, 6))
scatter = plt.scatter(data_pca[:, 0], data_pca[:, 1], c=best_labels, cmap='viridis', s=50, alpha=0.7, edgecolor='k')
plt.scatter(centroids_pca[:, 0], centroids_pca[:, 1], c='red', marker='X', s=200, label='质心', edgecolor='k')
plt.xlabel('第一主成分')
plt.ylabel('第二主成分')
plt.title(f'K-Means聚类结果 (K={best_k}) - PCA视图')
plt.colorbar(scatter, label='簇标签')
plt.legend()
plt.show()
# 返回结果
results = {
'best_k': best_k,
'model': best_kmeans,
'labels': best_labels,
'silhouette_score': final_silhouette,
'inertia': best_kmeans.inertia_,
'pca_data': data_pca,
'scaler': scaler
}
return results
# 使用示例(假设你有一个新的数据矩阵`my_data`)
# results = kmeans_analysis_pipeline(my_data, feature_names=my_feature_names)
```
这个函数自动化了从预处理、寻优到可视化的全过程。在实际项目中,我经常用它来快速探索一个新数据集的聚类潜力。当然,每个数据集都有其独特性,你可能需要根据具体情况调整参数或预处理步骤。比如,对于稀疏数据或文本数据(经过向量化后),余弦距离可能比欧氏距离更合适,这时可以考虑使用`KMeans`的变体或能够自定义距离度量的算法。
K-Means就像一把瑞士军刀,简单、实用,是数据科学家工具箱里的必备品。通过今天对Iris数据集的实战,我希望你不仅学会了如何调用`sklearn`的API,更重要的是理解了算法背后的思想、评估的方法以及解决问题的完整思路。下次当你面对一堆没有标签的数据时,不妨先用K-Means试试看,也许它能帮你发现意想不到的洞见。