在Matlab中有效地计算欧氏距离的成对平方
给定两组d维分。 我怎样才能最有效地计算Matlab中的平方欧氏距离matrix ？
符号：设置一个由(numA,d)matrixA给出，并且设置两个由(numB,d)matrixB 。 得到的距离matrix应该是格式(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
这里通常给出的答案是基于bsxfun （参见例如[1] ）。 我提出的方法是基于matrix乘法，结果比我能find的任何可比较的algorithm快得多：
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可以通过硬编码实现replacefor -loop，例如
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）。 事实certificate，不同的algorithm应该是在不同的设置首选。 在非时间紧急的情况下，只需使用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);
不过，这相当耗时。 用二维matrix替代较小的三维matrixhelp可能是有用的。 特别是对于d = 1它提供了一种通过简单的matrix乘法计算配对差异的方法：
pairDiffs = [ones(numA,1), A ] * [B'; -ones(1,numB)];
你有什么想法吗？
