Comparing many probability density functions

dist_one_one()

dist_one_one() is a lower level equivalent to distance(). However, instead of accepting a numeric data.frame or matrix, it expects two vectors representing probability density functions. In this example, we create two vectors, P and Q.

P <- 1:10 / sum(1:10)
Q <- 20:29 / sum(20:29)

To calculate the euclidean distance between them we can use several approaches - (a) build-in R dist() function, (b) philentropy::distance(), (c) philentropy::euclidean(), or the new dist_one_one().

# install.packages("microbenchmark")
microbenchmark::microbenchmark(
  dist(rbind(P, Q), method = "euclidean"),
  distance(rbind(P, Q), method = "euclidean", test.na = FALSE, mute.message = TRUE),
  euclidean(P, Q, FALSE),
  dist_one_one(P, Q, method = "euclidean", testNA = FALSE)
)

## Warning in microbenchmark::microbenchmark(dist(rbind(P, Q), method =
## "euclidean"), : less accurate nanosecond times to avoid potential integer
## overflows

## Unit: nanoseconds
##                                                                                    expr
##                                                 dist(rbind(P, Q), method = "euclidean")
##  distance(rbind(P, Q), method = "euclidean", test.na = FALSE,      mute.message = TRUE)
##                                                                  euclidean(P, Q, FALSE)
##                                dist_one_one(P, Q, method = "euclidean", testNA = FALSE)
##   min     lq     mean median     uq    max neval
##  5248 5514.5  6330.40   5740 5904.0  62730   100
##  7626 8036.0 12712.05   8200 8425.5 450713   100
##   492  533.0   954.89    574  656.0  37720   100
##   779  881.5  1258.70    943 1025.0  32472   100

All of them return the same, single value. However, as you can see in the benchmark above, some are more flexible, and others are faster.

dist_one_many()

The role of dist_one_many() is to calculate distances between one probability density function (in a form of a vector) and a set of probability density functions (as rows in a matrix).

Firstly, let’s create our example data.

set.seed(2020-08-20)
P <- 1:10 / sum(1:10)
M <- t(replicate(100, sample(1:10, size = 10) / 55))

P is our input vector and M is our input matrix.

Distances between the P vector and probability density functions in M can be calculated using several approaches. For example, we could write a for loop (adding a new code) or just use the existing distance() function and extract only one row (or column) from the results. The dist_one_many() allows for this calculation directly as it goes through each row in M and calculates a given distance measure between P and values in this row.

# install.packages("microbenchmark")
microbenchmark::microbenchmark(
  as.matrix(dist(rbind(P, M), method = "euclidean"))[1, ][-1],
  distance(rbind(P, M), method = "euclidean", test.na = FALSE, mute.message = TRUE)[1, ][-1],
  dist_one_many(P, M, method = "euclidean", testNA = FALSE)
)

## Unit: microseconds
##                                                                                             expr
##                                      as.matrix(dist(rbind(P, M), method = "euclidean"))[1, ][-1]
##  distance(rbind(P, M), method = "euclidean", test.na = FALSE,      mute.message = TRUE)[1, ][-1]
##                                        dist_one_many(P, M, method = "euclidean", testNA = FALSE)
##       min       lq       mean    median        uq      max neval
##    64.780   67.896   73.16450   70.9095   74.9685  150.142   100
##  6200.061 6293.172 6468.79058 6349.9775 6662.4180 8237.433   100
##     7.667    8.446    9.19343    8.7535    9.1225   38.007   100

The dist_one_many() returns a vector of values. It is, in this case, much faster than distance(), and visibly faster than dist() while allowing for more possible distance measures to be used.

dist_many_many()

dist_many_many() calculates distances between two sets of probability density functions (as rows in two matrix objects).

Let’s create two new matrix example data.

set.seed(2020-08-20)
M1 <- t(replicate(10, sample(1:10, size = 10) / 55))
M2 <- t(replicate(10, sample(1:10, size = 10) / 55))

M1 is our first input matrix and M2 is our second input matrix. I am not aware of any function build-in R that allows calculating distances between rows of two matrices, and thus, to solve this problem, we can create our own - many_dists()…

many_dists = function(m1, m2){
  r = matrix(nrow = nrow(m1), ncol = nrow(m2))
  for (i in seq_len(nrow(m1))){
    for (j in seq_len(nrow(m2))){
      x = rbind(m1[i, ], m2[j, ])
      r[i, j] = distance(x, method = "euclidean", mute.message = TRUE)
    }
  }
  r
}

… and compare it to dist_many_many().

# install.packages("microbenchmark")
microbenchmark::microbenchmark(
  many_dists(M1, M2),
  dist_many_many(M1, M2, method = "euclidean", testNA = FALSE)
)

## Unit: microseconds
##                                                          expr     min      lq
##                                            many_dists(M1, M2) 795.810 874.284
##  dist_many_many(M1, M2, method = "euclidean", testNA = FALSE)  12.505  13.940
##        mean   median       uq      max neval
##  1034.74775 934.6155 1000.216 3372.045   100
##    16.12899  15.4160   17.507   32.718   100

Both many_dists()and dist_many_many() return a matrix. The above benchmark concludes that dist_many_many() is about 30 times faster than our custom many_dists() approach.