performance - Calcule eficientemente la distancia euclidiana cuadrada por pares en Matlab

matrix distance

Question

Dado dos conjuntos de puntos d dimensional. ¿Cómo puedo calcular de manera más eficiente la matriz de distancia euclidiana cuadrada por pares en Matlab?

Notación: el conjunto uno está dado por una (numA,d) A y el conjunto dos está dado por una (numB,d) B La matriz de distancia resultante debe tener el formato (numA,numB) .

Puntos de ejemplo:

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

Answer 1

La respuesta usualmente dada aquí se basa en bsxfun (cf. ej. [1] ). Mi enfoque propuesto se basa en la multiplicación de matrices y resulta ser mucho más rápido que cualquier algoritmo comparable que pude encontrar:

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'';

Tenga en cuenta: para la constante d uno puede reemplazar el -loop por implementaciones codificadas, por ejemplo

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.

Evaluación:

%% 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;

Resultado:

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

Para una evaluación más detallada, la dimensión y el número de puntos de datos siguen la discusión a continuación (@comments). Resulta que diferentes algos deberían ser preferidos en diferentes configuraciones. En situaciones que no sean críticas, simplemente use la versión pdist2 .

Desarrollo adicional: se puede pensar en reemplazar el euclidiano cuadrado por cualquier otra métrica basada en el mismo principio:

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);

Sin embargo, esto es bastante lento. Podría ser útil reemplazar para d menor la help matriz tridimensional por d matrices bidimensionales. Especialmente para d = 1 proporciona un método para calcular la diferencia por pares mediante una simple multiplicación de matrices:

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

¿Tienes más ideas?