哈喽,大家好~
今儿和大家来聊聊中层次聚类,即:Hierarchical Clustering。
举个很简单的例子,比如有一堆小球散落在地上,目标是把它们分成若干簇(cluster)。
层次聚类就像你拿着放大镜观察这些小球:
-
自底向上:最开始每个小球都是一簇(最细粒度),你不断把“最近的两簇”合并,合并后再找下一个最近的,直到只剩一簇为止。整个过程会记录每一步合并是谁和合并时的距离,这就形成了一棵树(dendrogram,树状图)。你可以在树上任意“截断”得到不同数量的簇。 -
自顶向下:一开始把所有小球看成一簇,然后逐步切开(这在实践中少用)。
关键问题是:如何定义“最近”两个簇?
这就是 linkage(连接方式)和距离度量。常见的 linkage 有:
-
single(最小距离),两簇中任意两点的最小距离; -
complete(最大距离),任意两点的最大距离; -
average(平均距离),两簇间所有点对距离的平均; -
Ward,合并后导致簇内方差增加最小,更像在做贪心优化:每步选择使得总体 SSE 增量最小的合并。
要注意,除了 Ward 外,层次聚类通常是一种“贪心”过程,而不是对某个全局可微目标函数的直接全局优化(因此可能陷入局部最优)。但 Ward 方法可以看作是在每一步贪心地减少聚类内方差,这点与“优化算法”有天然联系~
直观表达
两点间常用距离
几种 linkage 的形式
对簇 :
-
single: -
complete: -
average: -
Ward(基于均值 ): Ward 的合并代价等价于聚合前后簇内平方和(SSE)的增量,属于每步贪心最小化的形式。
时间复杂度方面,朴素实现是 O(n^3),但用堆或更复杂数据结构可降低至 O(n^2 log n) 或 O(n^2)。
对于中小规模数据,几百到几千样本用朴素实现就够了。
案例实战
下面我们用 PyTorch 从零实现自底向上的层次聚类,支持 single/complete/average/Ward,最后通过一些可视化表达,让大家更好的理解,原始数据分布、不同截断下的聚类结果、dendrogram、距离矩阵热力图、轮廓系数曲线等。
代码与流程说明会一步步给出,大家可以自己操作起来~
1) 数据集
import numpy as np
import torch
import matplotlib.pyplot as plt
from scipy.cluster.hierarchy import dendrogram
from sklearn.metrics import silhouette_score
from matplotlib import cm
def make_complex_data(n1=1200, n2=1000, n3=800, n4=700, random_seed=42):
np.random.seed(random_seed)
# Gaussian blob 1
X1 = np.random.randn(n1,2)*0.25 + np.array([0.0, 0.0])
# Elongated Gaussian blob 2
A = np.array([[0.6, 0.3], [0.3, 0.2]]) # covariance-like transform
X2 = np.dot(np.random.randn(n2,2), A) + np.array([3.0, 1.0])
# Ring (circle) cluster 3
angles = np.random.rand(n3)*2*np.pi
r = 1.6 + 0.2*np.random.randn(n3)
X3 = np.stack([r*np.cos(angles), r*np.sin(angles)], axis=1) + np.array([-1.5, 2.0])
# Small tight blob 4
X4 = np.random.randn(n4,2)*0.12 + np.array([-2.5, -1.5])
X = np.vstack([X1, X2, X3, X4])
labels = np.array([0]*n1 + [1]*n2 + [2]*n3 + [3]*n4)
return X.astype(np.float32), labels
X_np, true_labels = make_complex_data()
X = torch.from_numpy(X_np) # PyTorch tensor
n_samples = X.shape[0]
print("samples:", n_samples)
我们合成了 4 类数据:紧密高斯、拉长高斯、环形与小紧簇,构成具有不同形状和尺度的复杂聚类场景,可以很好地展示层次聚类,尤其 Ward 能区分紧簇但对非球状簇可能表现不一。
2) 实现自底向上的层次聚类
实现思路:
-
维护一个距离矩阵 D(初始为点到点的距离),每个点开始为一个簇。 -
每一步找到两个簇 (i,j) 使得距离最小,合并它们,记录合并信息。 -
更新距离矩阵:对被合并簇与其他簇的距离按 linkage 更新。 -
对于 Ward,使用均值与簇大小代价公式,避免计算所有点对距离以提高效率。
def pairwise_distances_torch(X):
# X: (n, d)
# returns (n, n) distance matrix
with torch.no_grad():
D = torch.cdist(X, X, p=2)
return D
def hierarchical_clustering(X, method='average'):
"""
X: torch.tensor (n, d)
method: 'single', 'complete', 'average', 'ward'
returns:
- linkage_matrix: (n-1, 4) numpy array compatible with scipy dendrogram
- cluster_assignments(k): function to cut tree into k clusters
"""
n = X.shape[0]
D = pairwise_distances_torch(X).numpy()
# set diagonal to inf
np.fill_diagonal(D, np.inf)
# clusters: list of lists of indices
clusters = [[i] for i in range(n)]
sizes = [1]*n
# for ward we keep cluster means
means = X.numpy().copy() # each row initially
linkage = []
current_cluster_id = n # new cluster ids for scipy format
# maintain an active mask
active = [True]*n
for step in range(n-1):
# find pair with minimal distance among active clusters
# Build current cluster distance matrix:
m = len(clusters)
# compute pairwise distances between clusters
best_i, best_j, best_d = -1, -1, np.inf
for i in range(m):
for j in range(i+1, m):
Ci = clusters[i]
Cj = clusters[j]
if method == 'single':
d = np.min(D[np.ix_(Ci, Cj)])
elif method == 'complete':
d = np.max(D[np.ix_(Ci, Cj)])
elif method == 'average':
d = np.mean(D[np.ix_(Ci, Cj)])
elif method == 'ward':
# Ward distance via cluster means
mu_i = means[Ci].mean(axis=0)
mu_j = means[Cj].mean(axis=0)
ni = len(Ci); nj = len(Cj)
d = (ni*nj)/(ni+nj) * np.sum((mu_i-mu_j)**2)
else:
raise ValueError("unknown method")
if d < best_d:
best_d = d; best_i = i; best_j = j
# merge clusters best_i and best_j
Ci = clusters[best_i]; Cj = clusters[best_j]
new_cluster = Ci + Cj
# record for linkage: [idx1, idx2, dist, new_cluster_size]
# scipy expects cluster indices: for original points 0..n-1, for merged clusters use increasing ids
linkage.append([best_i, best_j, best_d, len(new_cluster)]) # temporary ids; we'll convert later
# replace clusters: remove larger idx first
if best_i < best_j:
first, second = best_i, best_j
else:
first, second = best_j, best_i
clusters.pop(second); clusters.pop(first)
clusters.append(new_cluster)
return linkage
上述朴素实现中,我们以 cluster lists 直接查询点对距离(借助初始点对距离矩阵 D),每步寻找最小合并对。
3) 使用 SciPy 快速生成 linkage,并计算不同截断的聚类结果
为了可视化和性能,我们用 SciPy 的 linkage (它内部做了优化并返回标准矩阵)来生成 dendrogram,但聚类思想与我们上述实现一致。
from scipy.cluster.hierarchy import linkage as scipy_linkage, fcluster
def compute_linkage_and_clusters(X_np, method='average'):
Z = scipy_linkage(X_np, method=method, metric='euclidean')
return Z
Z_avg = compute_linkage_and_clusters(X_np, method='average')
Z_ward = compute_linkage_and_clusters(X_np, method='ward')
4) 可视化分析
import seaborn as sns
# 1) 原始数据分布(真实标签)
plt.figure(figsize=(6,6))
colors = cm.tab10(true_labels % 10)
plt.scatter(X_np[:,0], X_np[:,1], c=colors, s=40, edgecolor='k')
plt.title("原始数据分布(真实标签)", fontsize=14)
plt.xlabel("x1"); plt.ylabel("x2")
plt.show()
# 2) dendrogram (Ward)
plt.figure(figsize=(12, 4))
dendrogram(Z_ward, truncate_mode='level', p=20, leaf_rotation=90., leaf_font_size=10.)
plt.title("Ward linkage 的 dendrogram(截断显示)", fontsize=14)
plt.xlabel("样本或簇"); plt.ylabel("距离")
plt.show()
# 3) 不同 K 下的聚类结果(使用 average linkage)
plt.figure(figsize=(12,10))
Ks = [2,3,4,5]
for i,k in enumerate(Ks):
labels_k = fcluster(Z_avg, k, criterion='maxclust') - 1 # 0..k-1
plt.subplot(2,2,i+1)
cmap = cm.get_cmap('tab10')
for label in range(k):
pts = X_np[labels_k==label]
plt.scatter(pts[:,0], pts[:,1], s=40, color=cmap(label), edgecolor='k', label=f"cluster {label}")
plt.title(f"Average linkage, K={k}")
plt.legend()
plt.tight_layout()
plt.show()
# 4) 距离矩阵热力图(按 avg linkage 聚类顺序重排)
order = np.argsort(fcluster(Z_avg, 4, criterion='maxclust')) # reorder by cluster id
D = np.linalg.norm(X_np[:,None,:]-X_np[None,:,:], axis=2)
D_ordered = D[np.ix_(order, order)]
plt.figure(figsize=(8,6))
plt.imshow(D_ordered, cmap='viridis')
plt.colorbar()
plt.title("按聚类顺序重排的距离矩阵热力图")
plt.xlabel("样本 index (reordered)"); plt.ylabel("样本 index (reordered)")
plt.show()
# 5) silhouette 随 K 变化(2..8)
scores = []
Ks = list(range(2,9))
for k in Ks:
labels_k = fcluster(Z_avg, k, criterion='maxclust') - 1
sc = silhouette_score(X_np, labels_k)
scores.append(sc)
plt.figure(figsize=(8,4))
plt.plot(Ks, scores, marker='o', color='crimson')
plt.xticks(Ks)
plt.xlabel("K (簇数)"); plt.ylabel("Silhouette Score")
plt.title("Average linkage 的 Silhouette Score 随 K 变化")
plt.grid(True)
plt.show()
原始数据分布:
可以直观看到 4 个簇的不同形状(环形、拉长、紧密等)。这是我们要从无标签数据恢复的目标。
dendrogram:
树状图展示了整个合并过程,横轴是合并步骤,纵轴是合并时的距离(或代价)。你可以在某个高度上“剪断”树得到不同数量的簇。高枝合并表明簇之间差异较大。
不同 K 下聚类结果:
展示平均连接(average linkage)在 K=2,3,4,5 下的分组情况。通过对比这些子图,可以看到在某些 K 下聚类更符合直觉(例如 K=4 时与真实簇比较接近),某些 K 会把环形簇与其他簇合并。
距离矩阵热力图:
先按聚类结果重排样本,再画样本间距离矩阵。块状(block)结构代表簇的存在——对角线上亮或暗的方块表示簇内样本彼此接近或远离,簇间区域则颜色不同。这个图很直观地显示了簇内一致性与簇间差异。
Silhouette Score 曲线:
轮廓系数越高代表簇分割效果越好(单个样本的轮廓在 -1 到 1 之间,越接近 1 越好)。这张图可以辅助确定合适的 K(通常取峰值位置或较平稳的高值)。
把 Ward 看成“优化算法”的角度
Ward linkage 每一步选择合并时,使得簇内总平方和(SSE)增加最小。
可写成:
因此 Ward 的层次聚类可以看作是一种贪心的、逐步优化的算法:每一步都选择“对全局 SSE 增量最小”的那对簇合并。虽然它并不是在一次性对一个全局可微目标做全局最小化(不是梯度下降意义上的连续优化),但其合并策略与优化思想密切相关,是“逐步贪心最小化某个评价指标”的代表。
-
优点:对于偏向球形且方差差异不大的簇,Ward 往往产生更均匀、紧凑的簇。 -
缺点:对非球形簇(如环形)表现不佳;层次聚类的贪心本质也可能导致局部次优的合并顺序。
总结
层次聚类通过不断合并或拆分簇得到一棵树,能在多个尺度上观察数据结构;
linkage 的定义决定聚类合并的“贪心策略”,Ward 可被视为每步最小化簇内方差的贪心优化;

