%GAUSSM Trainable mapping, mixture of Gaussians (MoG) density estimate
%
% W = GAUSSM(A,K,R,S,M)
% W = A*GAUSSM([],K,R,S,M)
% W = A*GAUSSM(K,R,S,M)
%
% INPUT
% A Dataset
% K Number of Gaussians per class
% R,S,M Regularization parameters, 0 <= R,S <= 1, see QDC
%
% OUTPUT
% W Mixture of Gaussians density estimate
%
% DESCRIPTION
% Estimation of a PDF by the dataset A by a mixture of Gaussians procedure.
% Use is made of EMCLUST(A,QDC,K). Unlabeled objects are neglected, unless
% A is entirely unlabeled or double. Then all objects are used. If A is a
% multi-class crisp labeled dataset the densities are estimated class by
% class and then weighted and combined according their prior probabilities.
% Use +A instead of A to obtain a single set of Gaussians. In all cases,
% just single density estimator W is returned.
%
% Note that it is necessary to set the label type of A to soft labels
% (A = LABTYPE(A,'soft') in order to use the traditional EM algorithm
% based on posterior probabilities instead of using crisp labels.
%
% The mapping W may be applied to a new dataset B using DENSITY = B*W.
%
% W = A*GAUSSM
%
% uses a single Gaussian per class (K=1) and no regularisation. If
% regulariisation is desired, also K should be supplied.
%
% EXAMPLE
% a = gendatb;
% w = a*gaussm(2);
% scatterd(a)
% plotm(w)
%
% SEE ALSO (<a href="http://37steps.com/prtools">PRTools Guide</a>)
% DATASETS, MAPPINGS, QDC, MOGC, EMCLUST, PLOTM, TESTC
% Copyright: R.P.W. Duin, r.p.w.duin@37steps.com
% Faculty EWI, Delft University of Technology
% P.O. Box 5031, 2600 GA Delft, The Netherlands
% $Id: gaussm.m,v 1.8 2009/02/02 21:57:39 duin Exp $
function w = gaussm(varargin)
mapname = 'Mixture of Gaussians';
argin = shiftargin(varargin,'scalar');
argin = setdefaults(argin,[],1,0,0,[]);
if mapping_task(argin,'definition')
w = define_mapping(argin,'untrained');
w = setname(w,mapname);
else
[a,n,r,s,dim] = deal(argin{:});
if isa(a,'prdataset')
labname = getname(a);
else
labname = '';
end
if ((~isdataset(a) && ~isdatafile(a)) || ...
(getsize(a,3) ~= 1 && islabtype(a,'crisp')))
% this handles the multiclass situation
w = mclassm(a,prmapping(mfilename,n),'weight');
w = setlabels(w,labname);
w = setname(w,mapname);
else
% here we have a single class
[m,k] = getsize(a);
if n == 1
% just one Gaussian, esitmate it
[U,G] = meancov(a);
res.mean = +U;
res.cov = G;
res.prior= 1;
w = normal_map(res,labname,k,1);
else
% multiple Gaussians, run EM
[e,v] = emclust(a,qdc([],r,s,dim),n);
ncomp0 = size(v.data.mean,1);
iter = 0;
while (ncomp0 ~= n & iter < 5) % repeat until exactly n components are found
[e,v1] = emclust(a,qdc([],r,s,m),n);
ncomp1 = size(v1.data.mean,1);
if ncomp1 > ncomp0
v = v1;
ncomp0 = ncomp1;
end
iter = iter + 1;
end
res = v.data;
res.nlab = ones(n,1); % defines that all Gaussian components have to be
% combined into a single class.
w = prmapping('normal_map','trained',res,labname,k,1);
end
w = setname(w,mapname);
end
end
return