Clustering

K-Means

Giả định: có $K$ cụm, mỗi cụm là một “đám” tròn quanh một tâm. K-Means tìm $K$ tâm $\boldsymbol{\mu}_1, \dots, \boldsymbol{\mu}_K$ cực tiểu tổng bình phương khoảng cách từ mỗi điểm đến tâm gần nhất (inertia):

Thuật toán Lloyd lặp hai bước đến khi không đổi:

1. Gán: mỗi điểm thuộc cụm có tâm gần nhất.
2. Cập nhật: mỗi tâm dời về trung bình các điểm của cụm.

Mỗi bước đều không làm $J$ tăng, và số cách chia là hữu hạn, nên thuật toán chắc chắn dừng. Nhưng nó chỉ dừng ở một cực tiểu địa phương: khởi tạo khác cho kết quả khác.

Hai cách giảm rủi ro kẹt:

• k-means++: chọn tâm đầu tiên ngẫu nhiên, rồi chọn mỗi tâm tiếp theo với xác suất tỉ lệ với bình phương khoảng cách từ điểm đến tâm gần nhất đã có. Nhờ vậy các tâm ban đầu thường nằm xa nhau.
• Chạy nhiều lần (n_init) và giữ kết quả có $J$ nhỏ nhất.

Chi phí mỗi vòng là $O(nKd)$, nên K-Means chạy được trên dữ liệu rất lớn; MiniBatchKMeans còn nhanh hơn nữa.

Chọn K. Vẽ $J$ theo $K$ rồi tìm “khuỷu tay” (elbow method), hoặc chọn $K$ có hệ số silhouette trung bình cao nhất (mục Đánh giá cụm). Elbow chỉ là gợi ý; nên kiểm thêm cách chia có ổn định khi đổi seed, bớt mẫu không. Thường thì $K$ còn phụ thuộc mục đích: chia khách hàng thành 4 nhóm hay 6 nhóm là quyết định kinh doanh, không phải định lý.

Giới hạn. K-Means thích cụm tròn, kích thước gần bằng nhau, cùng mật độ. Với cụm hình trăng lưỡi liềm hay cụm lồng nhau, nó chia sai một cách rất tự tin. Tâm lại là trung bình, nên một điểm ngoại lệ đủ xa có thể kéo tâm rời khỏi đám điểm của nó: cụm gồm 1, 2, 3 và 100 có tâm ở $26.5$.

Khởi tạo tệ: tự sửa hay kẹt?

Kết quả của Lloyd phụ thuộc tâm khởi tạo. Lấy dữ liệu một chiều $1, 2, 9, 10, 11$, $K = 2$, hai tâm đầu đặt tệ: $\mu_1 = 1$, $\mu_2 = 2$.

• Gán: $\{1\}$ về $\mu_1$; $\{2, 9, 10, 11\}$ về $\mu_2$.
• Cập nhật: $\mu_1 = 1$, $\mu_2 = 32/4 = 8$.
• Gán lại: điểm 2 giờ gần $\mu_1$ hơn ($1 < 6$), nên cụm thành $\{1, 2\}$ và $\{9, 10, 11\}$.
• Cập nhật: $\mu_1 = 1.5$, $\mu_2 = 10$. Gán lại không đổi, dừng với $J = 0.5 + 2 = 2.5$.

Lần này khởi tạo tệ vẫn tự sửa được, và $J$ giảm qua từng bước: $194 \to 50 \to 15 \to 2.5$. Nhưng chỉ cần dời một điểm là sự may mắn ấy biến mất (câu hỏi bên dưới).

Trong hai chiều cũng vậy: ba điểm $A(0, 2)$, $B(0, 1)$, $C(1, 0)$ với tâm đầu đặt ở $A$ và $B$ dừng ở $\{A\}$, $\{B, C\}$ với $J = 1$, trong khi $\{A, B\}$, $\{C\}$ chỉ có $J = 0.5$.

Cài một vòng Lloyd

