当前位置：首页 > news >正文

高斯混合模型参数估计使用EM算法对一维或者高维的高斯混合模型GMM进行密度函数参数估计

news 2026/8/11 5:34:56

高斯混合模型参数估计使用EM算法对一维或者高维的高斯混合模型GMM进行密度函数参数估计，算法使用kmeans参数初始化，使用EM算法进行迭代求解算法为Matlab编写，逻辑清晰，注释全面

直接上干货。咱们今天聊聊怎么用Matlab实现高斯混合模型（GMM）的参数估计，重点是用EM算法迭代求解。整个过程从kmeans初始化开始，手把手教你从零搭建一个能跑通的GMM模型。

先看参数初始化部分。这里用kmeans确定初始聚类中心，简单粗暴但有效：

function [mu, sigma, alpha] = init_params(data, K) [idx, centers] = kmeans(data, K); % K个聚类中心 mu = centers'; % 转置为列向量 sigma = zeros(1, 1, K); alpha = zeros(K, 1); for k = 1:K cluster_data = data(idx == k); sigma(:,:,k) = var(cluster_data); % 每个簇的方差 alpha(k) = sum(idx == k)/length(data); % 混合系数 end end

这段代码有意思的地方在于sigma的计算。注意三维数组的维度处理——虽然是一维数据，但保持高维扩展性。混合系数alpha直接按样本比例分配，这比随机初始化靠谱多了。

接下来是EM算法的核心循环。先看E步计算后验概率：

function gamma = e_step(x, mu, sigma, alpha) K = length(alpha); N = length(x); prob = zeros(N, K); for k = 1:K prob(:,k) = alpha(k) * normpdf(x, mu(k), sqrt(sigma(:,:,k))); end gamma = prob ./ sum(prob, 2); % 按行归一化 end

这里的gamma矩阵就是传说中的"责任值"。注意normpdf的参数处理——sigma需要开平方，很多新手会在这里踩坑。sum(prob,2)这个操作让每行的概率和为1，确保后验概率的合法性。

M步的参数更新才是重头戏：

function [mu_new, sigma_new, alpha_new] = m_step(x, gamma) [N, K] = size(gamma); mu_new = sum(gamma.*x) ./ sum(gamma); % 加权平均 sigma_new = zeros(1,1,K); for k = 1:K diff = x - mu_new(k); sigma_new(:,:,k) = (gamma(:,k)' * (diff.^2)) / sum(gamma(:,k)); end alpha_new = mean(gamma)'; % 直接取平均 end

mu的计算用了点乘技巧，比循环高效。sigma更新时注意保持维度一致性——这里故意用三维数组，方便后续扩展到多维情况。alpha更新简单到令人发指，就是各分量的平均责任值。

主函数把这几部分串起来：

function [mu, sigma, alpha, log_likelihood] = gmm_em(x, K, max_iter, tol) % 初始化 [mu, sigma, alpha] = init_params(x, K); prev_loglik = -inf; for iter = 1:max_iter % E步 gamma = e_step(x, mu, sigma, alpha); % M步 [mu, sigma, alpha] = m_step(x, gamma); % 计算似然值 log_likelihood = sum(log(sum(alpha' .* normpdf(x, mu', sqrt(sigma)), 2))); % 收敛判断 if abs(log_likelihood - prev_loglik) < tol break; end prev_loglik = log_likelihood; end end

这里有个小技巧：对数似然的计算避免了下溢问题。注意normpdf的参数排列方式，mu需要转置才能正确广播。收敛条件设为对数似然变化小于阈值，通常取1e-6左右就够用了。

实际跑代码时要注意几个坑：