EM算法的MATLAB代码

下面是一个简单的 MATLAB 实现 EM 算法的示例代码：

% EM算法的MATLAB代码
% 假设有一个二元高斯混合模型

% 初始化参数
mu1 = [1; 2];
mu2 = [3; 4];
sigma1 = [1, 0; 0, 1];
sigma2 = [2, 0; 0, 2];
p = [0.4, 0.6];

% 生成样本
rng(1); % 设置随机数种子
N = 1000; % 样本数
X = [mvnrnd(mu1, sigma1, round(N * p(1)));
     mvnrnd(mu2, sigma2, round(N * p(2)))];
 
% EM算法
maxIter = 100; % 最大迭代次数
tol = 1e-6; % 收敛阈值
L = zeros(maxIter, 1); % 记录每次迭代的对数似然值
gamma = zeros(N, 2); % 初始化后验概率
for iter = 1:maxIter
    % E步：计算后验概率
    gamma(:, 1) = p(1) * mvnpdf(X, mu1, sigma1);
    gamma(:, 2) = p(2) * mvnpdf(X, mu2, sigma2);
    gamma = gamma ./ sum(gamma, 2);

    % M步：更新参数
    Nk = sum(gamma, 1);
    mu1 = sum(gamma(:, 1) .* X, 1) / Nk(1);
    mu2 = sum(gamma(:, 2) .* X, 1) / Nk(2);
    sigma1 = (X - mu1)' * (gamma(:, 1) .* (X - mu1)) / Nk(1);
    sigma2 = (X - mu2)' * (gamma(:, 2) .* (X - mu2)) / Nk(2);
    p = Nk / N;

    % 计算对数似然值并判断是否收敛
    L(iter) = sum(log(sum(bsxfun(@times, gamma, p), 2)));
    if iter > 1 && abs(L(iter) - L(iter - 1)) < tol
        break;
    end
end

% 绘制结果
figure;
subplot(2, 2, 1);
plot(X(:, 1), X(:, 2), 'r.');
hold on;
ezcontour(@(x, y)mvnpdf([x, y], mu1, sigma1), [-6, 10], [-6, 10]);
title('True cluster 1');
subplot(2, 2, 2);
plot(X(:, 1), X(:, 2), 'r.');
hold on;
ezcontour(@(x, y)mvnpdf([x, y], mu2, sigma2), [-6, 10], [-6, 10]);
title('True cluster 2');
subplot(2, 2, 3);
plot(X(:, 1), X(:, 2), 'r.');
hold on;
ezcontour(@(x, y)mvnpdf([x, y], mu1, sigma1), [-6, 10], [-6, 10]);
ezcontour(@(x, y)mvnpdf([x, y], mu2, sigma2), [-6, 10], [-6, 10]);
title('True clusters');
subplot(2, 2, 4);
plot(X(gamma(:, 1) > 0.5, 1), X(gamma(:, 1) > 0.5, 2), 'r.');
hold on;
plot(X(gamma(:, 2) > 0.5, 1), X(gamma(:, 2) > 0.5, 2), 'b.');
ezcontour(@(x, y)mvnpdf([x, y], mu1, sigma1), [-6, 10], [-6, 10]);
ezcontour(@(x, y)mvnpdf([x, y], mu2, sigma2), [-6, 10], [-6, 10]);
title('Estimated clusters');

这个代码实现了一个简单的二元高斯混合模型，并使用 EM 算法进行参数估计，最后绘制出真实聚类和估计聚类的结果。注意，这里为了简化问题，只考虑了二维数据。如果需要处理更复杂的数据，需要根据实际情况进行修改。