两个逻辑向量中元素的最小距离（区间）快速计算方法（第二部分）

Question

我曾在这里提出一个相关问题，但意识到计算此复杂度量需要耗费太多时间（而目标是与随机化测试一起使用，因此速度是一个问题）。 因此，我决定放弃加权，并仅使用两个度量之间的最小距离。 所以在这里，我有两个向量（在数据框中进行演示，但实际上它们是两个向量）。

       x     y
1  FALSE  TRUE
2  FALSE FALSE
3   TRUE FALSE
4  FALSE FALSE
5  FALSE  TRUE
6  FALSE FALSE
7  FALSE FALSE
8   TRUE FALSE
9  FALSE  TRUE
10  TRUE  TRUE
11 FALSE FALSE
12 FALSE FALSE
13 FALSE FALSE
14 FALSE  TRUE
15  TRUE FALSE
16 FALSE FALSE
17  TRUE  TRUE
18 FALSE  TRUE
19 FALSE FALSE
20 FALSE  TRUE
21 FALSE FALSE
22 FALSE FALSE
23 FALSE FALSE
24 FALSE FALSE
25  TRUE FALSE

我已经编写了一些代码来找到最小距离，但是我需要更快的速度（减少不必要的调用和更好的向量化）。也许在基本的R中我无法更快。

## MWE EXAMPLE: THE DATA
x <- y <- rep(FALSE, 25)
x[c(3, 8, 10, 15, 17, 25)] <- TRUE
y[c(1, 5, 9, 10, 14, 17, 18, 20)] <- TRUE

## Code to Find Distances
xw <- which(x)
yw <- which(y)

min_dist <- function(xw, yw) {
    unlist(lapply(xw, function(x) {
        min(abs(x - yw))
    }))
}

min_dist(xw, yw)

有没有办法在基本的R中提高性能？使用dplyr或data.table吗？

我的向量要长得多（10,000+元素）。

根据flodel的测试结果进行编辑。flodel，我预料到了MWE中可能会出现的问题，但我也不知道该如何解决。如果任何x位置小于最小y位置，则会出现问题。

x <- y <- rep(FALSE, 25)
x[c(3, 8, 9, 15, 17, 25)] <- TRUE
y[c(5, 9, 10, 13, 15, 17, 19)] <- TRUE


xw <- which(x)
yw <- which(y)

flodel <- function(xw, yw) {
   i <- findInterval(xw, yw)
   pmin(xw - yw[i], yw[i+1L] - xw, na.rm = TRUE)
}

flodel(xw, yw)

## [1] -2 -1 -6 -2 -2 20
## Warning message:
## In xw - yw[i] :
##   longer object length is not a multiple of shorter object length

- Tyler Rinker

Answer 1

flodel <- function(x, y) {
  xw <- which(x)
  yw <- which(y)
  i <- findInterval(xw, yw, all.inside = TRUE)
  pmin(abs(xw - yw[i]), abs(xw - yw[i+1L]), na.rm = TRUE)
}

GG1 <- function(x, y) {
  require(zoo)
  yy <- ifelse(y, TRUE, NA) * seq_along(y)
  fwd <- na.locf(yy, fromLast = FALSE)[x]
  bck <- na.locf(yy, fromLast = TRUE)[x]
  wx <- which(x)
  pmin(wx - fwd, bck - wx, na.rm = TRUE)
}

GG2 <- function(x, y) {
  require(data.table)
  dtx <- data.table(x = which(x))
  dty <- data.table(y = which(y), key = "y")
  dty[dtx, abs(x - y), roll = "nearest"] 
}

示例数据：

x <- y <- rep(FALSE, 25)
x[c(3, 8, 10, 15, 17, 25)] <- TRUE
y[c(1, 5, 9, 10, 14, 17, 18, 20)] <- TRUE

X <- rep(x, 100)
Y <- rep(y, 100)

单元测试：

identical(flodel(X, Y), GG1(X, Y))
# [1] TRUE

