DAY 17 聚类-CSDN博客

本文链接：https://blog.csdn.net/Era_mo/article/details/147730723
知识点

聚类的指标
聚类常见算法：kmeans聚类、dbscan聚类、层次聚类
三种算法对应的流程
实际在论文中聚类的策略不一定是针对所有特征，可以针对其中几个可以解释的特征进行聚类，得到聚类后的类别，这样后续进行解释也更加符合逻辑。
聚类的流程
标准化数据
选择合适的算法，根据评估指标调参（）
KMeans 和层次聚类的参数是K值，选完k指标就确定
DBSCAN 的参数是 eps 和min_samples，选完他们出现k和评估指标
以及层次聚类的 linkage准则等都需要仔细调优。
除了经典的评估指标，还需要关注聚类出来每个簇对应的样本个数，避免太少没有意义。
将聚类后的特征添加到原数据中
原则t-sne或者pca进行2D或3D可视化
作业： 对心脏病数据集进行聚类。
输入：
## 预处理流程回顾
#1. 导入库
import pandas as pd
import pandas as pd    #用于数据处理和分析，可处理表格数据。
import numpy as np     #用于数值计算，提供了高效的数组操作。
import matplotlib.pyplot as plt    #用于绘制各种类型的图表
import seaborn as sns   #基于matplotlib的高级绘图库，能绘制更美观的统计图形。

# 设置中文字体（解决中文显示问题）
plt.rcParams['font.sans-serif'] = ['SimHei']  # 设置字体为黑体
plt.rcParams['axes.unicode_minus'] = False    # 正常显示负号

#2. 读取数据查看数据信息--理解数据
data = pd.read_csv(r'heart.csv')    #读取数据
print("数据基本信息：")
data.info()
print("\n数据前5行预览：")
print(data.head())

#3. 缺失值处理
# 先筛选字符串变量 
discrete_features = data.select_dtypes(include=['object']).columns.tolist()
print("\n离散变量：")
print(discrete_features)
# 依次查看内容
for feature in discrete_features:
    print(f"\n{feature}的唯一值：")
    print(data[feature].value_counts())
#本数据集中不纯在离散变量
# thal 标签编码
thal_mapping = {
    1: 1,
    2: 2,
    3: 3,
}
data['thal'] = data['thal'].map(thal_mapping)
# slope的独热编码，记得需要将bool类型转换为数值
data = pd.get_dummies(data, columns=['slope'])
data2 = pd.read_csv(r"heart.csv") # 重新读取数据，用来做列名对比
list_final = [] # 新建一个空列表，用于存放独热编码后新增的特征名
for i in data.columns:
    if i not in data2.columns:
       list_final.append(i) # 这里打印出来的就是独热编码后的特征名
for i in list_final:
    data[i] = data[i].astype(int) # 这里的i就是独热编码后的特征名
print(list_final) # 打印出来的就是独热编码后的特征名
# 布尔矩阵显示缺失值，这个方法返回一个布尔矩阵
print(data.isnull()) 
print(data.isnull().sum())   # 统计每一列缺失值的数量
#填补缺失值
continuous_features = data.select_dtypes(include=['int64', 'float64']).columns.tolist()  #把筛选出来的列名转换成列表
print("\n连续变量：")
print(continuous_features)

for feature in continuous_features:     
    mode_value = data[feature].mode()[0]            #获取该列的众数。
    data[feature] = data[feature].fillna(mode_value)          #用众数填充该列的缺失值，inplace=True表示直接在原数据上修改。
print(data.isnull().sum())   # 统计每一列缺失值的数量
# 4. 异常值处理
#异常值一般不处理，或者结合对照试验处理和不处理都尝试下，但是论文中要写这个，作为个工作量
data.info()  #查看数据基本信息
#此数据集无缺失值
#划分数据集
from sklearn.model_selection import train_test_split
X = data.drop('target', axis=1)  # 特征
y = data['target']  # 目标变量
X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.2, random_state=42)  # 划分训练集和测试集，test_size表示测试集占比，random_state表示随机种子，保证每次划分结果一致。

