function [M0, M1, acc, time] = traingmms(DATA, LABELS, NGAUSS)
% [M0, M1, acc, time] = traingmms(DATA, LABELS, NGAUSS)
%     Train two Gaussian mixture models and evaluate their accuracy 
%     on the training data.  DATA is rows of training feature vectors, 
%     with corresponding 0/1 labels in LABELS.  Train a pair of GMMs
%     with NGAUSS components, return in M1 (for label 0) and M2 (for 
%     label 1).  Report the % correct on the training data and return 
%     it in acc.  Return total elapsed time in time.
% 2003-06-30 dpwe@ee.columbia.edu muscontent practical

% Start execution timer
tic;

if nargin < 3
  NGAUSS = 5;
end

dd0 = DATA(LABELS == 0,:);
dd1 = DATA(LABELS == 1,:);

ndim = size(DATA, 2);

M0 = gmm(ndim, NGAUSS, 'diag');
M1 = gmm(ndim, NGAUSS, 'diag');

options = foptions;
options(14) = 5;  % 5 iterations of k-means
M0 = gmminit(M0, dd0, options);
M1 = gmminit(M1, dd1, options);

options = zeros(1,18);
options(14) = 20;  % 20 iterations of EM
M0 = gmmem(M0, dd0, options);
M1 = gmmem(M1, dd1, options);

% Data likelihood for 0-labelled frames
LR0 = log(gmmprob(M1, dd0)./gmmprob(M0, dd0));
% Data likelihood for 1-labelled frames
LR1 = log(gmmprob(M1, dd1)./gmmprob(M0, dd1));

% Overall classification accuracy on training data
acc = mean([(LR0' < 0),(LR1' > 0)]);
disp(['Accuracy on training data = ',num2str(round(1000*acc)/10), '%']);

% How long did it take?
time = toc;
disp(['Elapsed time = ', num2str(time),' secs']);