matlab中EM算法代码
时间: 2023-08-01 12:06:28 浏览: 98
EM算法是一种迭代算法,用于求解含有隐变量的概率模型参数的极大似然估计问题。下面是一个简单的MATLAB代码示例,演示了如何使用EM算法估计高斯混合模型的参数。
```matlab
% 生成高斯混合模型的数据
mu = [-5 0 5]; % 混合成分的均值
sigma = [1 1 1]; % 混合成分的标准差
p = [0.3 0.4 0.3]; % 混合系数
n = 1000; % 样本数量
data = zeros(n, 1);
for i = 1:n
% 按照混合系数随机选择成分
r = rand();
if r < p(1)
data(i) = mu(1) + randn() * sigma(1);
elseif r < p(1) + p(2)
data(i) = mu(2) + randn() * sigma(2);
else
data(i) = mu(3) + randn() * sigma(3);
end
end
% EM算法求解高斯混合模型参数
k = 3; % 成分数量
max_iter = 100; % 最大迭代次数
tol = 1e-6; % 收敛阈值
mu_init = [-4 1 3]; % 均值的初始值
sigma_init = [1 1 1]; % 标准差的初始值
p_init = [0.3 0.3 0.4]; % 混合系数的初始值
mu_hat = mu_init;
sigma_hat = sigma_init;
p_hat = p_init;
log_likelihood = -Inf;
for iter = 1:max_iter
% E步:计算每个样本属于每个成分的概率
gamma = zeros(n, k);
for j = 1:k
gamma(:, j) = p_hat(j) * normpdf(data, mu_hat(j), sigma_hat(j));
end
gamma = gamma ./ sum(gamma, 2);
% M步:更新模型参数
N = sum(gamma, 1);
mu_hat = sum(bsxfun(@times, gamma, data), 1) ./ N;
sigma_hat = sqrt(sum(bsxfun(@times, gamma, bsxfun(@minus, data, mu_hat).^2), 1) ./ N);
p_hat = N ./ n;
% 计算对数似然函数值
log_likelihood_new = sum(log(sum(bsxfun(@times, gamma, [p_hat(1)*normpdf(data, mu_hat(1), sigma_hat(1)); p_hat(2)*normpdf(data, mu_hat(2), sigma_hat(2)); p_hat(3)*normpdf(data, mu_hat(3), sigma_hat(3))]'), 2));
if log_likelihood_new - log_likelihood < tol
break;
end
log_likelihood = log_likelihood_new;
end
% 绘制真实数据分布和估计的混合模型分布
x = linspace(-10, 10, 100);
y_true = p(1)*normpdf(x, mu(1), sigma(1)) + p(2)*normpdf(x, mu(2), sigma(2)) + p(3)*normpdf(x, mu(3), sigma(3));
y_estimated = p_hat(1)*normpdf(x, mu_hat(1), sigma_hat(1)) + p_hat(2)*normpdf(x, mu_hat(2), sigma_hat(2)) + p_hat(3)*normpdf(x, mu_hat(3), sigma_hat(3));
figure;
plot(x, y_true);
hold on;
plot(x, y_estimated);
legend('True distribution', 'Estimated distribution');
```
阅读全文