import numpy as np
import pandas as pd
from sklearn.cluster import KMeans, DBSCAN, AgglomerativeClustering
from sklearn.preprocessing import StandardScaler
from sklearn.decomposition import PCA
import matplotlib.pyplot as plt
import seaborn as sns

# 标准化数据（聚类前通常需要标准化）
scaler = StandardScaler()
X_scaled = scaler.fit_transform(X)
# X_scaled

# 评估不同 k 值下的指标
k_range = range(2, 11)  # 测试 k 从 2 到 10
inertia_values = []
silhouette_scores = []
ch_scores = []
db_scores = []
agglo_scores = []
from sklearn.metrics import silhouette_score, calinski_harabasz_score, davies_bouldin_score
for k in k_range:
    # KMeans
    kmeans = KMeans(n_clusters=k, random_state=42)
    kmeans_labels = kmeans.fit_predict(X_scaled)
    inertia_values.append(kmeans.inertia_)  # 计算 KMeans 的轮廓系数
    silhouette = silhouette_score(X_scaled, kmeans_labels)  # 计算 KMeans 的轮廓系数
    silhouette_scores.append(silhouette)
    ch = calinski_harabasz_score(X_scaled, kmeans_labels)  # 计算 KMeans 的 CH 指数
    ch_scores.append(ch)
    db = davies_bouldin_score(X_scaled, kmeans_labels)  # 计算 KMeans 的 DB 指数
    db_scores.append(db)
    print(f"k={k}, 惯性: {kmeans.inertia_:.2f}, 轮廓系数: {silhouette:.3f}, CH 指数: {ch:.2f}, DB 指数: {db:.3f}")

# 绘制评估指标图
plt.figure(figsize=(15, 10))

# 肘部法则图（Inertia）
plt.subplot(2, 2, 1)
plt.plot(k_range, inertia_values, marker='o')
plt.title('肘部法则确定最优聚类数 k（惯性，越小越好）')
plt.xlabel('聚类数 (k)')
plt.ylabel('惯性')
plt.grid(True)

# 轮廓系数图
plt.subplot(2, 2, 2)
plt.plot(k_range, silhouette_scores, marker='o', color='orange')
plt.title('轮廓系数确定最优聚类数 k（越大越好）')
plt.xlabel('聚类数 (k)')
plt.ylabel('轮廓系数')
plt.grid(True)

# CH 指数图
plt.subplot(2, 2, 3)
plt.plot(k_range, ch_scores, marker='o', color='green')
plt.title('Calinski-Harabasz 指数确定最优聚类数 k（越大越好）')
plt.xlabel('聚类数 (k)')
plt.ylabel('CH 指数')
plt.grid(True)

# DB 指数图
plt.subplot(2, 2, 4)
plt.plot(k_range, db_scores, marker='o', color='red')
plt.title('Davies-Bouldin 指数确定最优聚类数 k（越小越好）')
plt.xlabel('聚类数 (k)')
plt.ylabel('DB 指数')
plt.grid(True)

plt.tight_layout()  # 调整子图间距
plt.show()

# 1. 肘部法则图： 找下降速率变慢的拐点，这里都差不多
# 2. 轮廓系数图：找局部最高点，这里选6
# 3. CH指数图： 找局部最高点，这里都差不多
# 4. DB指数图：找局部最低点，这里选4 5 6都行

# 综上，选择6比较合适。

# 提示用户选择 k 值
selected_k = 6

# 使用选择的 k 值进行 KMeans 聚类
kmeans = KMeans(n_clusters=selected_k, random_state=42)
kmeans_labels = kmeans.fit_predict(X_scaled)
X['KMeans_Cluster'] = kmeans_labels

# 使用 PCA 降维到 2D 进行可视化
pca = PCA(n_components=2)
X_pca = pca.fit_transform(X_scaled)
X_pca_df = pd.DataFrame(X_pca, columns=['PCA1', 'PCA2'])
X_pca_df['KMeans_Cluster'] = kmeans_labels

# KMeans 聚类结果可视化
plt.figure(figsize=(6, 5))
sns.scatterplot(x='PCA1', y='PCA2', hue='KMeans_Cluster', data=X_pca_df, palette='viridis')
plt.title(f'KMeans Clustering with k={selected_k} (PCA Visualization)')
plt.xlabel('PCA Component 1')
plt.ylabel('PCA Component 2')
plt.show()

