R 语言相关分析：cor()、cor.test() 与相关矩阵

两个变量之间的相关性：cor()

相关分析要回答的问题是：当一个变量变大时，另一个变量是倾向于也变大（正相关）、变小（负相关），还是各走各的（接近零）？在 R 中只需一次调用：

print(cor(mtcars$wt, mtcars$mpg))

结果约为 −0.87：车越重，油耗表现越差，而且这种关系很强。这个数就是 Pearson 相关系数，通常记作 r。

读懂 r

相关系数总是落在 −1 到 +1 之间。符号表示方向，绝对值大小表示强弱：

| |r| | 常见解读 | | --- | --- | | 0.0 至 0.2 | 可以忽略 | | 0.2 至 0.4 | 弱 | | 0.4 至 0.6 | 中等 | | 0.6 至 0.8 | 强 | | 0.8 至 1.0 | 非常强 |

把这些区间当作讨论的起点，而不是定律：在物理学中，0.6 的 r 令人失望；在心理学中，它却是职业生涯的高光时刻。有两个性质值得牢记：r 没有单位（用以吨为单位的重量和 mpg 计算相关系数，与用以千克为单位的重量计算，得到的 r 完全一样，因为 r 是在标准化后的值上计算的）；而且 r 只衡量线性关联：一个完美的 U 形关系，r 可能约等于 0。

还有一句必须说的话：相关不等于因果。在一年中的各个月份里，冰淇淋销量和溺水死亡人数高度相关，但这不是因为冰淇淋让人溺水，而是夏天同时推高了两者。相关性告诉你两个变量一起变化，却对原因只字不提。可能是 x 影响 y，可能是 y 影响 x，也可能是第三个因素（在这个例子中是季节）同时影响两者。要在这些可能性之间做出判断，需要实验或严谨的因果推理，而不是更大的 r。

Spearman 和 Kendall：当 Pearson 不适用时

Pearson 基于原始值计算，这让它对离群值很敏感，也看不出弯曲的关系。method 参数可以切换到基于秩的替代方法：

x <- c(1, 2, 3, 4, 5, 6, 7, 8)
y <- c(1, 4, 9, 16, 25, 36, 49, 300)  # increasing, with a wild last value

print(cor(x, y))                        # Pearson - dragged by the outlier
print(cor(x, y, method = "spearman"))   # ranks: perfectly monotonic
print(cor(x, y, method = "kendall"))

Spearman 把每个值替换成它的秩，然后对秩计算 Pearson 相关系数。由于 x 增大时 y 总是增大，所有的秩都对得上，Spearman 给出的结果正好是 1：离群值有多大不再重要，重要的只是它的位置。当数据是有序的（问卷量表）、严重偏斜，或者关系单调但不是直线时，就用 Spearman。Kendall 通过一致对和不一致对来回答类似的问题；它在小样本中更稳健，但速度更慢，而 Spearman 是更常用的默认选择。

相关矩阵

要一次性查看多个变量之间的关系，就把几个数值列交给 cor()：

m <- cor(mtcars[, c("mpg", "wt", "hp", "disp")])
print(round(m, 2))

每个变量都与其他每个变量配对，对角线上全是 1（任何变量与自身都完全相关），两侧互为镜像。保留两位小数比听上去重要得多：不取整的矩阵就是一堵数字墙，而矩阵的意义就在于能一眼扫过去。这里扫一眼就能看出，mpg 与另外三个变量都呈负相关（更重、马力更大、排量更大的车更费油），而 wt、hp 和 disp 两两之间都是很强的正相关。这组“大车”变量会在你学到线性回归及其令人头疼的多重共线性问题时派上用场。

缺失值：use 参数

有缺失数据时，cor() 默认返回 NA，而不是去猜：

d <- mtcars[, c("mpg", "wt", "hp")]
d$wt[c(3, 7)] <- NA  # knock two values out

print(cor(d)[1, 2])
print(cor(d, use = "complete.obs")[1, 2])
print(cor(d, use = "pairwise.complete.obs")[1, 2])

• use = "complete.obs" 删除所有含有任意 NA 的行，再用剩下的行计算整个矩阵。结果一致，但浪费数据：缺失一个 wt，这一行也会从 mpg 与 hp 这一对的计算中被删掉。
• use = "pairwise.complete.obs" 对每个单元格，使用这一对变量都有值的所有行来计算。它保留了更多数据，但不同单元格基于不同的数据子集，极少数情况下可能得到一个内部不一致的矩阵。

如果只有零星几个 NA，两种都可以；只要说明你用的是哪一种就行。

显著吗？cor.test()

cor() 只给出一个数，却无法告诉你它是不是噪声。cor.test() 会进行假设检验：

print(cor.test(mtcars$wt, mtcars$mpg))

逐块看一下输出：

• t = −9.56, df = 30：检验统计量。原假设是真实相关系数为零；观测到的 r 被转换成自由度为 n − 2（32 辆车 − 2）的 t 统计量。
• p-value = 1.29e-10：如果真实相关系数为零，在 32 个样本中观测到离零这么远的 r，概率约为 0.0000000001。这是关联真实存在的压倒性证据。不过要记住，p 值说明的是 r 是否不为零，而不是这种关系是否很大，或者是否存在因果。
• 95 percent confidence interval: −0.93 to −0.74：真实相关系数的合理范围。它往往比 p 值更有用：即使是这个区间中最乐观的一端，也是很强的负相关。
• sample estimates: cor = −0.87：和 cor() 给出的是同一个数。

对小样本要格外小心：当 n = 10 时，±0.5 的相关系数纯靠运气出现的频率高得惊人，而很宽的置信区间会如实告诉你这一点。报告结果时要给出区间，而不只是 p 值。

亲眼看看：一定要画图

相关系数把整个关系压缩成一个数，而这种压缩可能掩盖弯曲、聚类，或者某一个点包揽了全部作用的情况。在相信任何 r 之前，先看看散点图：

plot(mtcars$wt, mtcars$mpg)              # one pair
pairs(mtcars[, c("mpg", "wt", "hp", "disp")])  # every pair in the matrix

pairs() 会画出一个与相关矩阵对应的散点图网格，这是确认这些数字是否符合你的理解的最快方法。要绘制精美的热力图式相关矩阵，corrplot 包（先 install.packages("corrplot")，再 corrplot(cor(m))）是标准工具。

要点回顾

• cor(x, y) 给出 Pearson 的 r：符号表示方向，绝对值表示强弱，总是在 [−1, 1] 之间，没有单位。
• 相关性衡量的是线性的共同变化，对因果关系什么也没说：潜在的第三个变量永远是一个可能的解释。
• 用 method = "spearman" 基于秩计算：适用于有序数据、离群值、单调但弯曲的关系。
• 对数值列使用 cor(df) 得到相关矩阵；用 round(, 2) 取整，数据有缺失时注意 use = 参数。
• cor.test(x, y) 额外给出 p 值和置信区间：报告时要给出区间。
• 相信这个数之前，一定要先看散点图。

下一篇：当问题从“它们是否一起变化？”变成“这一组的均值和那一组的不同吗？”时，就要用到 t 检验。