gausspdf.m

模式识别工具箱,希望对大家有用!

function y = gausspdf(x,mu,sigma,lambda)
%
% Y = gausspdf(X,MU,SIGMA)
%
% High dimensional version of normpdf. Given the MU and SIGMA, the
% density at points X is computed. It is assumed that all objects are
% row objects.
% Per default, just the inverse of the covariance matrix is computed.
%
% Y = gausspdf(X,MU,SIGMA,lambda)
%
% When lambda>0, the covariance matrix is regularized by:
%
%  sigma' = sigma + lambda*eye(dim)
%
% For lambda<0 the pseudo-inverse is used.
%
% default 1D mu=0, sigma=1

% Copyright: D. Tax, R.P.W. Duin, [e-mail]
% Faculty of Applied Physics, Delft University of Technology
% P.O. Box 5046, 2600 GA Delft, The Netherlands


if nargin<4
  lambda = [];
end
if nargin<3
  sigma = 1;
end
if nargin<2
  mu = 0;
end
if nargin<1
  error('Requires at least on input argument.');
end

% check sizes:
[rx,cx] = size(x);
[rm,cm] = size(mu);
[rs,cs] = size(sigma);
if (cx ~= cm) | (cx ~= cs)
  error('Number of columns in X, MU and SIGMA should be equal');
end

% mean, inv-covariance and mahanalobis distance:
X = x - ones(rx,1)*mu;
if isempty(lambda)
  Sinv = inv(sigma);
else
  if (lambda<0)
    Sinv = pinv(sigma);
  else
    Sinv = inv(sigma + lambda*eye(cs));
  end
end
d2 = sum((X*Sinv).*X,2);

% normalize to pdf:
detSinv = det(Sinv);
if (detSinv<0)  % annoying when near-singular cov.matrix
  detSinv = -detSinv;  % hack hack hack
end

y = sqrt(detSinv/((2*pi)^cx))*exp(-d2/2);

return