c++ - Eigen 中两个矩阵之间的成对差异

Question

在 matlab/octave 矩阵之间的成对距离，例如k-means 是通过一个函数调用(参见 cvKmeans.m )计算的，以 distFunc(Codebook, X) 为参数，两个维度为 KxD 的矩阵。

在 Eigen 中，这可以通过使用广播对矩阵和一个 vector 完成，如 eigen.tuxfamily.org 中所述。 :

 (m.colwise() - v).colwise().squaredNorm().minCoeff(&index);

但是，在这种情况下，v 不仅仅是一个 vector ，而是一个矩阵。 Eigen 中用于计算两个矩阵之间所有条目的此类成对(欧几里德)距离的等效单线是什么？

Answer 1

我认为合适的解决方案是将此功能抽象为一个函数。该功能很可能是模板化的；而且它很可能会使用一个循环——毕竟这个循环真的很短。许多矩阵运算是使用循环实现的 - 这不是问题。

例如，给定你的例子......

MatrixXd p0(2, 4);
p0 <<
    1, 23, 6, 9,
    3, 11, 7, 2;

MatrixXd p1(2, 2);
p1 <<
    2, 20,
    3, 10;

然后我们可以构造一个矩阵 D 使得 D(i,j) = |p₀(i) - p₁(j)|²

MatrixXd D(p0.cols(), p0.rows());
for (int i = 0; i < p1.cols(); i++)
    D.col(i) = (p0.colwise() - p1.col(i)).colwise().squaredNorm().transpose();

我认为这很好 - 我们可以使用一些广播来避免 2 级嵌套:我们迭代 p₁ 的点，但不超过 p₀ 的点，也不超过它们的维度。

但是，如果您观察到 |p₀(i) - p_{1，则可以制作单行线}(j)|² = |p₀(i)|² + |p₁(j)|² - 2 p₀(i) ^T p₁(j)。特别是，最后一个分量只是矩阵乘法，所以 D = -2 p₀ sup>T p₁ + ...

剩下的待补空白由只依赖于行的组件组成；和一个仅依赖于列的组件:这些可以使用按行和按列的操作来表示。

那么最后的“oneliner”就是:

D = ( (p0.transpose() * p1 * -2
      ).colwise() + p0.colwise().squaredNorm().transpose()
    ).rowwise() + p1.colwise().squaredNorm();

您还可以将 rowwise/colwise 技巧替换为具有 1 vector 的(外)积。

这两种方法都会产生以下(平方)距离:

  1 410
505  10
 32 205
 50 185

您必须对最快的进行基准测试，但看到循环获胜我不会感到惊讶，而且我希望这也更具可读性。