第 1 章 相关、关联与尤尔–辛普森悖论
若说人类总在追问「为什么」，因果大概就站在那些问题的中心。古希腊人留下过两句，至今读来仍让人安静下来：
「我宁愿发现一条因果律，也不愿做波斯的国王。」
——德谟克利特
「我们尚未把握一物之因，便谈不上认识了它。」
——亚里士多德
可经典统计学教我们最多的，往往是关联，而不是因果。这一章，我们先温柔地认识几种常见的关联度量，再诚实地看看它们够不到的地方——先学会看影子，并不妨碍我们承认：影子还不是那个人。
1.1 统计学的传统看法
有一种很常见的看法：统计学的工作，是推断变量之间的相关或关联。照此说来，因果推断在统计学里几乎没有座位。两句听过许多遍的话，也常常跟着它一起出现：
- 「相关不等于因果。」
- 「你无法用统计学证明因果。」
本书想轻轻换一个立场：统计学对于理解因果，其实很要紧。 接下来我们会一起学习因果推断的形式语言，并在随机化实验与观察性研究里，试着去估计因果效应。
1.2 几种常用的关联度量
1.2.1 相关与回归
两个随机变量 \(Z\) 与 \(Y\) 的 Pearson 相关系数是
\[ \rho_{ZY} = \frac{\operatorname{cov}(Z,Y)}{\sqrt{\operatorname{var}(Z)\operatorname{var}(Y)}}, \]
它告诉我们的，主要是 \(Z\) 与 \(Y\) 之间线性牵绊的强弱。
若把 \(Y\) 对 \(Z\) 做线性回归：
\[ Y = \alpha + \beta Z + \varepsilon, \tag{1.1} \]
其中 \(\mathrm{E}(\varepsilon)=0\)，\(\mathrm{E}(\varepsilon Z)=0\)，则
\[ \beta = \frac{\operatorname{cov}(Z,Y)}{\operatorname{var}(Z)} = \rho_{ZY}\sqrt{\frac{\operatorname{var}(Y)}{\operatorname{var}(Z)}}. \]
所以 \(\beta\) 与 \(\rho_{ZY}\) 总是同号——它们朝同一个方向点头。
也可以把 \(Y\) 同时对 \(Z\) 与 \(X\) 回归：
\[ Y = \alpha + \beta Z + \gamma X + \varepsilon, \tag{1.2} \]
其中 \(\mathrm{E}(\varepsilon)=0\)，\(\mathrm{E}(\varepsilon Z)=0\)，\(\mathrm{E}(\varepsilon X)=0\)。这时人们常说：在「保持 \(X\) 不变」「条件于 \(X\)」或「控制了 \(X\)」之后，\(\beta\) 是 \(Z\) 对 \(Y\) 的「效应」。附录 B 会帮你把线性回归再温一遍。
有趣的是：上面两个回归里的 \(\beta\)，可以不一样，甚至可以一正一负。下面这段 R 代码，重新看了 Hainmueller (2012) 用过的 LaLonde 观察性数据。问题很朴素：就业培训项目，对收入究竟有没有帮助？把所有协变量都放进模型时，treat 的系数是 \(1067.5461\)；什么都不控制时，却变成了 \(-8506.4954\)。同一扇窗，擦亮与蒙尘，窗外的景色可以完全相反。
> dat <- read.table("cps1re74.csv", header = TRUE)
> dat$u74 <- as.numeric(dat$re74 == 0)
> dat$u75 <- as.numeric(dat$re75 == 0)
>
> ## 对结果做线性回归；. 表示用 dat 中其余全部变量
> lmoutcome = lm(re78 ~ ., data = dat)
> round(summary(lmoutcome)$coef[2, ], 3)
Estimate Std. Error t value Pr(>|t|)
1067.546 554.060 1.927 0.054
>
> lmoutcome = lm(re78 ~ treat, data = dat)
> round(summary(lmoutcome)$coef[2, ], 3)
Estimate Std. Error t value Pr(>|t|)
-8506.495 712.766 -11.934 0.0001.2.2 列联表
两个二值变量 \(Z\) 与 \(Y\) 的联合分布，可以放进一张小小的 \(2\times 2\) 表。记 \(p_{zy}=\operatorname{pr}(Z=z,Y=y)\)：
| \(Y=1\) | \(Y=0\) | |
|---|---|---|
| \(Z=1\) | \(p_{11}\) | \(p_{10}\) |
| \(Z=0\) | \(p_{01}\) | \(p_{00}\) |
若把 \(Z\) 看作处理（或暴露），\(Y\) 看作结果，我们可以把三种常见度量写清楚：
风险差（risk difference）
\[ \mathrm{rd} = \operatorname{pr}(Y=1\mid Z=1)-\operatorname{pr}(Y=1\mid Z=0) = \frac{p_{11}}{p_{11}+p_{10}}-\frac{p_{01}}{p_{01}+p_{00}}, \]
风险比（risk ratio）
\[ \mathrm{rr} = \frac{\operatorname{pr}(Y=1\mid Z=1)}{\operatorname{pr}(Y=1\mid Z=0)} = \frac{p_{11}/(p_{11}+p_{10})}{p_{01}/(p_{01}+p_{00})}, \]
优势比（odds ratio）1
\[ \mathrm{or} = \frac{\operatorname{pr}(Y=1\mid Z=1)/\operatorname{pr}(Y=0\mid Z=1)} {\operatorname{pr}(Y=1\mid Z=0)/\operatorname{pr}(Y=0\mid Z=0)} = \frac{p_{11}p_{00}}{p_{10}p_{01}}. \]
这些名字多半来自流行病学。那里的「结果」常常是疾病，用「风险」称呼发病的概率，听起来很自然，也带着一点对真实痛苦的体贴。
关于它们，有几件小事值得记住。
命题 1.1 (1) 下面几件事其实是一回事2：\(Z \perp\!\!\!\perp Y\)，\(\mathrm{rd}=0\)，\(\mathrm{rr}=1\)，以及 \(\mathrm{or}=1\)。(2) 若所有 \(p_{zy}\) 都为正，则 \(\mathrm{rd}>0\) 等价于 \(\mathrm{rr}>1\)，也等价于 \(\mathrm{or}>1\)。(3) 当 \(\operatorname{pr}(Y=1\mid Z=1)\) 与 \(\operatorname{pr}(Y=1\mid Z=0)\) 都很小时，\(\mathrm{or}\approx\mathrm{rr}\)。
(1)(2) 的证明留给习题 1.1。第 (3) 点说得不那么正式：对罕见病，\(p\approx 0\) 时，优势 \(p/(1-p)\) 几乎就是概率本身——泰勒展开告诉我们 \(p/(1-p)=p+p^2+\cdots\approx p\)。所以在流行病学里，若结果是罕见病，常常可以放心地假定这两个条件概率都很小。
若把概率换成「给定另一个变量 \(X\)」之后的条件概率，也能写出 \(\mathrm{rd}\)、\(\mathrm{rr}\)、\(\mathrm{or}\) 的条件版本。
用计数 \(n_{zy}=\#\{i:Z_i=z,Y_i=y\}\)，观测数据也可以排成同样的小表。再用样本比例 \(\hat p_{zy}=n_{zy}/n\) 去替换真实概率，就能估计这些量（见附录 A.3.2）。在 R 里，fisher.test 与 chisq.test 可以帮你检验 \(Z\) 与 \(Y\) 是否独立。
例 1.1 Bertrand and Mullainathan (2004) 做过一项关于简历的随机化实验：在虚构简历上随机写上「听起来像黑人」或「听起来像白人」的姓名，投向波士顿与芝加哥的招聘广告，看看面试回电有没有差别。数据长这样：
> resume = read.csv("resume.csv")
> Alltable = table(resume$race, resume$call)
> Alltable
0 1
black 2278 157
white 2200 235两行人数一样多，却能看见：听起来像白人的姓名，收到了更多回电。Fisher 精确检验告诉我们，这个差别很难用偶然解释：
> fisher.test(Alltable)
Fisher's Exact Test for Count Data
data: Alltable
p-value = 4.759e-05
alternative hypothesis: true odds ratio is not equal to 1
95 percent confidence interval:
1.249828 1.925573
sample estimates:
odds ratio
1.5497321.3 尤尔–辛普森悖论的一个例子
1.3.1 数据
经典的肾结石例子来自 Charig et al. (1986)。\(Z=1\) 是开放手术，\(Z=0\) 是小穿刺；\(Y=1\) 表示成功，\(Y=0\) 表示失败：
| \(Y=1\) | \(Y=0\) | |
|---|---|---|
| \(Z=1\) | 273 | 77 |
| \(Z=0\) | 289 | 61 |
估出来的风险差是
\[ \widehat{\mathrm{rd}} = \frac{273}{273+77}-\frac{289}{289+61} = 78\%-83\% = -5\% < 0. \]
看起来，小穿刺更好一点。
可这些数据并不是随机对照试验（RCT）3 得来的。接受两种手术的人，本就可能很不一样。研究里还藏着一个「躲在幕后的变量」：结石大小——有人小，有人大。我们按大小把数据轻轻拆开。
小结石：
| \(Y=1\) | \(Y=0\) | |
|---|---|---|
| \(Z=1\) | 81 | 6 |
| \(Z=0\) | 234 | 36 |
大结石：
| \(Y=1\) | \(Y=0\) | |
|---|---|---|
| \(Z=1\) | 192 | 71 |
| \(Z=0\) | 55 | 25 |
两张表加回去，应当正好还原第一张：
\[ 81+192=273,\quad 6+71=77,\quad 234+55=289,\quad 36+25=61. \]
小结石组：
\[ \widehat{\mathrm{rd}}_{\text{smaller}} = \frac{81}{81+6}-\frac{234}{234+36} = 93\%-87\% = 6\% > 0. \]
大结石组：
\[ \widehat{\mathrm{rd}}_{\text{larger}} = \frac{192}{192+71}-\frac{55}{55+25} = 73\%-69\% = 4\% > 0. \]
两边都在说：开放手术更好。可合在一起，却变成：开放手术更差。于是
\[ \widehat{\mathrm{rd}}<0,\qquad \widehat{\mathrm{rd}}_{\text{smaller}}>0,\qquad \widehat{\mathrm{rd}}_{\text{larger}}>0. \]
用人话说：无论小结石还是大结石，处理 \(1\) 都更好；可对「所有人」一汇总，处理 \(1\) 反而更差。若你真心想知道哪种手术更有帮助，这种矛盾会让人困惑，甚至有点委屈。统计学里，这叫尤尔–辛普森悖论：边际关联的方向，与每一层里的条件关联，正好相反。
1.3.2 解释
令 \(X=1\) 表示小结石，\(X=0\) 表示大结石。先看谁更常接受开放手术：
\[ \begin{aligned} &\widehat{\operatorname{pr}}(Z=1\mid X=1)-\widehat{\operatorname{pr}}(Z=1\mid X=0)\\ &= \frac{81+6}{81+6+234+36}-\frac{192+71}{192+71+55+25}\\ &= 24\%-77\% = -53\% < 0. \end{aligned} \]
大结石患者更常去做开放手术。也就是说，\(X\) 与 \(Z\) 是负向关联的。
再看成功率。在开放手术下：
\[ \begin{aligned} &\widehat{\operatorname{pr}}(Y=1\mid Z=1,X=1)-\widehat{\operatorname{pr}}(Y=1\mid Z=1,X=0)\\ &= \frac{81}{81+6}-\frac{192}{192+71} = 93\%-73\% = 20\% > 0; \end{aligned} \]
在小穿刺下：
\[ \begin{aligned} &\widehat{\operatorname{pr}}(Y=1\mid Z=0,X=1)-\widehat{\operatorname{pr}}(Y=1\mid Z=0,X=0)\\ &= \frac{234}{234+36}-\frac{55}{55+25} = 87\%-69\% = 18\% > 0. \end{aligned} \]
无论哪种手术，小结石患者都更容易成功。于是，在两个处理水平上，\(X\) 与 \(Y\) 都是正向关联的。
原书图 1.1 把这些关系画成一张安静的图：处理到结果有一条正向的小路 \(Z\to Y\)，同时还有一条更强的、绕经病情的负向小路 \(Z\leftarrow X\to Y\)。两相抵消之后，整体看起来，处理与结果竟是负相关。用人话说得更软一点：当效果较差的处理，更常被用在病情较轻的人身上时，它就会显得像更好的选择——不是世界故意捉弄我们，只是汇总的方式，把故事讲偏了。
一般地，只要 \(X\) 同时与 \(Z\)、\(Y\) 有关，\(Z\) 与 \(Y\) 的表面关联，就可能和「给定 \(X\) 之后」的关联，在方向上分道扬镳。这时我们说：\(X\) 是一个混杂变量（confounder），或者说 \(Z\) 与 \(Y\) 的关系被 \(X\) 混杂了。第 10 章和第 17 章会更仔细地陪你看这件事。
1.3.3 尤尔–辛普森悖论的几何图景
把全体数据、以及按 \(X=1\) 与 \(X=0\) 拆开的两张表放在一起，原书图 1.2 给出了一种几何直觉：纵轴是成功数，横轴是失败数。在每一层里，代表处理的那条线段斜率更大，于是处理显得有益；可把两层拼回去之后，斜率的大小关系却翻了过来——处理又显得有害。悖论，就这样从几何里慢慢长出来。
1.4 伯克利研究生院录取数据
Bickel et al. (1975) 看过伯克利研究生院男女学生的录取率。R 包 datasets 里有 UCBAdmissions。按六个最大院系整理后，跨院系汇总会得到：
> UCBAdmissions.sum = apply(UCBAdmissions, c(1, 2), sum)
> UCBAdmissions.sum
Admit
Gender Admitted Rejected
Male 1198 1493
Female 557 1278用下面这个小函数，可以得到风险差与 \(p\) 值：
> risk.difference = function(tb2)
+ {
+ p1 = tb2[1, 1]/(tb2[1, 1] + tb2[1, 2])
+ p2 = tb2[2, 1]/(tb2[2, 1] + tb2[2, 2])
+ testp = chisq.test(tb2)
+ return(list(p.diff = p1 - p2, pv = testp$p.value))
+ }合在一起看，男女录取率相差不小，而且很显著：
> risk.difference(UCBAdmissions.sum)
$p.diff
[1] 0.1416454
$pv
[1] 1.055797e-21可一旦按院系分开，差异往往变小，也多半不再显著；在院系 A，差异甚至显著为负：
> round(P.diff, 2)
[1] -0.20 -0.05 0.03 -0.02 0.04 -0.01
> round(PV, 2)
[1] 0.00 0.77 0.43 0.64 0.37 0.64远远望去，像一桩关于性别的沉重故事；走近每一间院系，故事的轮廓却柔和下来，甚至悄悄换了方向。悖论擅长的，正是这种事：它不总是恶意，有时只是提醒我们——整体与局部，可以同时诚实，又彼此不同。
1.5 习题
1.1 \(2\times 2\) 表中的独立性
证明命题 1.1 的 (1) 与 (2)。
1.2 更多尤尔–辛普森悖论的例子
请给出一个会出现该悖论的 \(2\times 2\times 2\) 数值例子；再找一个真实世界里的例子。
1.3 相关与偏相关
考虑三维正态随机向量：
\[ \begin{pmatrix} X \\ Y \\ Z \end{pmatrix} \sim \mathcal{N}\!\left( \begin{pmatrix} 0\\0\\0 \end{pmatrix}, \begin{pmatrix} 1 & \rho_{XY} & \rho_{XZ}\\ \rho_{XY} & 1 & \rho_{YZ}\\ \rho_{XZ} & \rho_{YZ} & 1 \end{pmatrix} \right). \]
证明偏相关系数
\[ \rho_{YZ\mid X} = \frac{\rho_{YZ}-\rho_{YX}\rho_{ZX}} {\sqrt{1-\rho_{YX}^2}\,\sqrt{1-\rho_{ZX}^2}}, \]
并给出一个数值例子，使 \(\rho_{YZ}>0\) 而 \(\rho_{YZ\mid X}<0\)。
注： 这是正态情形下的尤尔–辛普森悖论。证明时可借助附录 A.1.2。
1.4 设定搜索
回到 1.2.1 的 LaLonde 数据：结果是 re78，处理是 treat，另有 10 个协变量，于是有 \(2^{10}=1024\) 种回归设定。请把它们都跑一遍，看看处理系数有多少显著为正、显著为负、不显著；也欢迎写下你觉得有意思的发现。
1.5 再谈种族与回电
请把 Bertrand and Mullainathan (2004) 的数据按男性和女性分开再分析。你看见了什么？
1.6 推荐阅读
Bickel et al. (1975) 是 1.4 节的原始论文。Pearl and Mackenzie (2018) 后来又从因果推断的角度，重新探望了这项研究。