gmm_em_auto.m

高斯滤波器 matlab toolbox For GMMs and Gaussian kernels

function [g, w] = gmm_em_auto(s, K, N, type)
%function [g, w] = gmm_em_auto(s, K, N)
%
% INPUTS:
%   s - samples
%   K - number of components in resultant gmm
%   N - number of iterations of EM
%
% OUTPUT:
%   g - resultant gmm
%   w - negative log-likelihood of fit
%
% Tim Bailey 2005.
if nargin == 3, type = 1; end

NA = 10; % number of total search iterations
NK = 5; % number of kmeans iterations
NE = 5; % number of trial EM iterations

% Auto search for initial gmm components
switch type
case 1 % kmeans generates entire gmm
    Ps = repmat(eye(size(s,1))*eps, [1,1,K]);
    wb = 1e99;
    for i=1:NA
        g.x = s(:, ceil(stratified_random(K)*N));
        [g.x,g.P,g.w] = kmeans(s, g.x, NK);
        g.P = g.P + Ps; % protect against rank deficiency
        [gem, w] = gmm_em(s, g, NE);
        if w < wb
            gbest = gem;
            wb = w;
        end
    end
    
case 2 % kmeans generates means only, P and w are computed from ensemble statistics
    g.w = ones(1,K)/K;
    [xe,Pe] = sample_mean(s);
    g.P = repmat(Pe/K, [1,1,K]);
    
    wb = 1e99;
    for i=1:NA
        g.x = s(:, ceil(stratified_random(K)*N));
        [g.x,dmy,dmy] = kmeans(s, g.x, NK);
        [gem, w] = gmm_em(s, g, NE);
        if w < wb
            gbest = gem;
            wb = w;
        end
    end    
    
otherwise
    error('Invalid initialisation type')
end

% Optimise best gmm
[g, w] = gmm_em(s, gbest, N);