function result=libra_mahalanobis(x,locvct,varargin) %MAHALANOBIS computes the distance of each observation in x % from the location estimate (locvct) of the data, % relative to the shape of the data. % % Required input arguments: % x : data matrix (n observations in rows, p variables in columns) % locvct : location estimate of the data (p-dimensional vector) % cov or invcov : scatter estimate of the data or the inverse of the scatter estimate (pxp matrix) % % I/O: result=libra_mahalanobis(x,locvct,'cov',covmat,'n',size(x,1),'p',size(x,2)) % The user should only give the input arguments that have to change their default value. % The name of the input arguments needs to be followed by their value. % The order of the input arguments is of no importance. % % Examples: % result=libra_mahalanobis(x,loc,'cov',covx,'n',10) % result=libra_mahalanobis(x,loc,'p',2,'invcov',invcovx) % result=libra_mahalanobis(x,loc,'invcov',invcovx) % % Output: % A vector containing the distances of all the observations to locvct. % % This function is part of LIBRA: the Matlab Library for Robust Analysis, % available at: % http://wis.kuleuven.be/stat/robust.html % % Written by Katrien Van Driessen % Revisions by Sabine Verboven % Last update on 18/09/2003 % %Initialisation n = size(x,1); p = size(x,2); if nargin<3 error('Missing a required input variable') end counter=1; default=struct('cov',0,'invcov',NaN); list=fieldnames(default); options=default; IN=length(list); i=1; %reading the user's input if nargin>2 % %placing inputfields in array of strings % for j=1:nargin-2 if rem(j,2)~=0 chklist{i}=varargin{j}; i=i+1; end end % %Checking which default parameters have to be changed % and keep them in the structure 'options'. % while counter<=IN index=strmatch(list(counter,:),chklist,'exact'); if ~isempty(index) %in case of similarity for j=1:nargin-2 %searching the index of the accompanying field if rem(j,2)~=0 %fieldnames are placed on odd index if strcmp(chklist{index},varargin{j}) I=j; end end end options=setfield(options,chklist{index},varargin{I+1}); index=[]; end counter=counter+1; end end if size(locvct,2)==1 locvct=locvct'; %converting to a rowvector end if options.cov==0 & options.invcov==0 error('The scatter matrix or its inverse is a required input argument.') end %%%%%%MAIN%%%%%%%%% if ~isnan(options.invcov) covmat=options.invcov; if min(size(covmat))==1 covmat=diag(covmat); end else if min(size(options.cov))==1 options.cov=diag(options.cov); end covmat=pinv(options.cov); end hlp=x-repmat(locvct,n,1); dist=sum(hlp*covmat.*hlp,2)'; result=dist;