# 打印 KMeans 聚类标签的前几行
print(f"KMeans Cluster labels (k={selected_k}) added to X:")
print(X[['KMeans_Cluster']].value_counts())

# DBSCAN聚类
# 评估不同 eps 和 min_samples 下的指标
eps_range = np.arange(0.3, 0.8, 0.1)  # 测试 eps 从 0.3 到 0.7
min_samples_range = range(3, 8)  # 测试 min_samples 从 3 到 7
results = []

for eps in eps_range:
    for min_samples in min_samples_range:
        dbscan = DBSCAN(eps=eps, min_samples=min_samples)
        dbscan_labels = dbscan.fit_predict(X_scaled)
        # 计算簇的数量（排除噪声点 -1）
        n_clusters = len(np.unique(dbscan_labels)) - (1 if -1 in dbscan_labels else 0)
        # 计算噪声点数量
        n_noise = list(dbscan_labels).count(-1)
        # 只有当簇数量大于 1 且有有效簇时才计算评估指标
        if n_clusters > 1:
            # 排除噪声点后计算评估指标
            mask = dbscan_labels != -1
            if mask.sum() > 0:  # 确保有非噪声点
                silhouette = silhouette_score(X_scaled[mask], dbscan_labels[mask])
                ch = calinski_harabasz_score(X_scaled[mask], dbscan_labels[mask])
                db = davies_bouldin_score(X_scaled[mask], dbscan_labels[mask])
                results.append({
                    'eps': eps,
                    'min_samples': min_samples,
                    'n_clusters': n_clusters,
                    'n_noise': n_noise,
                    'silhouette': silhouette,
                    'ch_score': ch,
                    'db_score': db
                })
                print(f"eps={eps:.1f}, min_samples={min_samples}, 簇数: {n_clusters}, 噪声点: {n_noise}, "
                      f"轮廓系数: {silhouette:.3f}, CH 指数: {ch:.2f}, DB 指数: {db:.3f}")
        else:
            print(f"eps={eps:.1f}, min_samples={min_samples}, 簇数: {n_clusters}, 噪声点: {n_noise}, 无法计算评估指标")

# 将结果转为 DataFrame 以便可视化和选择参数
results_df = pd.DataFrame(results, columns=['eps', 'min_samples', 'n_clusters', 'n_noise', 'silhouette', 'ch', 'db'])
print(results_df)
#无少数簇
# # 绘制评估指标图，增加点论文中的工作量
# plt.figure(figsize=(15, 10))

# # 轮廓系数图
# plt.subplot(2, 2, 1)
# plt.plot(results_df['eps'], results_df['silhouette'], marker='o')
# plt.title('轮廓系数确定最优 eps（越大越好）')
# plt.xlabel('eps')
# plt.ylabel('轮廓系数')
# plt.grid(True)

# # CH 指数图
# plt.subplot(2, 2, 2)
# plt.plot(results_df['eps'], results_df['ch'], marker='o', color='orange')
# plt.title('Calinski-Harabasz 指数确定最优 eps（越大越好）')
# plt.xlabel('eps')
# plt.ylabel('CH 指数')
# plt.grid(True)

# # DB 指数图
# plt.subplot(2, 2, 3)
# plt.plot(results_df['eps'], results_df['db'], marker='o', color='green')
# plt.title('Davies-Bouldin 指数确定最优 eps（越小越好）')
# plt.xlabel('eps')
# plt.ylabel('DB 指数')
# plt.grid(True)

# # 簇数量图
# plt.subplot(2, 2, 4)
# plt.plot(results_df['eps'], results_df['n_clusters'], marker='o', color='red')
# plt.title('簇数量确定最优 eps（越多越好）')
# plt.xlabel('eps')
# plt.ylabel('簇数量')
# plt.grid(True)

# plt.tight_layout()  # 调整子图间距
# plt.show()