基准测试：

library(microbenchmark)
microbenchmark(flodel(X,Y), GG1(X,Y), GG2(X,Y))
# Unit: microseconds
#          expr       min         lq     median        uq        max neval
#  flodel(X, Y)   115.546   131.8085   168.2705   189.069   1980.316   100
#     GG1(X, Y)  2568.045  2828.4155  3009.2920  3376.742  63870.137   100
#     GG2(X, Y) 22210.708 22977.7340 24695.7225 28249.410 172074.881   100

[Matt Dowle编辑] 24695微秒=0.024秒。在使用微小数据进行的基准测试中得出的推论很少适用于有意义的数据大小。

[flodel编辑] 我的向量长度为2500，这相当有意义，考虑到Tyler的说法（10k），但好吧，让我们尝试一下长度为2.5e7的向量。在这种情况下，请原谅我使用system.time：

X <- rep(x, 1e6)
Y <- rep(y, 1e6)
system.time(flodel(X,Y))
#    user  system elapsed 
#   0.694   0.205   0.899 
system.time(GG1(X,Y))
#    user  system elapsed 
#  31.250  16.496 112.967 
system.time(GG2(X,Y))
# Error in `[.data.table`(dty, dtx, abs(x - y), roll = "nearest") : 
#   negative length vectors are not allowed

---

[Arun编辑] - 使用1.8.11版本进行2.5e7基准测试：
[Arun第二次编辑] - 在Matt最近更快的二分查找/合并后更新时间

require(data.table)
arun <- function(x, y) {
    dtx <- data.table(x=which(x))
    setattr(dtx, 'sorted', 'x')
    dty <- data.table(y=which(y))
    setattr(dty, 'sorted', 'y')
    dty[, y1 := y]
    dty[dtx, roll="nearest"][, abs(y-y1)]
}

# minimum of three consecutive runs
system.time(ans1 <- arun(X,Y))
#   user  system elapsed 
#  1.036   0.138   1.192 

# minimum of three consecutive runs
system.time(ans2 <- flodel(X,Y))
#   user  system elapsed 
#  0.983   0.197   1.221 

identical(ans1, ans2) # [1] TRUE

- flodel

---

有一个问题是MWE没有考虑到。我在我的解决方案中加入了编辑，但你的方法更快，我不想放弃它。 - Tyler Rinker

我认为我修好了，试一下吧。 - flodel

是的，确实。谢谢。我现在也明白了all.inside = TRUE的作用。感谢你向我介绍了一个新的函数。 - Tyler Rinker

@Tyler Rinker，你能发布一下10,000个数据集的相应时间吗？ - G. Grothendieck

@G.Grothendieck 是的，但我需要几天时间才能处理并测试它。 - Tyler Rinker

@G.Grothendieck 我一直在尝试对此进行基准测试，但是一旦数据的规模增大，解决方案就不再提供相同的价值。 - Ricardo Saporta

Answer 2

这里有两个解决方案，都不使用循环或应用函数。

1) 如果z为1，则第一个解决方案与我在您的prior question中发布的解决方案相同，但此处的简化假设使我们可以缩短它，并且相对于那个答案，我们已将答案减少了1。

library(zoo)

yy <- ifelse(y, TRUE, NA) * seq_along(y)
fwd <- na.locf(yy, fromLast = FALSE)[x]
bck <- na.locf(yy, fromLast = TRUE)[x]
wx <- which(x)
pmin(wx - fwd, bck - wx, na.rm = TRUE)

2) 第二种方法是使用 data.table。data.table 可以使用 roll="nearest" 参数，看起来正是你需要的：

library(data.table)

dtx <- data.table(x = which(x))
dty <- data.table(y = which(y), key = "y")
dty[dtx, abs(x - y), roll = "nearest"] 

我不确定这是否重要，但我正在使用data.table 1.8.11版本（目前CRAN版本是1.8.10）。

- G. Grothendieck