## Dissimilarity embedding for indefinite inner product space

November 29, 2010 Coded in Octave
``````function  [Xnk, lambda, calX, sig, T] = IIPSembedd(d)
% syntax: [Xnk, lambda, calX, sig, T] = IIPSembedd(d)
% creates indefinite inneer product space embedding for square dissimilarity matrix d
% Xnk will be the embedding/ mapping matrix, sig the signature of the
% resulting (in)definite inner product space.
% calX is the matrix of ordered eigenvectors, lambda the associated vector
% of (sortted) eigenvalues.
% T is the gramian in the indefinite inner product space.
% Reference:
%   Pekalska, Elzbieta et al: A Generalized Kernel Approach to Dissimilarity-based Classification
%		Journal of Machine Learning Research 2 (2001), pp. 175-211
%	Pekalska, Elzbieta et al: The Dissimilarity Representation for Pattern Recognition
%		World Scientific, Singapore, December 2005
%	Goldfarb, Lev: A New Approach to Pattern Recogntion, 1985
%
n = size(d, 1);
D2 = d.^2;                          % elementwise squaring
J = zeros(n) + 1/n;                 % 1/n *  - Matrix
J = eye(n) - J;                     % get centering matrix ...
T = -0.5 *J * D2 * J;               % ... and Gramian
T = 0.5 * (T + T');                 % symmetrize T (e.g. numerical errors)

% get Eigenvalue diagonal matrix and corresponding Eigenvector-matrix
[calX, lambda] = eig(T);            % todo: lambda=real(lambda)?
lambda=real(lambda);                % just in case. should always be real for a symmetric dissimilarity matrix 'd'.
evs = diag(lambda);                 % get vector of eigenvalues from diagonal matrix
[sev, idx] = sort(-evs);            % sort!
sgn = sign(sev);                    % get dimension information for sorting
ll = find(sgn>=0, 1);
if (isempty(ll))
ll = length(sgn);
endif
sig = ([length(sgn)-ll, ll]);

% Sort eigenvalues and associated eigenvectors:
% first positive EVs (descending), follwed by zero and negatives
% (descending w.r.t. absolute value).
evs = [-sev(1:ll-1); -sev(end:-1:ll)];
lambda=evs;
idx = [idx(1:ll-1); idx(end:-1:ll)];
calX = calX(:, idx);

% Selection of important EVs, for instance by threshold 0.00001, can be inserted here

% Compute final configuration
Xnk=calX * diag(sqrt(abs(evs)));
endfunction``````