# 层次聚类
# 标准化数据
scaler = StandardScaler()
X_scaled = scaler.fit_transform(X)
# 评估不同 n_clusters 下的指标
n_clusters_range = range(2, 11)  # 测试簇数量从 2 到 10
silhouette_scores = []
ch_scores = []
db_scores = []

for n_clusters in n_clusters_range:
    agglo = AgglomerativeClustering(n_clusters=n_clusters, linkage='ward')  # 使用 Ward 准则合并簇
    agglo_labels = agglo.fit_predict(X_scaled)
    
    # 计算评估指标
    silhouette = silhouette_score(X_scaled, agglo_labels)
    ch = calinski_harabasz_score(X_scaled, agglo_labels)
    db = davies_bouldin_score(X_scaled, agglo_labels)
    
    silhouette_scores.append(silhouette)
    ch_scores.append(ch)
    db_scores.append(db)
    
    print(f"n_clusters={n_clusters}, 轮廓系数: {silhouette:.3f}, CH 指数: {ch:.2f}, DB 指数: {db:.3f}")

# 绘制评估指标图
plt.figure(figsize=(15, 5))

# 轮廓系数图
plt.subplot(1, 3, 1)
plt.plot(n_clusters_range, silhouette_scores, marker='o')
plt.title('轮廓系数确定最优簇数（越大越好）')
plt.xlabel('簇数量 (n_clusters)')
plt.ylabel('轮廓系数')
plt.grid(True)

# CH 指数图
plt.subplot(1, 3, 2)
plt.plot(n_clusters_range, ch_scores, marker='o')
plt.title('Calinski-Harabasz 指数确定最优簇数（越大越好）')
plt.xlabel('簇数量 (n_clusters)')
plt.ylabel('CH 指数')
plt.grid(True)

# DB 指数图
plt.subplot(1, 3, 3)
plt.plot(n_clusters_range, db_scores, marker='o')
plt.title('Davies-Bouldin 指数确定最优簇数（越小越好）')
plt.xlabel('簇数量 (n_clusters)')
plt.ylabel('DB 指数')
plt.grid(True)

plt.tight_layout()
plt.show()

# 提示用户选择 n_clusters 值（这里可以根据图表选择最佳簇数）
selected_n_clusters = 4
# 使用选择的簇数进行 Agglomerative Clustering 聚类
agglo = AgglomerativeClustering(n_clusters=selected_n_clusters, linkage='ward')
agglo_labels = agglo.fit_predict(X_scaled)
X['Agglo_Cluster'] = agglo_labels

# 使用 PCA 降维到 2D 进行可视化
pca = PCA(n_components=2)
X_pca = pca.fit_transform(X_scaled)
X_pca_df = pd.DataFrame(X_pca, columns=['PCA1', 'PCA2'])
X_pca_df['Agglo_Cluster'] = agglo_labels

# Agglomerative Clustering 聚类结果可视化
plt.figure(figsize=(6, 5))
sns.scatterplot(x='PCA1', y='PCA2', hue='Agglo_Cluster', data=X_pca_df, palette='viridis')
plt.title(f'Agglomerative Clustering with n_clusters={selected_n_clusters} (PCA Visualization)')
plt.xlabel('PCA Component 1')
plt.ylabel('PCA Component 2')
plt.show()

# 打印 Agglomerative Clustering 聚类标签的分布
print(f"Agglomerative Cluster labels (n_clusters={selected_n_clusters}) added to X:")
print(X[['Agglo_Cluster']].value_counts())

# 层次聚类的树状图可视化
from scipy.cluster import hierarchy
import matplotlib.pyplot as plt

# 假设 X_scaled 是标准化后的数据
# 计算层次聚类的链接矩阵
Z = hierarchy.linkage(X_scaled, method='ward')  # 'ward' 是常用的合并准则

# 绘制树状图
plt.figure(figsize=(10, 6))
hierarchy.dendrogram(Z, truncate_mode='level', p=3)  # p 控制显示的层次深度
# hierarchy.dendrogram(Z, truncate_mode='level')  # 不用p这个参数，可以显示全部的深度
plt.title('Dendrogram for Agglomerative Clustering')
plt.xlabel('Cluster Size')
plt.ylabel('Distance')
plt.show()
输出：