Lấy bốn điểm một chiều 0, 1, 5, 6, hai tâm ban đầu 0 và 6. Mỗi điểm vào tâm gần nhất, nên nhóm trái gồm 0, 1 và nhóm phải gồm 5, 6. Tâm mới là 0.5 và 5.5. Inertia trước cập nhật là 2, sau cập nhật là 1. Code dưới cài đúng hai thao tác gán và cập nhật, và in ra các số vừa tính tay.

import numpy as np

values = np.array([[0.], [1.], [5.], [6.]])
centers = np.array([[0.], [6.]])
distance2 = np.sum((values[:, None, :] - centers[None, :, :]) ** 2, axis=2)
assigned = np.argmin(distance2, axis=1)
before = np.sum(distance2[np.arange(len(values)), assigned])
new_centers = np.array([values[assigned == k].mean(axis=0) for k in range(2)])
after = np.sum((values - new_centers[assigned]) ** 2)
print("gán:", assigned.tolist())
print("tâm mới:", new_centers.ravel().tolist())
print("inertia:", float(before), "->", float(after))

Kết quảgán: [0, 0, 1, 1]
tâm mới: [0.5, 5.5]
inertia: 2.0 -> 1.0

Nếu một cụm không còn điểm nào, mean của mảng rỗng là NaN, và NaN lặng lẽ lan sang mọi bước sau. Ví dụ nhỏ trên không gặp, nhưng bản cài tổng quát phải chọn trước một quy tắc cho trường hợp này.

DBSCAN

Giả định: cụm là vùng dày đặc, ngăn cách bởi vùng thưa. Hai tham số: bán kính eps và số điểm tối thiểu min_samples.

• Điểm lõi (core): có ít nhất min_samples điểm trong bán kính eps.
• Điểm biên (border): không phải lõi nhưng nằm trong bán kính của một điểm lõi.
• Điểm nhiễu (noise): không thuộc hai loại trên, gán nhãn $-1$.

Các điểm lõi nối với nhau thành cụm. DBSCAN không cần biết trước số cụm, tìm được cụm hình dạng bất kỳ và tự loại nhiễu. Điểm yếu: một eps duy nhất không hợp với dữ liệu có các cụm mật độ khác nhau. HDBSCAN khắc phục điều này và thường là lựa chọn tốt hơn khi có sẵn.

DBSCAN trên năm điểm

Trong scikit-learn, số điểm trong bán kính eps được đếm tính cả chính nó. Từ một điểm core, cụm lan qua các core lân cận rồi nhận thêm border ở rìa. Vì vậy điểm border ít hàng xóm vẫn thuộc cụm; bản cài chỉ giữ các điểm đủ ngưỡng sẽ biến mọi border thành noise, và đó đã là một thuật toán khác.

from sklearn.cluster import DBSCAN

density_points = np.array([[0.], [0.1], [0.2], [0.35], [3.]])
db = DBSCAN(eps=0.16, min_samples=3).fit(density_points)
print("nhãn cụm:", db.labels_.tolist())
print("core:", db.core_sample_indices_.tolist())
assert db.labels_[-1] == -1

Kết quảnhãn cụm: [0, 0, 0, 0, -1]
core: [1, 2]

Điểm 0.1 có ba điểm trong eps, còn 0 chỉ có hai; 0 vẫn vào cụm vì nó nối được với core 0.1. Điểm 3 quá xa nên nhãn −1. Số hiệu cụm 0, 1, … chỉ là tên tạm, có thể đổi thứ tự khi dữ liệu đổi. Một điểm border chạm được hai cụm thì thuộc cụm nào tùy thứ tự duyệt.

Hierarchical clustering

Agglomerative clustering bắt đầu với mỗi điểm là một cụm, rồi liên tục gộp hai cụm gần nhau nhất. Khoảng cách giữa hai cụm (linkage) có thể là:

