source: distools/lpdistm.m @ 25

Last change on this file since 25 was 10, checked in by bduin, 14 years ago
File size: 2.5 KB
RevLine 
[10]1%LPDISTM  l_p (p > 0) (Non)-Metric Distance Matrix
2%
3%     D = LPDISTM (A,B,P)
4%       OR
5%     D = LPDISTM (A,B)
6%       OR
7%     D = LPDISTM (A,P)
8%       OR
9%     D = LPDISTM (A)
10%
11% INPUT
12%   A   NxK Matrix or dataset
13%   B   MxK Matrix or dataset
14%   P   Parameter, P > 0
15%
16% OUTPUT
17%   D   NxM Dissimilarity matrix or dataset
18%
19% DESCRIPTION
20% Computation of the distance matrix D between two sets of vectors, A and B.
21% Distances between vectors X and Y are computed using the lp distance:
22%     d(X,Y) = (sum (|X_i - Y_i|.^P))^(1/P)
23%                i
24% If P = Inf, then the max-norm distance is computed:
25%     d(X,Y) = max (|X_i - Y_i|)
26%
27% If A and B are datasets, then D is a dataset as well with the labels defined
28% by the labels of A and the feature labels defined by the labels of B. If A is
29% not a dataset, but a matrix of doubles, then D is also a matrix of doubles.
30%
31% DEFAULT
32%   P = 1
33%   B = A
34%
35% REMARKS
36%   P >= 1      => D is metric
37%   P in (0,1)  => D is non-metric; D.^P is metric and l1-embeddable
38%   P = 1/2     => D is city block / Euclidean distance
39%
40% SEE ALSO
41%   FLPDISTM, EUDISTM
42%
43
44% Copyright: Elzbieta Pekalska, ela.pekalska@googlemail.com
45% Faculty EWI, Delft University of Technology and
46% School of Computer Science, University of Manchester
47
48
49function D = lpdistm (A,B,p)
50
51bisa = 0;
52if nargin < 2,
53  p = 1;
54  B = A;
55  bisa = 1;
56else
57  if nargin < 3,
58    if max (size(B)) == 1,
59      p = B;
60      B = A;
61      bisa = 1;
62    else
63      p = 1;
64    end
65  end
66end
67
68if p <= 0,
69  error ('The parameter p must be positive.');
70end
71
72isda = isdataset(A);
73isdb = isdataset(B);
74a    = +A;
75b    = +B;
76[ra,ca] = size(a);
77[rb,cb] = size(b);
78
79if ca ~= cb,
80  error ('The matrices should have the same number of columns.');
81end
82
83
84D = zeros(ra,rb);
85if p < Inf,
86  for i=1:rb
87    %if ~rem(i,50), fprintf('.'); end
88    D(:,i) = sum(abs(repmat(b(i,:),ra,1) - a).^p,2);
89  end
90else
91  for i=1:rb
92    %if ~rem(i,50), fprintf('.'); end
93    D(:,i) = max(abs(repmat(b(i,:),ra,1) - a),[],2);
94  end
95end
96
97
98%fprintf('\n');
99if p < Inf,
100  D = D.^(1/p);
101end
102
103% Check numerical inaccuracy
104D (find (D < eps)) = 0;   % Make sure that distances are nonnegative
105if bisa,
106  D = 0.5*(D+D');         % Make sure that distances are symmetric for D(A,A)
107end
108
109% Set object labels and feature labels
110if xor(isda, isdb),
111  prwarning(1,'One matrix is a dataset and the other not. ')
112end
113if isda,
114  if isdb,
115    D = setdata(A,D,getlab(B));
116  else
117    D = setdata(A,D);
118  end
119  D.name = 'Distance matrix';
120  if ~isempty(A.name)
121    D.name = [D.name ' for ' A.name];
122  end
123end
124return
Note: See TracBrowser for help on using the repository browser.