logo标签列表

高效地匹配一个向量中的所有值到另一个向量中。

rperformanceoptimizationsearchmatching
16

我正在寻找一种高效的方法，来匹配向量x中所有的值和向量y中的值，而不仅是像match()返回的第一个位置。本质上，我要的是pmatch()的默认行为，但不带有部分匹配：

x <- c(3L, 1L, 2L, 3L, 3L, 2L)
y <- c(3L, 3L, 3L, 3L, 1L, 3L)

期望输出：

pmatch(x, y)  
[1]  1  5 NA  2  3 NA

一种方法是使用ave()，但随着组数的增加，这种方法变得缓慢且非常浪费内存：
ave(x, x, FUN = \(v) which(y == v[1])[1:length(v)])
[1]  1  5 NA  2  3 NA

有人可以推荐一个高效的方法来实现这个，最好使用（但不强制）基本R吗？

用于基准测试的更大数据集：

set.seed(5)
x <- sample(5e3, 1e5, replace = TRUE)
y <- sample(x, replace = TRUE)
5个回答
9
使用split在base中的变量。
通过其值拆分两个向量的索引。使用第一个列表的名称对第二个列表进行子集，以使它们具有相同的顺序。将NULL更改为NA，并将第二个列表的长度调整为与第一个列表相同。按照第一个列表的索引重新排序第二个列表的索引。
x <- c(3L, 1L, 2L, 3L, 3L, 2L)
y <- c(3L, 3L, 3L, 3L, 1L, 3L)

a <- split(seq_along(x), x)
b <- split(seq_along(y), y)[names(a)]
b[lengths(b)==0] <- NA
b <- unlist(Map(`length<-`, b, lengths(a)), FALSE, FALSE)
`[<-`(b, unlist(a, FALSE, FALSE), b)
#[1]  1  5 NA  2  3 NA

我试图交换这部分

b <- split(seq_along(y), y)[names(a)]
b[lengths(b)==0] <- NA

使用

b <- list2env(split(seq_along(y), y))
b <- mget(names(a), b, ifnotfound = NA)

但它并不更快。

一个RCPP版本。
将第二个向量的索引存储在每个唯一值的unordered_map中的queue中。遍历第一个向量的所有值，并从queue中获取索引。

Rcpp::sourceCpp(code=r"(
#include <Rcpp.h>
#include <unordered_map>
#include <queue>

using namespace Rcpp;
// [[Rcpp::export]]
IntegerVector pm(const std::vector<int>& a, const std::vector<int>& b) {
  IntegerVector idx(no_init(a.size()));
  std::unordered_map<int, std::queue<int> > lut;
  for(int i = 0; i < b.size(); ++i) lut[b[i]].push(i);
  for(int i = 0; i < idx.size(); ++i) {
    auto search = lut.find(a[i]);
    if(search != lut.end() && search->second.size() > 0) {
      idx[i] = search->second.front() + 1;
      search->second.pop();
    } else {idx[i] = NA_INTEGER;}
  }
  return idx;
}
)")
pm(x, y)
#[1]  1  5 NA  2  3 NA

这是一个专门针对此用例的RCPP版本。
创建一个长度等于第一个向量最大值的向量，并计算每个值出现的次数。创建另一个同样长度的queue向量，并存储第二个向量的值的索引，直到达到第一个向量的数量为止。遍历第一个向量的所有值，并从queue中获取索引。

