如何高效计算超大数据集列对间的距离相关系数？

你的核心瓶颈在于重复调用na.omit——每次计算列对时都要重新检查NA、创建新矩阵，在2600列的场景下（总共有约340万对列），这种重复操作会产生巨大的性能开销。以下是针对性的优化方案：

方案1：预存非NA索引，直接提取有效行（最有效） #

提前记录每列的非NA行位置，两两计算时直接取两个列的非NA行交集，跳过na.omit的额外开销（它会做数据框转换、属性处理等冗余操作），同时只计算上三角并对称填充，减少一半计算量。

library(energy)

# 模拟带NA的数据集
set.seed(123)
N <- 100
K <- 2600
foo <- matrix(runif(N*K), nrow = N, ncol = K)
foo[sample(N*K, 5000)] <- NA  # 加入随机NA
colnames(foo) <- paste0("col", 1:K)

# 预计算每列的非NA行索引（一次性完成）
col_na_free <- lapply(1:K, function(col_idx) which(!is.na(foo[, col_idx])))

# 初始化结果矩阵
result <- matrix(NA, nrow = K, ncol = K)
dimnames(result) <- list(colnames(foo), colnames(foo))
diag(result) <- 1  # 自身的距离相关系数为1，直接赋值

# 只计算上三角，对称填充
for (i in 1:(K-1)) {
  idx_i <- col_na_free[[i]]
  for (j in (i+1):K) {
    idx_j <- col_na_free[[j]]
    # 取两列共同的非NA行
    common_idx <- intersect(idx_i, idx_j)
    if (length(common_idx) < 2) {  # dcor需要至少2个有效样本
      result[i,j] <- NA
      result[j,i] <- NA
      next
    }
    # 直接提取有效行计算，跳过na.omit
    result[i,j] <- energy::dcor(foo[common_idx, i], foo[common_idx, j])
    result[j,i] <- result[i,j]  # 距离相关系数对称
  }
}

方案2：并行计算加速（适合多核CPU） #

列对之间的计算完全独立，可以利用多核CPU并行处理，进一步缩短时间。结合方案1的预存索引，用foreach实现并行：

library(energy)
library(foreach)
library(doParallel)

# 初始化并行集群（留1个核心给系统）
cl <- makeCluster(detectCores() - 1)
registerDoParallel(cl)

set.seed(123)
N <- 100
K <- 2600
foo <- matrix(runif(N*K), nrow = N, ncol = K)
foo[sample(N*K, 5000)] <- NA
colnames(foo) <- paste0("col", 1:K)

col_na_free <- lapply(1:K, function(col_idx) which(!is.na(foo[, col_idx])))

# 并行计算上三角结果
result_list <- foreach(i = 1:(K-1), .combine = rbind) %dopar% {
  idx_i <- col_na_free[[i]]
  row_vals <- rep(NA, K)
  for (j in (i+1):K) {
    idx_j <- col_na_free[[j]]
    common_idx <- intersect(idx_i, idx_j)
    if (length(common_idx) >= 2) {
      row_vals[j] <- energy::dcor(foo[common_idx, i], foo[common_idx, j])
    }
  }
  row_vals
}

# 组装对称矩阵
result <- matrix(NA, nrow = K, ncol = K)
result[1:(K-1), ] <- result_list
result <- t(result) + result  # 填充下三角
diag(result) <- 1
dimnames(result) <- list(colnames(foo), colnames(foo))

# 关闭集群
stopCluster(cl)

方案3：替换为更高效的dcor实现 #

energy::dcor的性能不算最优，可以改用dcor包的同名函数，它内部用了更高效的距离矩阵计算逻辑，速度会快10%-30%左右（需确保结果一致性）：

library(dcor)
# 在方案1/2的代码中，将energy::dcor替换为dcor::dcor即可
result[i,j] <- dcor::dcor(foo[common_idx, i], foo[common_idx, j])

额外优化细节 #

• 优先用matrix而非data.frame存储数据，矩阵的索引操作比数据框快得多。
• 避免在循环内创建临时变量（比如单独赋值x/y），直接将索引后的向量传入dcor函数。
• 提前过滤掉有效样本数不足2的列对，跳过无意义的计算。

内容的提问来源于stack exchange，提问作者umbe1987