如何缩短欧氏距离计算的处理时间

How to improve processing time for euclidean distance calculation

我正在尝试计算具有相同列数(变量)和不同行数(观察值)的两个数据帧之间的加权欧氏距离(平方)。

计算公式如下:

DIST[m,i] <- sum(((DATA1[m,] - DATA2[i,]) ^ 2) * lambda[1,])

我特别需要用特定的权重 (lambda) 乘以每个 somatory 包裹。

下面 运行s 提供的代码正确,但如果我在数百次迭代中使用它,则需要大量处理时间。昨天,我花了 18 个小时才使用包含此计算的函数的多次迭代来创建图形。使用 library(profvis) profvis({ my code }) 我看到代码的这个特定部分占用了大约 80% 的处理时间。

我阅读了很多有关如何使用并行和矢量化操作减少处理时间的内容,但由于权重 lamb#,我不知道如何在这种特殊情况下实现它们。

有人可以帮助我减少使用此代码的处理时间吗?

有关代码和数据结构的更多信息可以在下面作为注释提供的代码中找到。

# Data frames used to calculate the euclidean distances between each observation 
#   from DATA1 and each observation from DATA2.
# The euclidean distance is between a [600x50] and a [8X50] dataframes, resulting 
#   in a [600X8] dataframe.
DATA1 <- matrix(rexp(30000, rate=.1), ncol=50) #[600x50]
DATA2 <- matrix(rexp(400, rate=.1), ncol=50) #[8X50]

# Weights used for each of the 50 variables to calculate the weighted 
#   euclidean distance.
# Can be a vector of different weights or a scalar of the same weight 
#   for all variables.
lambda <- runif(n=50, min=0, max=10)   ## length(lambda) > 1
# lambda=1   ## length(lambda) == 1

if (length(lambda) > 1) {
  as.numeric(unlist(lambda))
  lambda <- as.matrix(lambda)
  lambda <- t(lambda)
}

nrows1 <- nrow(DATA1)
nrows2 <- nrow(DATA2) 

# Euclidean Distance calculation
DIST <- matrix(NA, nrow=nrows1, ncol=nrows2 )  
for (m in 1:nrows1) {
  for (i in 1:nrows2) {
    if (length(lambda) == 1) { 
      DIST[m, i] <- sum((DATA1[m, ] - DATA2[i, ])^2) 
    }
    if (length(lambda) > 1){ 
      DIST[m, i] <- sum(((DATA1[m, ] - DATA2[i, ])^2) * lambda[1, ])
    }
    next
  }
  next
}

在所有建议之后,结合@MDWITT(对于 length(lambda > 1) 和@F.Privé(对于 length(lambda == 1))的答案,最终解决方案只用了一分钟 运行,而原来的代码花了我一个半小时 运行,在一个更大的代码中有那个计算。对于那些感兴趣的人,这个问题的最终代码是:

#Data frames used to calculate the euclidean distances between each observation from DATA1 and each observation from DATA2.
#The euclidean distance is between a [600x50] and a [8X50] dataframes, resulting in a [600X8] dataframe.
DATA1 <- matrix(rexp(30000, rate=.1), ncol=50) #[600x50]
DATA2 <- matrix(rexp(400, rate=.1), ncol=50) #[8X50]

#Weights used for each of the 50 variables to calculate the weighted euclidean distance.
#Can be a vector of different weights or a scalar of the same weight for all variables.
#lambda <- runif(n = 50, min = 0, max = 10)   ##length(lambda) > 1
lambda = 1   ##length(lambda) == 1

nrows1 <- nrow(DATA1)
nrows2 <- nrow(DATA2) 

#Euclidean Distance calculation
DIST <- matrix(NA, nrow = nrows1, ncol = nrows2)  

if (length(lambda) > 1){
  as.numeric(unlist(lambda))
  lambda <- as.matrix(lambda)
  lambda <- t(lambda)

  library(Rcpp)
  cppFunction('NumericMatrix weighted_distance (NumericMatrix x, NumericMatrix y, NumericVector lambda){

              int n_x = x.nrow();
              int n_y = y.nrow();


              NumericMatrix DIST(n_x, n_y);

              //begin the loop

              for (int i = 0 ; i < n_x; i++){
              for (int j = 0  ; j < n_y ; j ++) {
              double d = sum(pow(x.row(i) - y.row(j), 2)*lambda);
              DIST(i,j) = d;
              }
              }
              return (DIST) ;
  }')

    DIST <- weighted_distance(DATA1, DATA2, lambda = lambda)}


  if (length(lambda) == 1) { 
    DIST <- outer(rowSums(DATA1^2), rowSums(DATA2^2), '+') - tcrossprod(DATA1, 2 * DATA2)
  }

重写问题以使用线性代数和矢量化,这比循环快得多。

如果你没有lambda,这只是

outer(rowSums(DATA1^2), rowSums(DATA2^2), '+') - tcrossprod(DATA1, 2 * DATA2)

加上lambda,变成

outer(drop(DATA1^2 %*% lambda), drop(DATA2^2 %*% lambda), '+') -
    tcrossprod(DATA1, sweep(DATA2, 2, 2 * lambda, '*'))

这里使用 Rcpp 的另一种方法只是为了让这个概念文档。在一个名为 euclidean.cpp 的文件中,我有

#include <Rcpp.h>
#include <cmath>

using namespace Rcpp;

// [[Rcpp::export]]

NumericMatrix weighted_distance (NumericMatrix x, NumericMatrix y, NumericVector lambda){

  int n_x = x.nrow();
  int n_y = y.nrow();


  NumericMatrix out(n_x, n_y);

  //begin the loop

  for (int i = 0 ; i < n_x; i++){
    for (int j = 0  ; j < n_y ; j ++) {
      double d = sum(pow(x.row(i) - y.row(j), 2)*lambda);
      out(i,j) = d;
    }
  }
  return (out) ;
}

在R中,那么我有

library(Rcpp)
sourceCpp("libs/euclidean.cpp")

# Generate Data
DATA1 <- matrix(rexp(30000, rate=.1), ncol=50) #[600x50]
DATA2 <- matrix(rexp(400, rate=.1), ncol=50) #[8X50]
lambda <- runif(n=50, min=0, max=10)

# Run the program

out <- weighted_distance(DATA1, DATA2, lambda = lambda)

当我测试速度时使用:

microbenchmark(
  Rcpp_way = weighted_distance(DATA1, DATA2, lambda = lambda),
other = {DIST <- matrix(NA, nrow=nrows1, ncol=ncols)  
for (m in 1:nrows1) {
  for (i in 1:nrows2) {
    if (length(lambda) == 1) { 
      DIST[m, i] <- sum((DATA1[m, ] - DATA2[i, ])^2) 
    }
    if (length(lambda) > 1){ 
      DIST[m, i] <- sum(((DATA1[m, ] - DATA2[i, ])^2) * lambda[1, ])
    }
    next
  }
  next
}}, times = 100)

可以更快的看到是好剪辑:

Unit: microseconds
     expr       min        lq       mean    median         uq        max neval
 Rcpp_way   446.769   492.308   656.9849   562.667   846.9745   1169.231   100
    other 24688.821 30681.641 44153.5264 37511.385 50878.3585 200843.898   100