在模型选择的过程中，不仅需要确定模型中组件的数量，还需要确定协方差矩阵的类型。在本例中，将使用信息论标准来进行模型选择，具体包括Akaike信息准则（AIC）和贝叶斯信息准则（BIC）。尽管两者都能提供正确的结果，但将主要演示BIC的使用，因为它更适合在一组候选模型中识别出真实的模型。与贝叶斯方法不同，这种推断不需要先验知识。
通过随机抽样标准正态分布生成两个组件，每个组件包含n_samples个样本。其中一个组件保持球形但进行了位移和缩放处理，另一个组件则被变形以具有更一般的协方差矩阵。
import numpy as np
n_samples = 500
np.random.seed(0)
C = np.array([[0.0, -0.1], [1.7, 0.4]])
component_1 = np.dot(np.random.randn(n_samples, 2), C)
component_2 = 0.7 * np.random.randn(n_samples, 2) + np.array([-4, 1])
X = np.concatenate([component_1, component_2])
可以可视化这些不同的组件：
import matplotlib.pyplot as plt
plt.scatter(component_1[:, 0], component_1[:, 1], s=0.8)
plt.scatter(component_2[:, 0], component_2[:, 1], s=0.8)
plt.title("高斯混合模型组件")
plt.axis("equal")
plt.show()
改变组件的数量从1到6，并使用不同类型的协方差参数：“full”（每个组件都有自己的一般协方差矩阵）、“tied”（所有组件共享相同的一般协方差矩阵）、“diag”（每个组件都有自己的对角协方差矩阵）以及“spherical”（每个组件都有自己的单一方差）。使用GridSearchCV和用户定义的评分函数来评估不同的模型，该函数返回负的BIC分数，因为GridSearchCV旨在最大化分数（最大化负BIC等同于最小化BIC）。最佳参数集和估计器分别存储在best_parameters_和best_estimator_中。
from sklearn.mixture import GaussianMixture
from sklearn.model_selection import GridSearchCV
def gmm_bic_score(estimator, X):
""“Callable to pass to GridSearchCV that will use the BIC score.”""
return -estimator.bic(X)
param_grid = {
"n_components": range(1, 7),
"covariance_type": ["spherical", "tied", "diag", "full"],
}
grid_search = GridSearchCV(GaussianMixture(), param_grid=param_grid, scoring=gmm_bic_score)
grid_search.fit(X)
在Jupyter环境中，请重新运行此单元格以显示HTML表示或信任笔记本。在GitHub上，HTML表示无法渲染，请尝试使用nbviewer.org加载此页面。
为了便于绘制，可以从网格搜索的交叉验证结果中创建一个pandas.DataFrame。将BIC分数的符号反转以显示最小化它的效果。
import pandas as pd
df = pd.DataFrame(grid_search.cv_results_)[
["param_n_components", "param_covariance_type", "mean_test_score"]
]
df["mean_test_score"] = -df["mean_test_score"]
df = df.rename(columns={
"param_n_components": "组件数量",
"param_covariance_type": "协方差类型",
"mean_test_score": "BIC分数",
})
df.sort_values(by="BIC分数").head()
在当前情况下，具有2个组件和完整协方差（对应真实的生成模型）的模型具有最低的BIC分数，因此被网格搜索选中。
绘制一个椭圆来显示选定模型的每个高斯组件。为此，需要找到协方差矩阵的特征值，如covariances_属性返回的那样。这些矩阵的形状取决于covariance_type：“full”（n_components, n_features, n_features）、“tied”（n_features, n_features）、“diag”（n_components, n_features）以及“spherical”（n_components,）。
from matplotlib.patches import Ellipse
from scipy import linalg
color_iter = sns.color_palette("tab10", 2)[::-1]
Y_ = grid_search.predict(X)
fig, ax = plt.subplots()
for i, (mean, cov, color) in enumerate(zip(
grid_search.best_estimator_.means_,
grid_search.best_estimator_.covariances_,
color_iter,
)):
v, w = linalg.eigh(cov)
if not np.any(Y_ == i):
continue
plt.scatter(X[Y_ == i, 0], X[Y_ == i, 1], 0.8, color=color)
angle = np.arctan2(w[0][1], w[0][0])
angle = 180.0 * angle / np.pi # 转换为度
v = 2.0 * np.sqrt(2.0) * np.sqrt(v)
ellipse = Ellipse(mean, v[0], v[1], angle=180.0 + angle, color=color)
ellipse.set_clip_box(fig.bbox)
ellipse.set_alpha(0.5)
ax.add_artist(ellipse)
plt.title(f"选定的GMM: {grid_search.best_params_['covariance_type']}模型, {grid_search.best_params_['n_components']}组件")
plt.axis("equal")
plt.show()