Scikit-learn GaussianMixture 实战：3种协方差矩阵对比与聚类效果可视化
Scikit-learn GaussianMixture 实战：3种协方差矩阵对比与聚类效果可视化
1. 高斯混合模型的核心概念
高斯混合模型（Gaussian Mixture Model, GMM）是一种基于概率的聚类算法，它假设数据是由多个高斯分布组合生成的。与K-means等硬聚类算法不同，GMM属于软聚类方法，能够给出每个样本属于各个簇的概率。
关键数学概念：
- 每个高斯分布由均值向量μ和协方差矩阵Σ定义
- 混合系数α表示各高斯分布的权重，满足∑α=1
- 概率密度函数：p(x) = ∑αᵢN(x|μᵢ,Σᵢ)
在实际应用中，GMM通过EM算法迭代优化这些参数：
- E步：计算样本属于各分布的概率（责任值）
- M步：根据责任值更新分布参数
2. 协方差矩阵类型解析
Scikit-learn的GaussianMixture提供了4种协方差矩阵类型，直接影响聚类形状：
| 类型 | 参数数量 | 形状特点 | 适用场景 |
|---|---|---|---|
| 'full' | D(D+1)/2 | 每个簇有完全独立的协方差矩阵 | 各簇形状方向均不同 |
| 'tied' | D(D+1)/2 | 所有簇共享同一个协方差矩阵 | 簇形状相同但位置不同 |
| 'diag' | D | 对角矩阵，特征间独立 | 各轴缩放比例不同 |
| 'spherical' | 1 | 对角且元素相同（σ²I） | 各向同性分布 |
from sklearn.mixture import GaussianMixture # 创建不同协方差类型的GMM模型 gmm_full = GaussianMixture(n_components=3, covariance_type='full') gmm_tied = GaussianMixture(n_components=3, covariance_type='tied') gmm_diag = GaussianMixture(n_components=3, covariance_type='diag')3. 实战案例：鸢尾花数据集聚类
3.1 数据准备与可视化
import matplotlib.pyplot as plt from sklearn.datasets import load_iris import numpy as np # 加载数据 iris = load_iris() X = iris.data[:, :2] # 只使用前两个特征便于可视化 y = iris.target # 可视化原始数据 plt.scatter(X[:, 0], X[:, 1], c=y, s=40, cmap='viridis') plt.title("Iris True Labels") plt.xlabel("Sepal length") plt.ylabel("Sepal width")3.2 不同协方差矩阵效果对比
我们创建可视化函数来比较三种主要协方差类型：
def plot_gmm(gmm, X, title): # 训练模型 gmm.fit(X) # 预测聚类结果 labels = gmm.predict(X) # 创建网格用于绘制决策边界 x_min, x_max = X[:, 0].min() - 1, X[:, 0].max() + 1 y_min, y_max = X[:, 1].min() - 1, X[:, 1].max() + 1 xx, yy = np.meshgrid(np.linspace(x_min, x_max, 100), np.linspace(y_min, y_max, 100)) Z = gmm.predict(np.c_[xx.ravel(), yy.ravel()]) Z = Z.reshape(xx.shape) # 绘制结果 plt.figure(figsize=(8, 6)) plt.contourf(xx, yy, Z, alpha=0.3) plt.scatter(X[:, 0], X[:, 1], c=labels, s=40, cmap='viridis') # 绘制椭圆表示协方差 for i in range(gmm.n_components): # 计算2D协方差椭圆 if gmm.covariance_type == 'full': cov = gmm.covariances_[i][:2, :2] elif gmm.covariance_type == 'tied': cov = gmm.covariances_[:2, :2] elif gmm.covariance_type == 'diag': cov = np.diag(gmm.covariances_[i][:2]) elif gmm.covariance_type == 'spherical': cov = np.eye(2) * gmm.covariances_[i] v, w = np.linalg.eigh(cov) angle = np.arctan2(w[0][1], w[0][0]) angle = 180 * angle / np.pi # 转换为度 v = 2. * np.sqrt(2.) * np.sqrt(v) ell = plt.matplotlib.patches.Ellipse(gmm.means_[i, :2], v[0], v[1], 180 + angle, color='black') ell.set_clip_box(plt.gca().bbox) ell.set_alpha(0.5) plt.gca().add_artist(ell) plt.title(title) plt.xlabel("Sepal length") plt.ylabel("Sepal width") # 比较三种主要类型 plot_gmm(gmm_full, X, "GMM with full covariance") plot_gmm(gmm_tied, X, "GMM with tied covariance") plot_gmm(gmm_diag, X, "GMM with diagonal covariance")3.3 性能指标对比
我们使用轮廓系数和BIC（贝叶斯信息准则）评估模型：
from sklearn.metrics import silhouette_score models = { 'full': gmm_full, 'tied': gmm_tied, 'diag': gmm_diag, 'spherical': GaussianMixture(n_components=3, covariance_type='spherical') } results = [] for name, model in models.items(): model.fit(X) labels = model.predict(X) sil_score = silhouette_score(X, labels) bic = model.bic(X) results.append({ 'type': name, 'silhouette': sil_score, 'BIC': bic }) # 展示结果对比 import pandas as pd df_results = pd.DataFrame(results) print(df_results.sort_values('silhouette', ascending=False))典型输出结果：
| covariance_type | silhouette_score | BIC |
|---|---|---|
| full | 0.52 | 580 |
| diag | 0.48 | 610 |
| tied | 0.45 | 595 |
| spherical | 0.42 | 625 |
4. 高级应用技巧
4.1 确定最佳聚类数
使用BIC和轮廓系数确定最优组件数：
n_components = np.arange(1, 8) bics = [] silhouettes = [] for n in n_components: gmm = GaussianMixture(n_components=n, covariance_type='full') gmm.fit(X) bics.append(gmm.bic(X)) labels = gmm.predict(X) silhouettes.append(silhouette_score(X, labels)) # 绘制结果 fig, ax1 = plt.subplots() ax2 = ax1.twinx() ax1.plot(n_components, bics, 'b-', label='BIC') ax2.plot(n_components, silhouettes, 'r-', label='Silhouette') ax1.set_xlabel('Number of components') ax1.set_ylabel('BIC', color='b') ax2.set_ylabel('Silhouette', color='r') plt.title("Optimal Number of Components")4.2 处理高维数据
当特征维度较高时，建议：
- 先使用PCA降维
- 选择'diag'或'spherical'减少参数数量
- 使用BayesianGaussianMixture自动确定组件数
from sklearn.decomposition import PCA from sklearn.mixture import BayesianGaussianMixture # 高维数据处理示例 X_highdim = iris.data # 使用全部4个特征 # 方法1：PCA降维 pca = PCA(n_components=2) X_pca = pca.fit_transform(X_highdim) gmm_pca = GaussianMixture(n_components=3, covariance_type='full') gmm_pca.fit(X_pca) # 方法2：自动确定组件数 bgmm = BayesianGaussianMixture(n_components=10, covariance_type='diag') bgmm.fit(X_highdim) print(f"实际使用的组件数: {np.sum(bgmm.weights_ > 0.01)}")5. 常见问题解决方案
问题1：模型收敛警告
解决方案：增加max_iter或tol参数，或尝试不同初始化
gmm = GaussianMixture(n_components=3, covariance_type='full', max_iter=500, tol=1e-4, init_params='kmeans')问题2：奇异矩阵错误
解决方案：使用更简单的协方差类型，或添加正则化
gmm = GaussianMixture(n_components=3, covariance_type='diag', reg_covar=1e-6) # 添加小的正则化项问题3：聚类结果不稳定
解决方案：固定random_state，或使用BayesianGaussianMixture
gmm = GaussianMixture(n_components=3, covariance_type='full', random_state=42)