data.table join 中的 Rcpp 函数

Question

问题

为什么我的 data.table 连接内的 Rcpp 函数与在连接外部使用时产生不同(且不正确)的结果？

示例

我有两个 data.table，我想找到两个表中每对坐标之间的欧几里得距离。

为了进行距离计算，我定义了两个函数，一个以 R 为基础，另一个使用 Rcpp。

library(Rcpp)
library(data.table)

rEucDist <- function(x1, y1, x2, y2) return(sqrt((x2 - x1)^2 + (y2 - y1)^2))

cppFunction('NumericVector cppEucDistance(NumericVector x1, NumericVector y1,
                                          NumericVector x2, NumericVector y2){

            int n = x1.size();
            NumericVector distance(n);

            for(int i = 0; i < n; i++){
              distance[i] = sqrt(pow((x2[i] - x1[i]), 2) + pow((y2[i] - y1[i]), 2));
            }
            return distance;
}')


dt1 <- data.table(id = rep(1, 6),
                  seq1 = 1:6,
                  x = c(1:6),
                  y = c(1:6))

dt2 <- data.table(id = rep(1, 6),
                  seq2 = 7:12,
                  x = c(6:1),
                  y = c(6:1))

当先进行连接，然后计算距离时，两个函数产生相同的结果

dt_cpp <- dt1[ dt2, on = "id", allow.cartesian = T]
dt_cpp[, dist := cppEucDistance(x, y, i.x, i.y)]

dt_r <- dt1[ dt2, on = "id", allow.cartesian = T]
dt_r[, dist := rEucDist(x, y, i.x, i.y)]

all.equal(dt_cpp$dist, dt_r$dist)
# [1] TRUE

但是，如果我在连接内进行计算，结果会有所不同； cpp版本不正确。

dt_cppJoin <- dt1[
    dt2,
    { (cppEucDistance(x, y, i.x, i.y)) },
    on = "id",
    by = .EACHI
    ]

dt_rJoin <- dt1[
    dt2,
    { (rEucDist(x, y, i.x, i.y)) },
    on = "id",
    by = .EACHI
    ]

all.equal(dt_cppJoin$V1, dt_rJoin$V1)
# "Mean relative difference: 0.6173913"

## note that the R version of the join is correct
all.equal(dt_r$dist, dt_rJoin$V1)
# [1] TRUE

Rcpp 实现导致连接版本给出不同结果的原因是什么？

Answer 1

tl;博士

1. 最终，这归因于我对 EACHI 的工作原理缺乏了解。

2. 我的 C++ 循环目前依赖于 x1 和 x2 的长度相同。使用 EACHI 时情况并非如此，因此我在 C++ 函数中的向量取子集将不正确。因此我需要一个更复杂的 C++ 函数。

<小时/>

我相信我看到的“问题”归结为 EACHI 正在做什么。一遍又一遍地重新阅读 Arun's answer 后，我想我明白了。

特别是，EACHI“评估 Y 中每一行的匹配行上的 j 表达式”。因此，它采用 Y 表(在我的例子中为 dt2)并一次对每一行进行评估。如果您在 Rcpp 函数内使用一些prints，就可以看到这一点

## I'm reducing the data sets so the 'prints' are readible
dt1 <- data.table(id = rep(1, 2),
                  seq1 = 1:2,
                  x = c(1:2),
                  y = c(1:2))

dt2 <- data.table(id = rep(1, 2),
                  seq2 = 7:8,
                  x = c(6:5),
                  y = c(6:5))

cppFunction('NumericVector cppEucDistance(NumericVector x1, NumericVector y1,
                                          NumericVector x2, NumericVector y2){

  int n = x1.size();
  NumericVector distance(n);

  Rcpp::Rcout << "x1 size: " << x1.size() << std::endl;
  Rcpp::Rcout << "x2 size: " << x2.size() << std::endl;
  Rcpp::Rcout << "n size: " << n << std::endl;

  for(int i = 0; i < n; i++){
    distance[i] = sqrt(pow((x2[i] - x1[i]), 2) + pow((y2[i] - y1[i]), 2));
  }
  return distance;
}')

dt_cppJoin <- dt1[
  dt2,
  {print(i.x); (cppEucDistance(x, y, i.x, i.y)) },
  on = "id",
  by = .EACHI
  ]
# [1] 6
# x1 size: 2
# x2 size: 1
# n size: 2
# [1] 5
# x1 size: 2
# x2 size: 1
# n size: 2

请注意，在打印语句的输出中，x2 每次迭代的长度仅为 1。而 n 是 2。因此，当 i 到达第二个元素时，我的 x2[i] 显然会“失败”(请记住，C++ 数组是0-索引)

而删除 EACHI 与完全连接方法具有相同的效果

dt_cppJoin <- dt1[
  dt2,
  {print(i.x); (cppEucDistance(x, y, i.x, i.y)) },
  on = "id"
  ]
# [1] 6 6 5 5
# x1 size: 4
# x2 size: 4
# n size: 4

<小时/>