Rcpp::sourceCpp(code=r"(
#include <Rcpp.h>
#include <vector>
#include <array>
#include <queue>
#include <algorithm>

using namespace Rcpp;
// [[Rcpp::export]]
IntegerVector pm2(const std::vector<int>& a, const std::vector<int>& b) {
  IntegerVector idx(no_init(a.size()));
  int max = 1 + *std::max_element(a.begin(), a.end());
  std::vector<int> n(max);
  for(int i = 0; i < a.size(); ++i) ++n[a[i]];
  std::vector<std::queue<int> > lut(max);
  for(int i = 0; i < b.size(); ++i) {
    if(b[i] < max && n[b[i]] > 0) {
      --n[b[i]];
      lut[b[i]].push(i);
    }
  }
  for(int i = 0; i < idx.size(); ++i) {
    auto & P = lut[a[i]];
    if(P.size() > 0) {
      idx[i] = P.front() + 1;
      P.pop();
    } else {idx[i] = NA_INTEGER;}
  }
  return idx;
}
)")
pm2(x,y)
#[1]  1  5 NA  2  3 NA

基准测试

set.seed(5)
x <- sample(5e3, 1e5, replace = TRUE)
y <- sample(x, replace = TRUE)

library(data.table)

matchall <- function(x, y) {
  data.table(y, rowid(y))[
    data.table(x, rowid(x)), on = .(y = x, V2), which = TRUE
  ]
}

rmatch <- function(x, y) {
  xp <- cbind(seq_along(x), x)[order(x),]
  yp <- cbind(seq_along(y), y)[order(y),]
  result <- numeric(length(x))
  
  xi <- yi <- 1
  Nx <- length(x)
  Ny <- length(y)
  while (xi <= Nx) {
    if (yi > Ny) {
      result[xp[xi,1]] <- NA
      xi <- xi + 1
    } else if (xp[xi,2] == yp[yi,2]) {
      result[xp[xi,1]] = yp[yi,1]
      xi <- xi + 1
      yi <- yi + 1
    } else if (xp[xi,2] < yp[yi,2]) {
      result[xp[xi,1]] <- NA
      xi <- xi + 1
    } else if (xp[xi,2] > yp[yi,2]) {
      yi <- yi + 1
    }
  }
  result  
}

bench::mark(
ave = ave(x, x, FUN = \(v) which(y == v[1])[1:length(v)]),
rmatch = rmatch(x, y),
make.name = match(make.names(x, TRUE), make.names(y, TRUE)),
paste = do.call(match, lapply(list(x, y), \(v) paste(v, ave(v, v, FUN = seq_along)))),
make.unique = match(make.unique(as.character(x)), make.unique(as.character(y))),
split = {a <- split(seq_along(x), x)
  b <- split(seq_along(y), y)[names(a)]
  b[lengths(b)==0] <- NA
  b <- unlist(Map(`length<-`, b, lengths(a)), FALSE, FALSE)
  `[<-`(b, unlist(a, FALSE, FALSE), b)},
data.table = matchall(x, y),
RCPP = pm(x, y),
RCPP2 = pm2(x, y)
)

结果

  expression       min   median `itr/sec` mem_alloc `gc/sec` n_itr  n_gc
  <bch:expr>  <bch:tm> <bch:tm>     <dbl> <bch:byt>    <dbl> <int> <dbl>
1 ave            1.66s    1.66s     0.603    3.73GB    68.7      1   114
2 rmatch      258.29ms 259.35ms     3.86     5.34MB    30.8      2    16
3 make.name   155.69ms 156.82ms     6.37    14.06MB     1.59     4     1
4 paste         93.8ms 102.06ms     9.74    18.13MB     7.79     5     4
5 make.unique  81.67ms   92.8ms    10.4      9.49MB     5.22     6     3
6 split        12.66ms  13.16ms    65.8      7.18MB    16.0     33     8
7 data.table    6.22ms   6.89ms   114.       5.13MB    28.0     57    14
8 RCPP          3.06ms    3.2ms   301.     393.16KB     3.98   151     2
9 RCPP2         1.64ms   1.82ms   514.     393.16KB     8.00   257     4

在这种情况下，C++ 版本是最快的，并且分配的内存最少。如果使用 base，那么 splitB 变体是最快的，rmatch 分配的内存最少。
谢谢@GKi。我很乐意接受任何提供的答案，但你的split选项是最有效的基本解决方案。 - Ritchie Sacramento
通过一些更改，它可以比我的第一篇文章快得多。也许还有其他在base中更快的方法。 - GKi
7

