基于基因表达数据，如何从Scikit-learn的GaussianMixture获取类mclust的分类属性？

嘿，刚好我之前也处理过类似的基因表达数据聚类需求，用Scikit-learn的GaussianMixture完全能实现和R中mclust一样的分类效果，下面给你一步步拆解怎么做：

1. 先做好准备工作：导入库+整理数据 #

首先得导入需要的工具库，同时把你的基因表达数据调整成GMM要求的格式——它需要输入二维数组（形状为(n_samples, n_features)），哪怕是单特征的基因表达数据，也要转成列向量形式。

示例代码：

import numpy as np
from sklearn.mixture import GaussianMixture

# 假设你的基因表达数据是一维数组（比如从文件读取后的结果）
gene_expr = np.random.normal(loc=[1, 5], scale=[1, 1], size=(100, 1)).ravel()
# 转成GMM要求的二维格式
X = gene_expr.reshape(-1, 1)

2. 拟合双峰高斯混合模型 #

因为你要定位双峰分布，所以直接指定n_components=2（对应两个高斯组分），然后拟合模型即可：

# 初始化GMM模型，固定random_state保证结果可复现
gmm = GaussianMixture(n_components=2, random_state=42)
gmm.fit(X)

3. 获取核心分类结果（对应mclust的$classification） #

在Scikit-learn的GaussianMixture里，有两个关键方式可以拿到分类信息：

方法1：直接获取样本的组分标签 #

用predict()方法就能得到每个样本被分配到的组分（标签为0或1），这和R中mclust返回的$classification完全对应：

# 每个样本对应的分类标签，形状为(n_samples,)
cluster_labels = gmm.predict(X)

cluster_labels就是你要的核心分类属性，每个元素对应一个样本所属的双峰簇。

方法2：查看样本的组分归属概率（可选） #

如果你想知道每个样本被分到不同组分的概率（类似mclust的$z矩阵），可以用predict_proba()方法：

# 每个样本属于两个组分的概率，形状为(n_samples, 2)
cluster_probs = gmm.predict_proba(X)

比如cluster_probs[i, 0]就是第i个样本属于第0个组分的概率，cluster_probs[i, 1]是属于第1个组分的概率。

4. 额外验证：确认双峰分布的拟合效果 #

如果你想验证模型是否真的捕捉到了双峰分布，可以绘制原始数据的直方图和拟合的高斯曲线：

import matplotlib.pyplot as plt

# 绘制原始数据的直方图
plt.hist(X, bins=20, density=True, alpha=0.5, label='Gene Expression')

# 生成整体拟合的GMM概率密度曲线
x_range = np.linspace(X.min(), X.max(), 1000).reshape(-1, 1)
pdf = np.exp(gmm.score_samples(x_range))
plt.plot(x_range, pdf, label='Fitted GMM')

# 分别绘制两个组分的高斯曲线（带权重）
for i in range(gmm.n_components):
    mean = gmm.means_[i][0]
    cov = gmm.covariances_[i][0][0]
    # 单个组分的概率密度
    component_pdf = (1 / np.sqrt(2 * np.pi * cov)) * np.exp(-(x_range - mean)**2 / (2 * cov))
    # 乘以组分权重，得到在混合模型中的贡献
    plt.plot(x_range, component_pdf * gmm.weights_[i], label=f'Component {i+1}')

plt.legend()
plt.xlabel('Gene Expression Level')
plt.ylabel('Density')
plt.show()

小提醒 #

• 如果不确定数据到底是双峰还是多峰，可以用GMM的bic()或aic()方法计算准则值，选择值最小的n_components来确定最优组分数量。
• 固定random_state参数可以让每次运行的分类结果一致，方便后续复现分析。

内容的提问来源于stack exchange，提问作者Jeff The Liu