许多矩阵对之间的相似度/距离

Question

我希望通过计算每对（多维）点集之间的距离平均值来量化群体相似度。

我可以像这样轻松地为每个组合手动完成这项任务：

library(dplyr)
library(tibble)
library(proxy)

# dummy data
set.seed(123)
df1 <- data.frame(x = rnorm(100,0,4), 
                  y = rnorm(100,1,5), 
                  z = rbinom(100, 1, 0.1))
df2 <- data.frame(x = rnorm(100,-1,3), 
                  y = rnorm(100,0,6), 
                  z = rbinom(100, 1, 0.1))
df3 <- data.frame(x = rnorm(100,-30,4), 
                  y = rnorm(100,10,2), 
                  z = rbinom(100, 1, 0.9))

# compute distance (unscaled, uncentred data)
dist(df1, df2, method = "gower") %>% mean
dist(df1, df3, method = "gower") %>% mean
dist(df2, df3, method = "gower") %>% mean

但我希望以某种方式将其向量化，因为我的实际数据有30多个组。可以使用简单的for循环来实现此目的，如下所示：

# combine data and scale, centre
df <- rbind(df1, df2, df3) %>% 
  mutate(id = rep(1:3, each = 100))
df <- df %>% 
  select(-id) %>%
  transmute_all(scale) %>% 
  add_column(id = df$id)

# create empty matrix for comparisons
n <- df$id %>% unique %>% length
m <- matrix(nrow = n, ncol = n)

# loop through each pair once
for(i in 1:n) {
  for(j in 1:i) { #omit top right corner
    if(i == j) {
      m[i,j] <- NA #omit diagonal
    } else {
      m[i,j] <- dist(df[df$id == i,1:3], df[df$id == j,1:3], method = "gower") %>% mean
    }
  }
}

m
          [,1]      [,2] [,3]
[1,]        NA        NA   NA
[2,] 0.2217443        NA   NA
[3,] 0.8446070 0.8233932   NA

然而，这种方法的可伸缩性非常差；快速基准测试表明，使用我的实际数据（每组1000多行，共30多组），需要花费90小时以上。

有没有人能够建议一种更有效的解决方案，或者可能是我忽略的根本不同的问题框架？

- jogall

Answer 1

我不确定这种方法是否有效，但这是另一种方法。您可以使用ls获取矩阵的名称，使用combn生成两个矩阵对，然后使用get获取计算dist所需的矩阵。"最初的回答"

do.call(rbind,
        combn(ls(pattern = "df\\d+"), 2, FUN = function(x)
            data.frame(pair = toString(x),
                       dist = mean(dist(get(x[1]), get(x[2]), method = "gower")),
                       stringsAsFactors = FALSE),
            simplify = FALSE
        ))
#      pair      dist
#1 df1, df2 0.2139304
#2 df1, df3 0.8315169
#3 df2, df3 0.8320911

- d.b

---

感谢@d.b的回复。不幸的是，这并没有快多少；使用上面的示例，组大小为500，在我的机器上需要32.16秒，而for循环需要32.84秒。 - jogall

@jogall，看起来最耗时的步骤是 dist。除非有其他包含更快实现的 dist 的包，否则可能无法做太多来提高速度。 - d.b

是的，看起来是这样--我猜我希望有人能够建议一种在概念上不同的方法来解决问题，例如可以采用某种有效的快捷方式来计算组平均值，而不是计算集合中每对向量之间的距离。 - jogall

@jogall，请查看gower::gower_dist()。它似乎更快，但会给出不同的值。 - d.b

Answer 2

您可以将每对组合并，然后仅在该组内计算不相似度矩阵。显然，这意味着您在某种程度上正在将一组与自身进行比较，但它可能仍适用于您的用例，并且使用可以在数据规模较大时快速完成。

library(cluster)

n <- 30
groups <- vector("list", 30)

# dummy data
set.seed(123)
for(i in 1:30) {
  groups[[i]] = data.frame(x = rnorm(1000,ceiling(runif(1, -10, 10)),ceiling(runif(1, 2, 4))), 
                           y = rnorm(1000,ceiling(runif(1, -10, 10)),ceiling(runif(1, 2, 4))), 
                           z = rbinom(1000,1,runif(1,0.1,0.9)))
}

m <- matrix(nrow = n, ncol = n)

