在Matlab中高效计算成对平方欧几里得距离

Question

在Matlab中高效计算成对平方欧几里得距离

performancematlabmatrixdistanceeuclidean-distance

13

给定两组 d 维点，如何在 Matlab 中最有效地计算 成对平方欧几里得距离矩阵？

符号表示： 第一组由一个 (numA,d) 矩阵 A 给出，第二组由一个 (numB,d) 矩阵 B 给出。结果距离矩阵的格式应为 (numA,numB)。

示例点：

d = 4;            % dimension
numA = 100;       % number of set 1 points
numB = 200;       % number of set 2 points
A = rand(numA,d); % set 1 given as matrix A
B = rand(numB,d); % set 2 given as matrix B

- matheburg

1

你看过pdist2函数吗？http://www.mathworks.com/help/stats/pdist2.html - rayryeng

@rayryeng，可以看一下我回答中的评估部分 :) - matheburg

2个回答

1

对于欧几里得距离的平方，也可以使用以下公式

||a-b||^2 = ||a||^2 + ||b||^2 - 2<a,b>

其中<a,b>是向量a和b的点积。

nA = sum( A.^2, 2 ); %// norm of A's elements
nB = sum( B.^2, 2 ); %// norm of B's elements
distMat = bsxfun( @plus, nA, nB' ) - 2 * A * B' ;

最近我被告知，在R2016b版本中，计算平方欧几里得距离的这种方法比现有的方法更快。

- Shai

网页内容由stack overflow 提供, 点击上面的

可以查看英文原文，
原文链接

- matheburg · Accepted Answer

通常给出的答案是基于bsxfun（参见例如[1]）。我提出的方法基于矩阵乘法，结果比我找到的任何可比较算法都要快得多。

helpA = zeros(numA,3*d);
helpB = zeros(numB,3*d);
for idx = 1:d
    helpA(:,3*idx-2:3*idx) = [ones(numA,1), -2*A(:,idx), A(:,idx).^2 ];
    helpB(:,3*idx-2:3*idx) = [B(:,idx).^2 ,    B(:,idx), ones(numB,1)];
end
distMat = helpA * helpB';

请注意：对于常量d，可以通过硬编码实现替换for-循环，例如：

helpA(:,3*idx-2:3*idx) = [ones(numA,1), -2*A(:,1), A(:,1).^2, ... % d == 2
                          ones(numA,1), -2*A(:,2), A(:,2).^2 ];   % etc.

评估：

%% create some points
d = 2; % dimension
numA = 20000;
numB = 20000;
A = rand(numA,d);
B = rand(numB,d);

%% pairwise distance matrix
% proposed method:
tic;
helpA = zeros(numA,3*d);
helpB = zeros(numB,3*d);
for idx = 1:d
    helpA(:,3*idx-2:3*idx) = [ones(numA,1), -2*A(:,idx), A(:,idx).^2 ];
    helpB(:,3*idx-2:3*idx) = [B(:,idx).^2 ,    B(:,idx), ones(numB,1)];
end
distMat = helpA * helpB';
toc;

% compare to pdist2:
tic;
pdist2(A,B).^2;
toc;

% compare to [1]:
tic;
bsxfun(@plus,dot(A,A,2),dot(B,B,2)')-2*(A*B');
toc;

% Another method: added 07/2014
% compare to ndgrid method (cf. Dan's comment)
tic;
[idxA,idxB] = ndgrid(1:numA,1:numB);
distMat = zeros(numA,numB);
distMat(:) = sum((A(idxA,:) - B(idxB,:)).^2,2);
toc;

结果：

Elapsed time is 1.796201 seconds.
Elapsed time is 5.653246 seconds.
Elapsed time is 3.551636 seconds.
Elapsed time is 22.461185 seconds.

针对维度和数据点数量的更详细评估，请参见下面的讨论（@comments）。结果表明，在不同的情境中应该优先选择不同的算法。在非时间关键的情况下，只需使用pdist2版本。

进一步发展： 可以考虑用基于相同原理的任何其他度量替换平方欧几里得距离：

help = zeros(numA,numB,d);
for idx = 1:d
    help(:,:,idx) = [ones(numA,1), A(:,idx)     ] * ...
                    [B(:,idx)'   ; -ones(1,numB)];
end
distMat = sum(ANYFUNCTION(help),3);

然而，这样做非常耗时。将 3 维矩阵 help 替换为 d 个 2 维矩阵可能会很有用，特别是对于较小的 d。特别是当 d = 1 时，可以通过简单的矩阵乘法计算成对差异。

pairDiffs = [ones(numA,1), A ] * [B'; -ones(1,numB)];

你有其他的想法吗？