计算一组间隔之间的成对距离

Question

假设我有一组由此矩阵表示的闭合线性区间:

interval.mat = matrix(c(1,2,3,5,4,6,8,9), byrow = TRUE, ncol = 2)

interval.mat[,1]区间起点在哪里,interval.mat[,2]是它们对应的终点.

我正在寻找一种有效的(因为这个例子矩阵是一个玩具,实际上我的矩阵包含几千个间隔)的方式来产生一个矩阵,它将保持间隔之间的所有成对正距离.一对间隔之间的距离应该是间隔的开始,两者之间的较大端减去间隔的末端,两者之间的末端较小.例如间隔之间的距离c(1,2)和c(3,5)应3 - 2 = 1,由于第二间隔中的第一一个之后结束.在情况下,时间间隔重叠的距离应为0.因此,例如,在的情况下,c(3,5)和c(4,6)该距离将是0.

因此,上述间隔的成对距离矩阵将是:

> matrix(c(0,1,2,6,1,0,0,3,2,0,0,2,6,3,2,0), byrow = TRUE, nrow = 4, ncol = 4)
     [,1] [,2] [,3] [,4]
[1,]    0    1    2    6
[2,]    1    0    0    3
[3,]    2    0    0    2
[4,]    6    3    2    0

Answer 1

这是一个Rcpp解决方案.它将具有快速和内存效率(详见下文).

首先让我们定义一个辅助函数来计算所有成对距离.如果n是要考虑的间隔数,我们有n*(n-1)/2唯一的向量对(当然,我们不会考虑相同的间隔,因为它们之间的距离是0).

library('Rcpp')
library('inline')    
cppFunction("
   NumericVector distint_help(NumericMatrix x) {
      int n = x.nrow(); // number of rows
      NumericVector out(n*(n-1)/2); // result numeric vector
      int k = 0;
      for (int i=0; i<n-1; ++i) {
         for (int j=i+1; j<n; ++j) {
            if (x(i,0) >= x(j,1))
               out[k++] = x(i,0)-x(j,1);
            else if (x(j,0) > x(i,1))
               out[k++] = x(j,0)-x(i,1);
            else
               out[k++] = 0.0;
         }
      }
      return out;
   }
")

上面的函数返回一个带有计算距离的数值向量.让我们尝试模仿内置dist函数的输出(检查结果x <- dist(interval.mat); unclass(x)).

现在主要功能:

distint <- function(interval) {
   stopifnot(is.numeric(interval), is.matrix(interval), ncol(interval) == 2)

   res <- distint_help(interval) # use Rcpp to calculate the distances

   # return the result similar to the one of dist()
   structure(res, class='dist', Size=nrow(interval), Diag=FALSE, Upper=FALSE)
}
distint(interval.mat)
##   1 2 3
## 2 1    
## 3 2 0  
## 4 6 3 2

上述对象可以转换为"普通"方阵:

as.matrix(distint(interval.mat))
##   1 2 3 4
## 1 0 1 2 6
## 2 1 0 0 3
## 3 2 0 0 2
## 4 6 3 2 0

除非距离矩阵稀疏(存在许多零),否则上述解决方案是存储有效的.

基准:

test <- matrix(runif(1000), ncol=2)
library('microbenchmark')
library(proxy)
f <- function(x,y) max(min(x)-max(y),0)
microbenchmark(distint(test), as.matrix(dist(test, method=f)), times=10)

## Unit: milliseconds
##                               expr        min         lq     median         uq        max neval
##                      distint(test)   1.584548   1.615146   1.650645   3.071433   3.164231    10
##  as.matrix(dist(test, method = f)) 455.300974 546.438875 551.596582 599.977164 609.418194    10