如何获取cor()函数计算相关系数时的稀疏矩阵配对样本数量

在R中，当用cor()计算稀疏矩阵的相关系数时，要获取每对变量的有效配对样本数（即两个变量均无缺失值的样本行数），可以通过以下步骤实现，同时保证配对数与cor()使用的样本完全对应：

步骤1：准备稀疏矩阵与依赖包 #

首先确保你加载了Matrix包（处理稀疏矩阵的常用包），并准备好你的稀疏矩阵（比如dgCMatrix类型）。

library(Matrix)

# 模拟一个带缺失值的稀疏矩阵示例
set.seed(123)
# 生成10行10列的稀疏矩阵，70%的概率为1，其余为0
sp_mat <- Matrix(rbinom(100, 1, 0.7), nrow=10, ncol=10, sparse=TRUE)
# 随机设置20个位置为NA（模拟缺失值）
sp_mat[sample(length(sp_mat), 20)] <- NA

步骤2：标记非缺失位置 #

根据你的稀疏矩阵定义，标记每个元素是否为有效样本：

• 如果NA代表缺失值，0是有效数值：

non_missing <- !is.na(sp_mat)

• 如果0代表缺失值，只有非0且非NA的元素是有效样本：

non_missing <- sp_mat != 0 & !is.na(sp_mat)

步骤3：计算有效配对数矩阵 #

利用tcrossprod()函数计算每对变量的共同有效样本数——这个函数会对标记矩阵做转置乘积，每个位置(i,j)的结果就是第i列和第j列同时为有效样本的行数：

# 计算配对数矩阵，结果为稠密矩阵（因为配对数是整数，无稀疏性）
n_pairs <- tcrossprod(non_missing)

# 查看结果
print(n_pairs)

示例输出（简化版）：

[,1] [,2] [,3]
[1,]   8    7    6
[2,]   7    9    8
[3,]   6    8    10

步骤4：计算对应相关系数矩阵 #

使用cor()时，一定要指定use="pairwise.complete.obs"，这样才能保证相关系数的计算与我们的配对数完全对应：

cor_mat <- cor(sp_mat, use="pairwise.complete.obs")

# 查看相关系数矩阵
print(cor_mat)

示例输出（简化版）：

[,1]      [,2]      [,3]
[1,]  1.0000000 0.1212121 -0.340909
[2,]  0.1212121 1.0000000  0.560606
[3,] -0.3409091 0.5606061  1.000000

可选：合并相关系数与配对数 #

如果你想要一个同时显示相关系数和配对数的矩阵，可以用如下代码生成：

# 生成带配对数的结果矩阵，保留两位小数
cor_with_n <- matrix(
  paste0(round(cor_mat, 2), " (", as.matrix(n_pairs), ")"),
  nrow = nrow(cor_mat),
  ncol = ncol(cor_mat),
  dimnames = dimnames(cor_mat)
)

print(cor_with_n)

示例输出（简化版）：

[,1]          [,2]          [,3]         
[1,] "1 (8)"     "0.12 (7)"   "-0.34 (6)"  
[2,] "0.12 (7)"  "1 (9)"      "0.56 (8)"   
[3,] "-0.34 (6)" "0.56 (8)"   "1 (10)"     

注意事项 #

• tcrossprod()对稀疏矩阵的效率很高，它会利用稀疏结构避免不必要的计算，即使矩阵很大也能快速运行。
• 如果你的变量数量极多（比如上万列），配对数矩阵会变成一个巨大的稠密矩阵，可能占用过多内存——这种情况下建议只计算你需要的变量对，而不是全矩阵。
• 一定要保证cor()的use参数与你标记有效样本的逻辑一致，否则配对数和相关系数会不匹配。

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