-2

我正在尝试计算

1）欧几里得距离，和

2) 马氏距离

对于 r 中的一组矩阵。我一直在这样做：

v1 <- structure(c(0.508, 0.454, 0, 2.156, 0.468, 0.488, 0.682, 1, 1.832, 
            0.44, 0.928, 0.358, 1, 1.624, 0.484, 0.516, 0.378, 1, 1.512, 
            0.514, 0.492, 0.344, 0, 1.424, 0.508, 0.56, 0.36, 1, 1.384, 0.776, 
            1.888, 0.388, 0, 1.464, 0.952, 0.252, 0.498, 1, 1.484, 0.594, 
            0.256, 0.54, 2, 2.144, 0.402, 0.656, 2.202, 1, 1.696, 0.252), 
          .Dim = c(5L, 10L), 
          .Dimnames = list(NULL, c("KW_1", "KW_2", "KW_3", "KW_4", "KW_5", "KW_6", "KW_7", "KW_8", "KW_9", "KW_10")))

v2 <- structure(c(1.864, 1.864, 1.864, 1.864, 1.864, 1.6, 1.6, 1.6, 
            1.6, 1.6, 1.536, 1.536, 1.536, 1.536, 1.536, 1.384, 1.384, 1.384, 
            1.384, 1.384, 6.368, 6.368, 6.368, 6.368, 6.368, 2.792, 2.792, 
            2.792, 2.792, 2.792, 2.352, 2.352, 2.352, 2.352, 2.352, 2.624, 
            2.624, 2.624, 2.624, 2.624, 1.256, 1.256, 1.256, 1.256, 1.256, 
            1.224, 1.224, 1.224, 1.224, 1.224), 
          .Dim = c(5L, 10L), 
          .Dimnames = list(NULL, c("KW_1", "KW_2", "KW_3", "KW_4", "KW_5", "KW_6", "KW_7", "KW_8", "KW_9", "KW_10")))

L2 <- sqrt(rowSums((v1-v2)^2))  # Euclidean distance for each row

它提供：

[1] 7.132452 7.568359 7.536904 5.448696 7.163580

那很完美！但我听说您也可以使用以下形式计算欧几里得/L2 距离：

在此处输入图像描述

我想以这种方式计算我的距离，因为马氏距离就是这个和协方差矩阵。看到这个。

但是，我还没有弄清楚如何在 r 中对此进行编码。我试过了：

sqrt(crossprod((t(v1)-t(v2))))

和

sqrt((v1-v2) %*% t(v1-v2))

但他们就是不给我我想要的。建议？

注意 - 我希望将其作为单个操作来执行，而不是在任何类型的循环中。它必须非常快，因为我要多次执行数百万行。也许这是不可能的。我愿意改变v1and的格式v2。

4

1 回答 1

1

您需要将公式分别应用于每一行，例如：

> sapply(1:nrow(v1), function(i) {
+     q = v1[i, ] - v2[i, ]
+     d = sqrt(t(q) %*% q)
+     d
+ })
[1] 7.132452 7.568359 7.536904 5.448696 7.163580

如果你需要更快的东西，你总是可以在 C++ 中尝试同样的事情（代码改编自这里）：

#include <Rcpp.h>

using namespace Rcpp;

double dist2(NumericVector x, NumericVector y){
    double d = sqrt( sum( pow(x - y, 2) ) );
    return d;
}

// [[Rcpp::export]]
NumericVector calc_l2 (NumericMatrix x, NumericMatrix y){
    int out_length = x.nrow();
    NumericVector out(out_length);

    for (int i = 0 ; i < out_length; i++){
        NumericVector v1 = x.row(i);
        NumericVector v2 = y.row(i);
        double d = dist2(v1, v2);
        out(i) = d;
    }
    return (out) ;
}

在 R 中运行：

library(Rcpp)

sourceCpp("calc_L2.cpp")
calc_l2(v1, v2)
于 2017-05-23T23:04:46.067 回答