如何在Python中基于距离矩阵评估自定义层次聚类的质量

问题背景 #

使用变量间的相关距离构建自定义层次聚类流程时，无法直接用scikit-learn、SciPy等库的内置评估指标（Davies–Bouldin Index、Calinski-Harabasz Index、Silhouette Score）评估聚类质量——这些指标默认要求输入原始特征矩阵，而非预计算的距离矩阵。

解决方案：手动实现基于距离矩阵的评估指标 #

这三个指标的核心逻辑均基于样本间的距离，完全可以直接利用预计算的距离矩阵手动实现，无需依赖原始特征。以下是具体实现和改造后的聚类流程：

1. 手动实现三个评估指标 #

import numpy as np

def silhouette_score_from_distances(dist_matrix, labels):
    """基于距离矩阵计算轮廓系数"""
    n_samples = len(labels)
    silhouette_vals = []
    
    for i in range(n_samples):
        # 计算当前样本到同簇其他样本的平均距离(a)
        same_cluster_mask = labels == labels[i]
        same_cluster_mask[i] = False  # 排除自身
        a = dist_matrix[i, same_cluster_mask].mean() if np.any(same_cluster_mask) else 0
        
        # 计算当前样本到最近异簇的平均距离(b)
        unique_clusters = np.unique(labels)
        b_vals = []
        for c in unique_clusters:
            if c != labels[i]:
                other_cluster_mask = labels == c
                b_val = dist_matrix[i, other_cluster_mask].mean()
                b_vals.append(b_val)
        b = min(b_vals) if b_vals else 0
        
        # 计算当前样本的轮廓值
        if a == 0 and b == 0:
            silhouette_val = 0
        else:
            silhouette_val = (b - a) / max(a, b)
        silhouette_vals.append(silhouette_val)
    
    return np.mean(silhouette_vals)

def davies_bouldin_score_from_distances(dist_matrix, labels):
    """基于距离矩阵计算Davies-Bouldin Index"""
    unique_clusters = np.unique(labels)
    n_clusters = len(unique_clusters)
    cluster_intra_means = []
    
    # 计算每个簇的内部平均距离
    for c in unique_clusters:
        cluster_mask = labels == c
        cluster_dist = dist_matrix[cluster_mask][:, cluster_mask]
        # 取上三角（排除重复和自身）计算平均
        triu_indices = np.triu_indices(cluster_dist.shape[0], k=1)
        mean_dist = cluster_dist[triu_indices].mean() if len(triu_indices[0]) > 0 else 0
        cluster_intra_means.append(mean_dist)
    
    # 计算每个簇与其他簇的DB相似度
    db_vals = []
    for i in range(n_clusters):
        max_sim = 0
        for j in range(n_clusters):
            if i != j:
                # 计算簇i和j之间的平均距离
                cluster_i_mask = labels == unique_clusters[i]
                cluster_j_mask = labels == unique_clusters[j]
                inter_dist = dist_matrix[cluster_i_mask][:, cluster_j_mask].mean()
                if inter_dist == 0:
                    sim = np.inf
                else:
                    sim = (cluster_intra_means[i] + cluster_intra_means[j]) / inter_dist
                if sim > max_sim:
                    max_sim = sim
        db_vals.append(max_sim)
    
    return np.mean(db_vals)

def calinski_harabasz_score_from_distances(dist_matrix, labels):
    """基于距离矩阵计算Calinski-Harabasz Index"""
    n_samples = len(labels)
    unique_clusters = np.unique(labels)
    n_clusters = len(unique_clusters)
    
    # 计算总离散度（所有样本对距离平方和的1/2，利用距离矩阵对称性）
    triu_indices = np.triu_indices(n_samples, k=1)
    total_ss = np.sum(dist_matrix[triu_indices] ** 2)
    
    # 计算簇内离散度之和
    within_ss = 0
    for c in unique_clusters:
        cluster_mask = labels == c
        cluster_dist = dist_matrix[cluster_mask][:, cluster_mask]
        triu_cluster = np.triu_indices(cluster_dist.shape[0], k=1)
        within_ss += np.sum(cluster_dist[triu_cluster] ** 2)
    
    # 簇间离散度 = 总离散度 - 簇内离散度
    between_ss = total_ss - within_ss
    
    # CH指数公式
    if within_ss == 0:
        return np.inf
    ch_score = (between_ss / (n_clusters - 1)) / (within_ss / (n_samples - n_clusters))
    return ch_score

2. 改造后的层次聚类评估流程 #

将原代码替换为调用上述自定义函数即可：

# 假设已导入必要库：from sklearn import cluster as skclus
DistanceMatrix = 1 - abs(pipe6.Stats())
N = DistanceMatrix.shape[0]
scores = {
    'davies_bouldin': [],
    'calinski_harabasz': [],
    'silhouette': []
}

for K in range(2, N-1):
    hclus = skclus.FeatureAgglomeration(n_clusters=K,
                                        linkage="complete",
                                        metric="precomputed")
    hclus.fit(DistanceMatrix)
    clusters = hclus.labels_
    
    # 调用自定义评估函数
    db_score = davies_bouldin_score_from_distances(DistanceMatrix, clusters)
    ch_score = calinski_harabasz_score_from_distances(DistanceMatrix, clusters)
    sil_score = silhouette_score_from_distances(DistanceMatrix, clusters)
    
    scores['davies_bouldin'].append(db_score)
    scores['calinski_harabasz'].append(ch_score)
    scores['silhouette'].append(sil_score)

指标解读 #

• 轮廓系数：数值越接近1表示聚类效果越好，接近-1表示样本被分到错误簇，0表示簇之间重叠。
• Davies-Bouldin Index：数值越小表示聚类效果越好，0是理论最优值。
• Calinski-Harabasz Index：数值越大表示簇间差异越大、簇内越紧凑，聚类效果越好。

内容的提问来源于stack exchange，提问作者Eduardo Augusto Pereira Gomes