【问题标题】:Is there a way to calculate the following specified matrix by avoiding loops? in R or Matlab有没有办法通过避免循环来计算以下指定的矩阵?在 R 或 Matlab 中
【发布时间】:2012-06-19 00:45:05
【问题描述】:

我有一个 N×M 矩阵 X,我需要计算一个 N×N 矩阵 Y

Y[i, j] = sum((X[i,] - X[j,]) ^ 2)   0 <= i,j <= N

现在,我必须使用嵌套循环来处理 O(n2)。我想知道是否有更好的方法,比如使用矩阵运算。

更一般地说,sum(....) 可以是一个函数,fun(x1,x 2) 其中x1x2 是 M×1 向量。

【问题讨论】:

标签: r matlab vectorization euclidean-distance


【解决方案1】:

您可以使用expand.grid 获取可能对的data.frame:

X <- matrix(sample(1:5, 50, replace=TRUE), nrow=10)

row.ind <- expand.grid(1:dim(X)[1], 1:dim(X)[2])

然后apply 沿着每一对使用一个函数:

myfun <- function(n) {
  sum((X[row.ind[n, 1],] - X[row.ind[n, 2],])^2)
}

Y <- matrix(unlist(lapply(1:nrow(row.ind), myfun)), byrow=TRUE, nrow=nrow(X))


> Y
      [,1] [,2] [,3] [,4] [,5]
 [1,]    0   28   15   31   41
 [2,]   31   28   33   30   33
 [3,]   28    0   15    7   19
 [4,]   33   30   19   34   11
 [5,]   15   15    0   12   22
 [6,]   10   19   10   21   20
 [7,]   31    7   12    0    4
 [8,]   16   17   16   13    2
 [9,]   41   19   22    4    0
[10,]   14   11   28    9    2
> 

我敢打赌有更好的方法,但现在是星期五,我累了!

【讨论】:

  • @ralu 我不想考虑复杂性,但你可以在函数中替换任何你想要的东西来简化矩阵数学,但这避免了循环结构并且相当直接地多如果需要,线程 w/mclapply。而且,就像我说的,今天是星期五!
  • 也许更好的方法是使用dist(),因为我们只是在计算行对之间的欧几里得距离。这将是 Ansari 的“作弊”Matlab 答案的 R 等价物。
  • 但是 OP 要求提供一个通用解决方案,其中 sum(...) 可以是任何东西。
  • 是的,我想要一个通用的解决方案。我正在编写一个 SVM 函数,矩阵是 R 中 quadprog 的 Dmat,总和...片段是核函数。
【解决方案2】:
(x[i]-x[j])^2 = x[i]² - 2*x[i]*x[j] + x[j]²

中间部分只是矩阵乘法-2*X*tran(X)(矩阵),其他部分只是 vetrors,你必须在每个元素上运行它

这有 O(n^2.7) 或任何矩阵乘法复杂度

伪代码:

vec=sum(X,rows).^2
Y=X * tran(X) * -2 
for index [i,j] in Y:
   Y[i,j] =  Y[i,j] + vec[i]+vec[y]

【讨论】:

    【解决方案3】:

    在 MATLAB 中,对于您的特定 f,您可以这样做:

    Y = pdist(X).^2;
    

    对于非“作弊”版本,试试这样的(MATLAB):

    [N, M] = size(X);
    f = @(u, v) sum((u-v).^2);
    helpf = @(i, j) f(X(i, :), X(j, :))
    Y = arrayfun(helpf, meshgrid(1:N, 1:N), meshgrid(1:N, 1:N)');
    

    使用特定函数sum(...) 有更有效的方法,但您的问题是您想要通用函数f 的通用方法。通常,此操作将是每个向量对操作复杂度的 O(n^2) 倍,因为这就是需要完成的操作数。如果f是特殊形式,部分计算结果可以重复使用。

    【讨论】:

      猜你喜欢
      • 1970-01-01
      • 1970-01-01
      • 1970-01-01
      • 2013-01-11
      • 1970-01-01
      • 1970-01-01
      • 2013-02-07
      • 1970-01-01
      • 1970-01-01
      相关资源
      最近更新 更多