• single: hai điểm gần nhất, tạo cụm dạng chuỗi dài;
• complete: hai điểm xa nhất, tạo cụm gọn;
• average: trung bình mọi cặp;
• Ward: mức tăng tổng phương sai khi gộp, gần với tinh thần của K-Means.

Kết quả là một cây gộp (dendrogram); cắt cây ở độ cao nào thì được số cụm tương ứng. Cách này cho cái nhìn nhiều tầng, nhưng tốn $O(n^2)$ bộ nhớ nên chỉ hợp với vài chục nghìn điểm.

Spectral clustering

Spectral clustering làm ba bước:

1. Dựng đồ thị tương đồng giữa các điểm, ví dụ nối mỗi điểm với các láng giềng gần.
2. Lấy vài vector riêng của ma trận Laplacian của đồ thị làm tọa độ mới.
3. Chạy K-Means trên tọa độ đó.

Hai điểm cùng nằm trên một “dải” dữ liệu sẽ gần nhau trong tọa độ mới dù xa nhau theo đường thẳng. Cách này mạnh với cụm không lồi, nhưng tốn kém với $n$ lớn.

Laplacian là $L = D - A$, với $A$ ma trận trọng số cạnh và $D$ ma trận đường chéo chứa bậc của mỗi đỉnh. Số trị riêng bằng 0 đúng bằng số thành phần liên thông, nên các vector riêng ứng với trị riêng nhỏ nhất chỉ ra cách cắt đồ thị.

Với hai vòng đồng tâm còn một cách thủ công: thêm đặc trưng $x^2 + y^2$, bình phương khoảng cách tới tâm, rồi chạy K-Means.

Gaussian mixture: gán mềm

K-Means gán mỗi điểm cứng vào đúng một cụm, và chỉ nhìn khoảng cách tới tâm. Gaussian mixture model (GMM) coi dữ liệu được sinh từ $K$ phân phối chuẩn, mỗi thành phần $k$ có trọng số $\pi_k$, trung bình $\boldsymbol\mu_k$ và hiệp phương sai $\Sigma_k$ riêng. Mỗi điểm nhận một xác suất thuộc từng thành phần (gán mềm), gọi là responsibility $\gamma_{ik}$.

GMM học bằng thuật toán EM, lặp hai bước giống Lloyd:

1. E-step: với tham số hiện tại, tính $\gamma_{ik}$, xác suất điểm $i$ thuộc thành phần $k$.
2. M-step: cập nhật $\pi_k$, $\boldsymbol\mu_k$, $\Sigma_k$ bằng trung bình có trọng số $\gamma_{ik}$.

Khởi tạo thường lấy tâm của K-Means; lặp tới khi log-likelihood gần như không tăng. Cần một cụm cứng thì gán mỗi điểm vào thành phần có $\gamma_{ik}$ lớn nhất. K-Means là trường hợp đặc biệt khi mọi thành phần có cùng hiệp phương sai tròn và rất nhỏ.

Khác biệt hiện rõ khi hai cụm có độ rộng khác nhau. Dưới đây là 300 điểm quanh 0 với độ lệch chuẩn 0.5, và 300 điểm quanh 5 với độ lệch chuẩn 2.

from sklearn.cluster import KMeans
from sklearn.metrics import adjusted_rand_score
from sklearn.mixture import GaussianMixture

gen = np.random.default_rng(0)
x1 = np.concatenate([gen.normal(0, 0.5, 300), gen.normal(5, 2.0, 300)]).reshape(-1, 1)
truth = np.r_[np.zeros(300), np.ones(300)]
km = KMeans(n_clusters=2, n_init=10, random_state=0).fit(x1)
gm = GaussianMixture(n_components=2, random_state=0).fit(x1)
narrow = int(np.argmin(gm.covariances_.ravel()))
print("Tâm K-Means:", np.sort(km.cluster_centers_.ravel()).round(2),
      " ranh giới ở", round(km.cluster_centers_.mean(), 2))
print("GMM: trung bình", gm.means_.ravel().round(2),
      " độ lệch chuẩn", np.sqrt(gm.covariances_.ravel()).round(2))