# loop through each pair once
for(i in 1:n) {
  for(j in 1:i) { #omit top right corner
    if(i == j) {
      m[i,j] <- NA #omit diagonal
    } else {
      # concatenate groups
      dat <- rbind(df_list[[i]], df_list[[j]])

      # compute all distances (between groups and within groups), return matrix
      mm <- dat %>% 
        daisy(metric = "gower") %>%
        as.matrix

      # retain only distances between groups
      mm <- mm[(nrow(df_list[[i]])+1):nrow(dat) , 1:nrow(df_list[[i]])]

      # write mean distance to global comparison matrix
      m[i,j] <- mean(mm)
    }
  }
}

- Thom

---

谢谢@Thom，daisy对Gower的实现与proxy相比速度惊人地快；对于1000行的3个组，您的方法需要2.4秒，而我的方法需要146.5秒，即使它也计算了组内距离。我只会在您的代码中添加一点内容，展示如何在计算平均值之前仅保留组间距离，因为这些是我在这里感兴趣的所有内容。 - jogall

Answer 3

proxy 可以处理矩阵列表作为输入，你只需要定义一个包装函数来实现你想要的功能：

nested_gower <- function(x, y, ...) {
  mean(proxy::dist(x, y, ..., method = "gower"))
}

proxy::pr_DB$set_entry(
  FUN = nested_gower,
  names = c("ngower"),
  distance = TRUE,
  loop = TRUE
)

df_list <- list(df1, df2, df3)
proxy::dist(df_list, df_list, method = "ngower")
     [,1]      [,2]      [,3]     
[1,] 0.1978306 0.2139304 0.8315169
[2,] 0.2139304 0.2245903 0.8320911
[3,] 0.8315169 0.8320911 0.2139049

这仍然会很慢，但它应该比纯R中的for循环更快（proxy在后台使用C语言）。

重要提示：请注意，所得交叉距离矩阵的对角线不为零。 如果您像proxy::dist(df_list, method = "ngower")那样调用dist, proxy将假定distance(x, y) = distance(y, x)（对称性）， 以及distance(x, x) = 0，而后者在这种情况下并不成立。 向dist传递两个参数可以避免这种假设。 如果您真的不关心对角线， 只传递一个参数可以节省一些额外的时间，避免计算上三角。 或者，如果您关心对角线但仍想避免计算上三角， 先使用一个参数调用dist，然后再调用proxy::dist(df_list, df_list, method = "ngower", pairwise = TRUE)。

附注：如果您想模仿d.b建议的gower包的行为， 您可以将包装函数定义为：

nested_gower <- function(x, y, ...) {
  distmat <- sapply(seq_len(nrow(y)), function(y_row) {
      gower::gower_dist(x, y[y_row, , drop = FALSE], ...)
  })

  mean(distmat)
}

然而，返回的值似乎会根据传递给函数的记录数量而发生变化， 因此很难确定最佳方法是什么。

*如果您想重新定义proxy中的函数，请先使用proxy::pr_DB$delete_entry("ngower")。

---

如果您更喜欢使用代理的Gower交叉距离矩阵版本，那么您可以利用我的dtwclust包的一些功能并行计算：

library(dtwclust)
library(doParallel)

custom_dist <- new("tsclustFamily", dist = "ngower", control = list(symmetric = TRUE))@dist

workers <- makeCluster(detectCores())
registerDoParallel(workers)

distmat <- custom_dist(df_list)

stopCluster(workers); registerDoSEQ()

对于实际使用情况来说，这可能会更快（但对于这里的小样本数据而言并不是太明显）。同样需要注意对角线（因此请使用custom_dist(df_list, df_list) 或 custom_dist(df_list, pairwise = TRUE)）。如果您想了解更多信息，请参见第3.2节here和tsclustFamily文档。

- Alexis

---

感谢您的详细回答，@Alexis。我不知道可以使用proxy来定义包装函数，这真的很有用。然而，我接受了@Thom的答案（稍作修改），仅仅因为它更快，总计算时间是我面临的最大障碍--对于每个有1000行的3组样本，该方法需要2.4秒，而此方法需要89.4秒。 - jogall

@jogall，你仍然可以在 nested_gower 中嵌入建议的 daisy 过程，并让 proxy 进行循环。 - Alexis