只是想指出，你可以使用match + make.unique来实现相同的功能。从速度上来说，它可能比data.table方法慢：

match(make.unique(as.character(x)), make.unique(as.character(y)))

[1]  1  5 NA  2  3 NA

match(make.names(x, TRUE), make.names(y, TRUE))
[1]  1  5 NA  2  3 NA
6

使用 data.table 连接，受到 this 问答的启发。

library(data.table)

matchall <- function(x, y) {
  data.table(y, rowid(y))[
    data.table(x, rowid(x)), on = .(y = x, V2), which = TRUE
  ]
}

检查行为

x <- c(3L, 1L, 2L, 3L, 3L, 2L)
y <- c(3L, 3L, 3L, 3L, 1L, 3L)

matchall(x, y)
#> [1]  1  5 NA  2  3 NA

处理更大的向量的时间：

set.seed(5)
x <- sample(5e3, 1e5, replace = TRUE)
y <- sample(x, replace = TRUE)

system.time(z1 <- matchall(x, y))
#>    user  system elapsed 
#>    0.06    0.00    0.01

system.time(z2 <- ave(x, x, FUN = \(v) which(y == v[1])[1:length(v)]))
#>    user  system elapsed 
#>    0.88    0.43    1.31

identical(z1, z2)
#> [1] TRUE
4
如果您有一些额外的内存可以使用，您可以通过对值进行排序并基本上执行两个指针遍历来加快进程以匹配数据。以下是示例：
rmatch <- function(x, y) {
  xp <- cbind(seq_along(x), x)[order(x),]
  yp <- cbind(seq_along(y), y)[order(y),]
  result <- numeric(length(x))
  
  xi <- yi <- 1
  Nx <- length(x)
  Ny <- length(y)
  while (xi <= Nx) {
    if (yi > Ny) {
      result[xp[xi,1]] <- NA
      xi <- xi + 1
    } else if (xp[xi,2] == yp[yi,2]) {
      result[xp[xi,1]] = yp[yi,1]
      xi <- xi + 1
      yi <- yi + 1
    } else if (xp[xi,2] < yp[yi,2]) {
      result[xp[xi,1]] <- NA
      xi <- xi + 1
    } else if (xp[xi,2] > yp[yi,2]) {
      yi <- yi + 1
    }
  }
  result  
}

我用这里发布的其他一些基本R选项进行了测试

mbm <- microbenchmark::microbenchmark(
  ave = ave(x, x, FUN = \(v) which(y == v[1])[1:length(v)]),
  rmatch = rmatch(x, y),
  pmatch = pmatch(x, y),
  times = 20
)

并且发现它似乎表现良好

Unit: milliseconds
   expr        min         lq       mean     median         uq        max neval
    ave  1227.6743  1247.6980  1283.1024  1264.1485  1324.1569  1349.3276    20
 rmatch   198.1744   201.1058   208.3158   204.5933   209.4863   247.7279    20
 pmatch 39514.4227 39595.9720 39717.5887 39628.0892 39805.2405 40105.4337    20

它们都返回相同的值向量。

1
使用C++实现，也许可以与data.table具有竞争力？ - s_baldur
2
当然可以。但我的目标是只使用基本的R，没有依赖项（包括编译C++所需的系统工具）。已经采用data.table方法，因为大部分工作都在C++后端完成，所以速度更快。 - MrFlick
1
@RitchieSacramento 谢谢你的测试用例。我发现了一个差一的错误并进行了修复。但我同意 split() 方法是更好的选择。 - MrFlick
2

您可以简单地运行match + paste + ave

> do.call(match, lapply(list(x, y), \(v) paste(v, ave(v, v, FUN = seq_along))))
[1]  1  5 NA  2  3 NA
网页内容由stack overflow 提供, 点击上面的
可以查看英文原文，
原文链接