print("Điểm 1.5 thuộc cụm hẹp với xác suất", gm.predict_proba([[1.5]])[0, narrow].round(2))
print("ARI: K-Means", round(adjusted_rand_score(truth, km.labels_), 2),
      " GMM", round(adjusted_rand_score(truth, gm.predict(x1)), 2))

Kết quảTâm K-Means: [0.14 5.4 ]  ranh giới ở 2.77
GMM: trung bình [ 4.88 -0.04]  độ lệch chuẩn [2.04 0.5 ]
Điểm 1.5 thuộc cụm hẹp với xác suất 0.11
ARI: K-Means 0.79  GMM 0.96

K-Means đặt ranh giới ở giữa hai tâm, 2.77, nên các điểm từ 1.5 tới 2.7 bị xếp vào cụm hẹp. GMM biết cụm quanh 0 rất hẹp: điểm 1.5 đã cách tâm 3 độ lệch chuẩn của cụm đó, nên chỉ có xác suất 0.11 thuộc về nó. Nhờ học cả độ rộng của từng cụm, GMM đạt ARI 0.96 so với 0.79.

Đánh giá cụm

Không có nhãn thì đo độ “gọn và tách biệt”. Với điểm $i$, gọi $a(i)$ là khoảng cách trung bình đến các điểm cùng cụm, $b(i)$ là khoảng cách trung bình nhỏ nhất đến một cụm khác. Hệ số silhouette

$s$ gần 1 là điểm nằm đúng cụm, gần 0 là nằm giữa hai cụm, âm là có lẽ xếp nhầm. $b(i)$ lấy cụm khác gần nhất, không phải xa nhất: $a = 0.4$, khoảng cách trung bình tới hai cụm khác là 0.9 và 0.6 thì $s = (0.6 - 0.4)/0.6 \approx 0.33$.

Silhouette cũng “thích” cụm tròn, nên đánh giá thấp kết quả tốt của DBSCAN trên dữ liệu hình trăng (mục Dùng thư viện). Đừng chọn số cụm chỉ vì silhouette cao.

Inertia lại càng không dùng được để chọn $K$: inertia tối ưu không tăng khi $K$ tăng, nên “chọn $K$ có inertia nhỏ nhất” luôn đẩy về rất nhiều cụm.

Nếu có nhãn thật (đề thi chấm phân cụm thường có), dùng Adjusted Rand Index (ARI) hoặc Normalized Mutual Information (NMI). Nhãn chỉ dùng để chấm, không dùng khi fit. Cả hai không phụ thuộc cách đặt tên cụm: cụm “0” của mô hình có là cụm “3” của đáp án cũng không sao. Accuracy so thẳng thì khác: chỉ cần đổi tên cụm là điểm sụp.

Dùng thư viện

from sklearn.cluster import DBSCAN, KMeans
from sklearn.datasets import make_moons
from sklearn.metrics import adjusted_rand_score, silhouette_score
from sklearn.preprocessing import StandardScaler

X, y = make_moons(n_samples=400, noise=0.06, random_state=0)
X = StandardScaler().fit_transform(X)
for name, model in [("K-Means", KMeans(n_clusters=2, n_init=10, random_state=0)),
                    ("DBSCAN ", DBSCAN(eps=0.3, min_samples=5))]:
    labels = model.fit_predict(X)
    print(name, "ARI =", round(adjusted_rand_score(y, labels), 2),
          " silhouette =", round(silhouette_score(X, labels), 2))

Kết quảK-Means ARI = 0.47  silhouette = 0.5
DBSCAN  ARI = 1.0  silhouette = 0.39

DBSCAN khớp nhãn hoàn toàn theo ARI, nhưng silhouette lại chấm K-Means cao hơn: đúng điều mục Đánh giá cụm đã cảnh báo. Con số eps=0.3 cũng chỉ có nghĩa sau StandardScaler: bán kính luôn đi kèm thang đo.