第3章 图形化因果模型
第 3 章 图形因果模型
在第 \(1\) 章中，你了解到因果推断可拆解为两个问题：识别和估计。本章中，你将更深入探究识别部分，而识别可以说是最具挑战性的环节。本章理论性较强，因为你会接触图形模型，且不一定会用数据估计其参数。但别被这唬住，识别是因果推断的核心，学习其理论对于解决现实生活中的因果问题至关重要。在本章中，你将：
初步认识图形模型，了解什么是因果关系的图形模型、关联在图中如何传递，以及如何使用现成软件查询图形。
透过图形模型的视角，重新审视识别这一概念。
了解阻碍识别的两种极为常见的偏差来源、它们的因果图结构，以及应对这些偏差的方法 。
思考因果关系
你有没有注意到，\(YouTube\) 视频里的那些厨师对食物的描述多么精妙？“把酱汁收浓，直到呈现出如天鹅绒般的质地。” 要是你刚开始学做饭，根本不明白这是什么意思。直接告诉我这东西该在炉子上煮多久不就行了嘛！在因果关系的学习上，情况也类似。假设你走进一家酒吧，听到人们讨论因果关系（这酒吧很可能挨着经济系 ）。这种情况下，你会听到他们说收入方面的混杂因素如何让识别移民对社区失业率的影响变得困难，所以他们得用工具变量。而现在，你可能压根不懂他们在说啥。没关系。对于因果推断这门语言，你才刚触及皮毛。你已经学了一点关于反事实结果和偏差的知识，足以让你理解因果推断试图解决的关键问题，也足以让你明白因果推断最有力工具（随机对照试验 ）背后的原理。但这种工具不会总有可用之机，或者单纯就不管用（你很快会在 “选择偏差” 一节看到 ）。当你遇到更具挑战性的因果推断问题时，你也需要对因果推断的语言有更广泛的理解，这样才能正确理解自己面临的情况以及如何应对。
一套表述清晰的语言能让你思路清晰。本章旨在拓展你的因果推断词汇量 。你可以把图形模型视为因果关系的基础语言之一。它们是构建因果推断问题、明确（甚至可视化 ）识别假设的有力方式。图形模型能让你的想法对他人和你自己都清晰透明 。
结构因果模型
一些科学家用 “结构因果模型\(（SCM, Structural \space Causal \space Model）\)” 这一术语来指代因果推断的统一语言。这些模型由图形和因果方程构成。在这里，我会主要聚焦于结构因果模型的图形方面 。
作为进入奇妙的图形世界的起点，让我们以之前估算电子邮件对转化率影响的例子为例。在那个例子中，处理 \(T\) 是交叉销售邮件，结果 \(Y\) 是客户是否转化为购买新产品 。
import warnings
warnings.filterwarnings('ignore')
import pandas as pd
import numpy as np
import graphviz as gr
import seaborn as sns
from matplotlib import pyplot as plt
color=['0.3', '0.5', '0.7', '0.9']
linestyle=['-', '--', ':', '-.']
marker=['o', 'v', 'd', 'p']
pd.set_option('display.max_rows', 6)
gr.set_default_format("png")## 'pdf'import pandas as pd # for data manipulation
import numpy as np # for numerical computation
data = pd.read_csv("./data/cross_sell_email.csv")
data## gender cross_sell_email age conversion
## 0 0 short 15 0
## 1 1 short 27 0
## 2 1 long 17 0
## .. ... ... ... ...
## 320 0 no_email 15 0
## 321 1 no_email 16 0
## 322 1 long 24 1
##
## [323 rows x 4 columns]
让我们回顾上一章的内容，在这个问题中，\(T\) 是随机分配的。因此，我们可以认为处理与潜在结果相互独立，即 \((Y_0, Y_1) \perp T\)，这使得相关性等同于因果关系：
\(E[Y_1 - Y_0] = E[Y|T = 1] - E[Y|T = 0]\)
重要的是，仅通过查看数据，绝对无法确定独立性假设是否成立。你只能依据关于处理分配机制的信息来判定该假设成立。也就是说，你知道邮件是随机发送的 。
可视化因果关系
你可以把这些知识编码到一个图形中，该图形能体现你对 “什么导致什么” 的看法。在这个简单例子里，假设你认为交叉销售邮件会导致转化。你还认为自己测量过的其他变量，比如年龄\(（age）\)和性别\(（gender）\)，也会导致转化。此外，你也可以把未测量的变量添加到图形里。我们通常用字母 \(U\) 来表示它们，因为这些变量是未观测到的\(（unobserved）\)。可能有很多未观测到的变量会导致转化（比如客户收入、社会背景），也会影响年龄（比如你的产品对不同人群的吸引力、公司运营所在的城市）。但由于你没测量这些变量，你可以把所有这些因素归到一个 \(U\) 节点里，代表所有这些未测量的变量。最后，你可以添加一个指向 \(T\) 的随机化节点，以此体现你所知道的 “交叉销售邮件是随机发送的” 这一事实。
有向无环图\(（DAG）\)
你可能会发现人们把因果图称为 \(DAGs\) 。这个缩写代表 “有向无环图\(（directed \space acyclic \space graph）\)” 。“有向\(（directed）\)” 意味着图中的边有方向，这和无向图（比如社交网络）不同。“无环\(（acyclic）\)” 意味着图中没有环路或循环。因果图通常是有向且无环的，因为因果关系是不可逆的 。
要把你的那些看法添加到图形中并直观地呈现出来，你可以使用 \(graphviz\)（图形可视化工具） ：
import graphviz as gr
g_cross_sell = gr.Digraph()
g_cross_sell.edge("U", "conversion")
g_cross_sell.edge("U", "age")
g_cross_sell.edge("U", "gender")
g_cross_sell.edge("rnd", "cross_sell_email")
g_cross_sell.edge("cross_sell_email", "conversion")
g_cross_sell.edge("age", "conversion")
g_cross_sell.edge("gender", "conversion")
# 生成图像文件（默认格式为 PNG，也可指定 svg/pdf 等） filename 为输出文件名（不含后缀），format 为图像格式,view=True 会自动打开系统默认图像查看器
g_cross_sell.render(filename="cross_sell_graph", format="png", view=True)## 'cross_sell_graph.png'
图中的每个节点都是一个随机变量。你可以用箭头（或边）来表示一个变量是否导致另一个变量。在这个图形模型中，你是在表明 \(cross\_sell\_email\)会导致\(conversion\)，\(U\)会导致\(age\)、\(conversion\)和\(gender\)等等。图形模型这种表达方式会帮助你理清对因果关系的思考，因为它把你对世界运行方式的看法明确呈现出来。要是你觉得这不切实际 —— 毕竟，在如今的数据应用中，通常存在成百上千个变量，你不可能把所有变量都编码进去 —— 别担心，你不需要这么做。实际上，你可以通过合并节点来大幅简化问题，同时保留你想要传达的大致因果逻辑。比如，你可以把前面图中的可观测变量合并成一个\(X\)节点。由于这些变量都由 \(U\) 导致，且都会导致\(conversion\)，把它们合并后，你的因果逻辑依然完整。
另外，当你要表示已被随机化或干预过的变量时，你可以直接移除指向该变量的所有入向箭头 ：
# rankdir:LR layers the graph from left to right
g_cross_sell = gr.Digraph(graph_attr={"rankdir": "LR"})
g_cross_sell.edge("U", "conversion")
g_cross_sell.edge("U", "X")
g_cross_sell.edge("cross_sell_email", "conversion")
g_cross_sell.edge("X", "conversion")
g_cross_sell.render(filename="cross_sell_graph", format="png", view=True)## 'cross_sell_graph.png'
有意思的是，要意识到在有向无环图中，或许最重要的信息恰恰是图中没有的部分：从一个变量指向另一个变量的边缺失，意味着假设这两个变量之间不存在直接因果联系。比如，在前面的图中，你是在假设不存在同时导致处理\(（treatment）\)和结果\(（outcome）\)的因素 。
就像学习任何一种语言时一样，你看着这些内容，可能会觉得不是所有内容都能完全理解。这很正常。我大可以给你抛出一大堆规则和最佳实践，教你如何在图中表示变量之间的因果关系。但这很可能是效率最低的学习方式。相反，我的计划是让你接触大量例子。随着时间推移，你会掌握窍门的。眼下，我只是希望你记住，图形是理解 “相关性并非因果性” 的极为有力的工具 。
咨询顾问值得聘请吗？
为了展现有向无环图的作用，咱们来看一个更有意思的例子，其中处理并非随机分配的。假设你是一家公司的经理，正在考虑是否聘请一些顶尖咨询顾问。你知道他们收费高昂，但也清楚，他们因与商界最优秀的公司合作过，掌握着专业知识。
而情况更复杂的是，你不确定这些顶尖咨询顾问会助力你的公司业务改善，还是仅仅因为只有盈利极为可观的公司才雇得起他们，所以他们的存在才与强劲的业务表现相关。要是有人能随机安排咨询顾问的介入情况，那就太棒了，因为这样一来，回答这个问题会很轻松。但显然，你没这样的便利条件，所以得另想办法。想必你现在已经看出来了，这是一个要从关联中梳理出因果关系的问题。为理解它，你可以把自己对其因果机制的看法编码到一个图形中：
g_consultancy = gr.Digraph(graph_attr={"rankdir": "LR"})
g_consultancy.edge("U1", "profits_next_6m")
g_consultancy.edge("U2", "consultancy")
g_consultancy.edge("U3", "profits_prev_6m")
g_consultancy.edge("consultancy", "profits_next_6m")
g_consultancy.edge("profits_prev_6m", "consultancy")
g_consultancy.edge("profits_prev_6m", "profits_next_6m")
g_consultancy.render(filename="cross_sell_graph", format="png", view=True)## 'cross_sell_graph.png'
注意我是如何给每个这些变量添加 \(U\) 节点的，以此表示存在其他我们无法测量、却会对它们产生影响的因素这一事实。由于图形通常表示随机变量，预计会有一个随机成分导致所有变量，而那些 \(U\) 就是用来表示这个随机成分的。不过，它们不会对我接下来要讲的因果逻辑有任何补充，所以我也完全可以省略它们：
g_consultancy = gr.Digraph(graph_attr={"rankdir":"LR"})
g_consultancy.edge("consultancy", "profits_next_6m")
g_consultancy.edge("profits_prev_6m", "consultancy")
g_consultancy.edge("profits_prev_6m", "profits_next_6m")
g_consultancy.render(filename="cross_sell_graph", format="png", view=True)## 'cross_sell_graph.png'
在这里，我想说的是，一家公司的过往业绩会促使该公司聘请顶尖咨询顾问。要是公司经营出色，就有能力支付这项昂贵的服务；要是经营得不太好，就无力承担。因此，过往业绩（此处以过往利润衡量 ）决定了一家公司聘请咨询顾问的可能性。要记住，这种关系不一定是确定性的。我只是说，经营状况良好的公司更有可能聘请顶尖咨询顾问。
不仅如此，过去 \(6\) 个月经营良好的公司，在接下来 \(6\) 个月也很可能继续表现出色。当然，这并非绝对会发生，但平均而言是这样，这就是为什么还有一条从过往业绩指向未来业绩的边。最后，我添加了一条从咨询服务指向公司未来业绩的边。你的目标是了解这种关联的强度，这是你关心的因果关系：咨询服务真的会促使公司业绩提升吗？
回答这个问题并不简单，因为咨询服务和未来业绩之间的关联存在两个来源，一个是因果性的，另一个不是。要理解并梳理清楚它们，你首先需要快速了解一下关联在因果图中是如何传递的 。
图形模型速成课
学校会开设一整个学期的图形模型课程。当然了，要是你想深入研究图形模型，这对你理解因果推断会大有益处。但就本书而言，至关重要的是你要明白图形模型隐含着何种独立性以及条件独立性假设。你会看到，关联在图形模型中的传递，就如同水流经溪流一般。你可以阻止或让这种传递发生，这取决于你如何处理图中的变量。
为理解这一点，咱们来研究一些常见的图形结构和例子。它们会相当直观，但却是理解图形模型中关联、独立性和条件独立性传递的充分基石 。
链式结构
首先，看这个非常简单的图。它被称为链式结构。在这里，\(T\) 导致 \(M\)，\(M\) 导致 \(Y\)。你有时可以把中间节点称为中介变量，因为它在 \(T\) 和 \(Y\) 之间的关系中起中介作用 ：
g = gr.Digraph(graph_attr={"rankdir": "LR"})
g.edge("T", "M")
g.edge("M", "Y")
g.node("M", "M")
g.edge("causal knowledge", "solve problems")
g.edge("solve problems", "job promotion")
g.render(filename="chains", format="png", view=True)## 'chains.png'
在第一个图中，尽管因果关系仅沿箭头方向传递，但相关性是双向流动的。举个更具体的例子，假设了解因果推断知识能提升你的问题解决能力，而问题解决能力会增加你获得晋升的几率。所以，因果知识会使你的问题解决能力提升，进而让你获得工作晋升。在这里，你可以说工作晋升依赖于因果知识。因果专业知识越丰富，你获得晋升的几率就越大。同样，获得晋升的几率越大，你具备因果知识的可能性也就越大。否则，就很难获得晋升。换句话说，工作晋升与因果推断专业知识相关联，就如同因果推断专业知识与工作晋升相关联一样，即便只有一个方向存在因果关系。当两个变量相互关联时，你可以说它们是相关联的\(（dependent）\)或不独立的\(（not \space independent ）\)：
\(T \not\perp Y\)
现在，咱们把中间变量固定住。你可以这么做，比如在例子中只观察具备相同 \(M\)（即相同问题解决能力）的人。正式来说，你可以称这是在对 \(M\) 进行条件限定。这种情况下，相关性就被阻断了。所以，在给定 \(M\) 的条件下，\(T\) 和 \(Y\) 是相互独立的。你可以用数学方式表示为：
\(T \perp Y \mid M\)
为了表明我们在对一个节点进行条件限定，我会把它着色：
g = gr.Digraph(graph_attr={"rankdir": "LR"})
g.edge("T", "M")
g.edge("M", "Y")
g.node("M", "M")
g.node("M", color="lightgrey", style="filled")
g.edge("causal knowledge", "solve problems")
g.edge("solve problems", "job promotion")
g.node("solve problems", color="lightgrey", style="filled")
g.render(filename="chains_m", format="png", view=True)## 'chains_m.png'
要理解这在我们的例子中意味着什么，想想以人们的问题解决能力为条件。如果你观察一群问题解决能力相同的人，知道其中哪些人擅长因果推断，并不会为你提供关于他们获得工作晋升几率的更多信息。用数学术语来说：
\(E[\text{Promotion} \mid \text{Solve problems}, \text{Causal knowledge}] = E[\text{Promotion} \mid \text{Solve problems}]\)
反过来也成立；一旦我知道你解决问题的能力有多强，知道你的工作晋升状态，并不会为我提供关于你了解因果推断的可能性的更多信息。
一般规则是，如果你有一个像前面图中那样的链式结构，当你以中间变量 \(M\) 为条件时，从 \(T\) 到 \(Y\) 路径上的关联就会被阻断。或者说：
\(T \not\perp Y\)
但
\(T \perp Y \mid M\)
叉形结构
接着，咱们来看叉形结构。在这种结构里，存在一个共同原因：同一个变量会导致图中另外两个变量。在叉形结构中，相关性会沿着箭头反向传递：
g = gr.Digraph()
g.edge("X", "Y")
g.edge("X", "T")
g.node("X", "X")
g.edge("statistics", "causal inference")
g.edge("statistics", "machine learning")
g.render(filename="forks", format="png", view=True)## 'forks.png'
比如，假设你对统计学的知识会让你对因果推断和机器学习都了解更多。然而，懂因果推断对你学习机器学习没帮助，反之亦然，所以这两个变量之间没有边。
这个图告诉你，要是你不知道某人的统计知识水平，那么了解到他们擅长因果推断，会让你更倾向于认为他们也擅长机器学习，即便因果推断对机器学习没帮助。这是因为，就算你不知道某人的统计知识水平，你也能从他们的因果推断知识推断出来。要是他们擅长因果推断，他们很可能统计学得好，也就更有可能擅长机器学习。叉形结构顶端的变量，即便彼此之间没有因果关系，也会一同变化，原因很简单，它们都是由同一个因素导致的。在因果推断文献中，当处理和结果之间存在共同原因时，我们把这个共同原因称为混杂因素\(（confounder）\) 。
叉形结构在因果推断中极为重要，值得再举一个例子。你知道科技行业的招聘人员有时会让你解决那些在你申请的工作中可能根本不会遇到的问题吗？比如，他们让你翻转二叉树，或者用 \(Python\) 统计重复元素。其实，他们本质上是在利用这样一个事实：在下面这个图中，相关性会通过叉形结构传递：
g = gr.Digraph()
g.edge("good programmer", "can invert a binary tree")
g.edge("good programmer", "good employee")
g.render(filename="good_programmer", format="png", view=True)## 'good_programmer.png'
招聘人员知道，优秀的程序员往往是表现出色的员工。但面试你时，他们不知道你是不是优秀的程序员。所以，他们会问你一个只有优秀程序员才能回答的问题。这个问题不一定与你申请的工作中会遇到的问题相关，它只是用来表明你是否是优秀的程序员。要是你能回答这个问题，你就很可能是优秀的程序员，这意味着你也很可能会是优秀的员工。
现在，假设招聘人员已经知道你是一名优秀的程序员。也许他们在之前的公司就认识你，或者你有令人印象深刻的学历。在这种情况下，知道你能否回答招聘流程中的问题，并不会为你是否会是一名优秀员工提供更多信息。用专业术语来说，一旦以 “是优秀程序员” 为条件，回答问题和成为优秀员工就是相互独立的。
更一般地讲，如果你有一个叉形结构，两个共享共同原因的变量是相关联的，但在以该共同原因为条件时，它们是相互独立的。或者说：
\(T \not\perp Y\)
但
\(T \perp Y \mid X\)
碰撞结构
唯一缺失的结构是 “碰撞结构”（没错，这是个专业术语 ）。当两个节点共享一个子节点，且它们之间没有直接关系时，就形成了不道德结构。换种说法，就是两个变量存在一个共同效应。这个共同效应常被称作 “碰撞器”，因为有两条箭头在此交汇：
g = gr.Digraph()
g.edge("Y", "X")
g.edge("T", "X")
g.edge("statistics", "job promotion")
g.edge("flatter", "job promotion")
g.render(filename="job_promotion", format="png", view=True)## 'job_promotion.png'
在碰撞结构中，两个父节点彼此独立。但要是以共同效应为条件，它们就会变得相关。比如，假设获得工作晋升有两种途径：要么擅长统计学，要么讨好上司。要是我不以你的工作晋升为条件（也就是不知道你会不会得到晋升），那么你的统计学水平和讨好上司的能力是相互独立的。换句话说，知道你统计学有多好，并不会让我了解到你讨好上司的本事如何。但另一方面，要是你确实得到了晋升，突然之间，知道你的统计学水平，就能让我推断出你讨好上司的水平。要是你统计学不好但却得到了晋升，那你很可能很擅长讨好上司。否则，你就不太可能获得晋升。反过来，要是你统计学很好，那么你不擅长讨好上司的可能性更大，因为统计学好本身就解释了你的晋升原因。这种现象有时被称为 “解释消除”，因为一个原因已经解释了结果，使得另一个原因的可能性降低。
一般规则是，以碰撞器为条件会开启关联路径，让变量变得相关；不以它为条件，路径就是关闭的。或者说：
\(T \perp Y\)
并且
\(T \not\perp Y \mid X\)
重要的是，如果你不以碰撞器为条件，而是以该碰撞器的某个效应（直接或间接的）为条件，也能开启相同的依赖路径。继续用我们的例子，假设获得工作晋升会大幅增加你的薪水，这样就有了下面这个图：
g = gr.Digraph()
g.edge("Y", "X1")
g.edge("T", "X1")
g.edge("X1", "X2")
g.node("X2", color="lightgrey", style="filled")
g.edge("statistics", "job promotion")
g.edge("flatter", "job promotion")
g.edge("job promotion", "high salary")
g.node("high salary", color="lightgrey", style="filled")
g.render(filename="high_salary", format="png", view=True)## 'high_salary.png'
在这个图中，即便你不以碰撞器为条件，而是以它的一个原因为条件，碰撞器的原因也会变得相关。比如，就算我不知道你是否得到晋升，但我知道你薪水很高，那么你的统计学知识和讨好上司的能力就会变得相关：具备其中一项，会降低你具备另一项的可能性。
关联传递速查表
了解了链式结构、叉形结构和不道德结构这三种结构后，你就能推导出关于图中独立性和关联传递的更通用规则。
一条路径被阻断，当且仅当：
它包含一个已被条件限定的非碰撞器。
它包含一个未被条件限定、且其所有被条件限定的后代也不存在的碰撞器。
\(图3-1\) 是一张关于图中依赖关系如何传递的速查表。
要是这些规则看起来有点晦涩难懂、难以把握，那现在我可以告诉你，幸运的是，你可以使用现成的算法来检查图中的两个变量是相互关联还是相互独立。为了把你所学的知识融会贯通，咱们来看最后一个例子，这样我就能向你展示如何编写相关代码了。
在 Python 中查询图形
看下面这个图形：
g = gr.Digraph(graph_attr={"rankdir": "LR"})
g.edge("C", "A")
g.edge("C", "B")
g.edge("D", "A")
g.edge("B", "E")
g.edge("F", "E")
g.edge("A", "G")
g.render(filename="check", format="png", view=True)## 'check.png'
很快，你会把这个图输入到一个 \(python\) 库中，这样就能轻松回答关于它的问题。但在这之前，为了巩固你刚学到的概念，先自己试着回答以下问题，当作练习：
\(D\) 和 \(C\) 是相关联的吗？
在给定 \(A\) 的条件下，\(D\) 和 \(C\) 是相关联的吗？
在给定 \(G\) 的条件下，\(D\) 和 \(C\) 是相关联的吗？
\(A\) 和 \(B\) 是相关联的吗？
在给定 \(C\) 的条件下，\(A\) 和 \(B\) 是相关联的吗？
\(G\) 和 \(F\) 是相关联的吗？
在给定 \(E\) 的条件下，\(G\) 和 \(F\) 是相关联的吗？
现在，要验证你回答得是否正确，你可以把那个图输入到 \(networkx\) 库的 \(DiGraph\) 中。\(networkx\) 是一个用于处理图形模型的库，有一系列实用算法可帮助你检查这个图：
import networkx as nx
model = nx.DiGraph([
("C", "A"),
("C", "B"),
("D", "A"),
("B", "E"),
("F", "E"),
("A", "G"),
])首先，来看 \(D\) 和 \(C\) 。它们构成了你之前见过的碰撞结构，其中 \(A\) 是碰撞因子。根据碰撞结构中独立性的规则，你知道 \(D\) 和 \(C\) 是相互独立的。你也知道，要是以碰撞因子 \(A\) 为条件，关联就会在它们之间开始传递。\(is\_d\_separator\) 方法能告诉你图中两个变量之间是否有关联传递（\(d\) 分离是表达图中两个变量独立性的另一种方式 ）。要以某个变量为条件，你可以把该变量添加到观测集合中。比如，要检查在给定 \(A\) 的条件下 \(D\) 和 \(C\) 是否相关，你可以使用 \(is\_d\_separator\) ，并传入第四个参数 \(z = \{ "A" \}\) ，像这样：
print("Are D and C dependent?")## Are D and C dependent?print(not(nx.is_d_separator(model, {"D"}, {"C"}, set())))## Falseprint("Are D and C dependent given A?")## Are D and C dependent given A?print(not(nx.is_d_separator(model, {"D"}, {"C"}, {"A"})))## Trueprint("Are D and C dependent given G?")## Are D and C dependent given G?print(not(nx.is_d_separator(model, {"D"}, {"C"}, {"G"})))## True接下来，注意 \(D\)、\(A\) 和 \(G\) 构成一条链。你知道关联性会在链中传递，所以 \(D\) 和 \(G\) 不独立。不过，要是以中间变量A为条件，就能阻断关联性的传递：
print("Are G and D dependent?")## Are G and D dependent?print(not(nx.is_d_separator(model, {"G"}, {"D"}, set())))## Trueprint("Are G and D dependent given A?")## Are G and D dependent given A?print(not(nx.is_d_separator(model, {"G"}, {"D"}, {"A"})))## False你需要复习的最后一种结构是叉形结构。你可以看到，\(A\)、\(B\) 和 \(C\) 构成一个叉形结构，其中 \(C\) 是 \(A\) 和 \(B\) 的共同原因。你知道关联性会通过叉形结构传递，所以 \(A\) 和 \(B\) 不相互独立。不过，要是以这个共同原因为条件，关联性的传递路径就会被阻断：
print("Are A and B dependent?")## Are A and B dependent?print(not(nx.is_d_separator(model, {"A"}, {"B"}, set())))## Trueprint("Are A and B dependent given C?")## Are A and B dependent given C?print(not(nx.is_d_separator(model, {"A"}, {"B"}, {"C"})))## False最后，咱们把所有内容整合起来，聊聊 \(G\) 和 \(F\)。关联性会在它们之间传递吗？从 \(G\) 开始看。你知道关联性会在 \(G\) 和 \(E\) 之间传递，因为它们处于一个叉形结构中。不过，关联性在碰撞因子 \(E\) 处会中断，这意味着 \(G\) 和 \(F\) 是相互独立的。但要是以 \(E\) 为条件，关联性就会开始通过碰撞器传递，路径开启，把 \(G\) 和 \(F\) 关联起来：
print("Are G and F dependent?")## Are G and F dependent?print(not(nx.is_d_separator(model, {"G"}, {"F"}, set())))## Falseprint("Are G and F dependent given E?")## Are G and F dependent given E?print(not(nx.is_d_separator(model, {"G"}, {"F"}, {"E"})))## True这太棒了。你不仅学习了图中的三种基本结构，还了解了如何使用现成的算法来检查图中的独立性。但这和因果推断有什么关系呢？是时候回到本章开头我们探讨的问题了。
回想一下，我们当时试图理解聘请昂贵的顶尖顾问对业务绩效的影响，我们将其描绘成了如下图形：
你可以用新学到的技能，来弄清楚在这个图里为什么关联不等于因果。注意，这个图里有一个叉形结构。所以，在聘请顾问和公司未来业绩之间，存在两种关联传递路径：一条是直接的因果路径，另一条是受共同原因干扰的非因果路径。后一种路径被称为 “后门路径” 。这个图中存在混淆性的后门路径，这表明观察到的聘请顾问和公司业绩之间的关联，不能单纯归因于因果关系。
理解关联如何通过非因果路径在图中传递，会让你在谈论关联和因果的区别时，表述得更精准。基于这个原因，结合图形模型的新视角，重新审视识别这一概念是很有价值的。
再谈识别
到目前为止，在没有随机化的情况下，我用来解释为何难以找到因果效应的理由是，处理组和未处理组无法相互比较。比如，聘请顾问的公司，通常比不聘请昂贵顾问的公司过往业绩更好。这就会产生你之前见过的那种偏差：
\(E[Y \mid T = 1] - E[Y \mid T = 0] = \underbrace{E[Y_1 - Y_0 \mid T = 1]}_{ATT} + \underbrace{\{E[Y_0 \mid T = 1] - E[Y_0 \mid T = 0]\}}_{BIAS}\)
既然你已经学习了因果图，就能更精准地描述这种偏差的性质，更重要的是，你能明白如何消除它。识别与图形模型中的独立性密切相关。要是你有一个图，描绘了处理、结果和其他相关变量之间的因果关系，你可以把识别看作是在这个图中分离出处理和结果之间因果关系的过程。在识别阶段，本质上你要关闭所有不想要的关联传递路径。
以顾问咨询的因果图为例。正如你之前所见，处理和结果之间有两条关联路径，但其中只有一条是因果路径。
要检查偏差，你可以构建一个和原图类似的因果图，但移除其中的因果关系。如果在这个图中，处理和结果仍有关联，那一定是因为存在非因果路径，这就表明存在偏差：
consultancy_sev = gr.Digraph(graph_attr={"rankdir": "LR"})
consultancy_sev.edge("profits_prev_6m", "profits_next_6m")
consultancy_sev.edge("profits_prev_6m", "consultancy")
consultancy_sev.render(filename="profits_prev_6m", format="png", view=True)## 'profits_prev_6m.png'
consultancy_model_severed = nx.DiGraph([
("profits_prev_6m", "profits_next_6m"),
("profits_prev_6m", "consultancy"),
# ("consultancy", "profits_next_6m"), # causal relationship removed
])
not(nx.is_d_separator(consultancy_model_severed, {"consultancy"}, {"profits_next_6m"}, set()))## True这些非因果的关联传递被称为 “后门路径” 。要识别T和Y之间的因果关系，你需要关闭这些路径，让仅剩下因果路径。在顾问咨询的例子中，你知道以共同原因（也就是公司的过往业绩）为条件，就能关闭那条路径，具体如下：
g_consultancy = gr.Digraph(graph_attr={"rankdir":"LR"})
g_consultancy.edge("profits_prev_6m", "profits_next_6m")
g_consultancy.edge("profits_prev_6m", "consultancy")
g_consultancy.edge("consultancy", "profits_next_6m")
g_consultancy.node("profits_prev_6m", color="lightgrey", style="filled")
g_consultancy.render(filename="profits_next_6m", format="png", view=True)## 'profits_next_6m.png'
CIA 与调整公式
你刚看到，以 \(profits\_prev\_6m\) 为条件，会阻断处理（聘请顾问）和结果（公司未来业绩）之间的非因果关联传递。因此，要是你观察一组过往业绩相似的公司，在组内对比聘请顾问和未聘请顾问的公司的未来业绩，差异就能完全归因于顾问。这在直觉上说得通，对吧？处理组（聘请顾问的公司）和未处理组未来业绩的差异，一是源于处理（聘请顾问）本身，二是源于聘请顾问的公司原本往往运营良好这一事实。要是你只对比运营状况同样良好的处理组和未处理组公司，第二种差异来源就会消失。
当然，和因果推断中的所有情况一样，你在这里做了一个假设。具体来说，你假设处理组和结果之间所有非因果关联的来源，都源于你能够测量并以之为条件的共同原因。这和你之前看到的独立性假设非常相似，但形式更弱：
\((Y_0, Y_1) \perp T \mid X\)
这个条件独立性假设\(（CIA）\)表明，如果你对比协变量 \(X\) 水平相同的单元（比如公司），它们的潜在结果平均而言是相同的。换句话说，要是观察协变量 \(X\) 取值相同的单元，处理就好像是随机分配的。
CIA 的其他名称
CIA 在很多因果推断研究中都存在，并且有很多其他名称，比如可忽略性\(（ignorability）\)、外生性 \(（exogeneity）\) 或可交换性 \(（exchangeability）\) 。
条件独立性假设还给出了一种从数据中的可观测数量识别因果效应的非常简单的方法。如果在由协变量 \(X\) 定义的组内，处理看起来就像随机分配的一样，那么你要做的就是在每个由 \(X\) 定义的组内对比处理组和未处理组，然后以组的规模为权重对结果取平均：
\(ATE = E_X [E[Y \mid T = 1] - E[Y \mid T = 0]]\)
\(\begin{align*} ATE &= \sum_{x} \{(E[Y \mid T = 1, X = x] - E[Y \mid T = 0, X = x]) P(X = x)\} \\ &= \sum_{x} \{E[Y \mid T = 1, X = x] P(X = x) - E[Y \mid T = 0, X = x] P(X = x)\} \end{align*}\)
这被称为调整公式或条件性原理。它表明，如果你对 \(X\) 进行条件设定（或控制），平均处理效应可以识别为处理组和控制组在组内差异的加权平均。同样，如果对 \(X\) 进行条件设定能阻断图中非因果路径的关联传递，像平均处理效应这样的因果量就变得可识别了，这意味着你可以从可观测数据中计算它。通过对混杂因素进行调整来关闭后门路径的过程，有个极具创意的名字叫后门调整。
正性假设
调整公式也凸显了正性的重要性。由于你要在协变量 \(X\) 的各个分组上，对处理和结果之间的差异取平均，所以必须确保，对于 \(X\) 的所有分组，处理组和对照组中都有一定数量的单元，否则差异无法定义。更正式地说，处理的条件概率需要严格大于 \(0\) 且小于 \(1\)，即 \(1 > P(T \mid X) > 0\) 。当正性假设不成立时，识别仍有可能实现，但这需要你进行有风险的外推。
正性假设的其他名称
由于正性假设在因果推断中也很常用，它也有很多其他名称，比如共同支撑\(（common \space support）\)或重叠 \(（overlap）\) 。
一个带数据的识别示例
由于这可能有点抽象，咱们结合一些数据看看实际情况。为简化示例，假设你收集了\(6\)家公司的数据，其中\(3\)家在过去\(6\)个月利润较低（\(100\) 万美元），另外\(3\)家利润较高。正如你推测的，高利润公司更有可能聘请顾问。\(3\)家高利润公司中有\(2\)家聘请了顾问，而\(3\)家低利润公司中只有\(1\)家聘请了顾问（要是样本量小让你困扰，就假设这里每个数据点实际代表\(10000\)家公司 ）：
df = pd.DataFrame(dict(
profits_prev_6m=[1.0, 1.0, 1.0, 5.0, 5.0, 5.0],
consultancy=[0, 0, 1, 0, 1, 1],
profits_next_6m=[1, 1.1, 1.2, 5.5, 5.7, 5.7],
))
df## profits_prev_6m consultancy profits_next_6m
## 0 1.0 0 1.0
## 1 1.0 0 1.1
## 2 1.0 1 1.2
## 3 5.0 0 5.5
## 4 5.0 1 5.7
## 5 5.0 1 5.7要是你直接对比聘请顾问的公司和没聘请顾问的公司的 \(profits\_next\_6m\) （未来六个月利润），你会发现利润差异为 \(166\) 万美元：
(df.query("consultancy==1")["profits_next_6m"].mean()
- df.query("consultancy==0")["profits_next_6m"].mean())## 1.666666666666667但你更清楚。这并非顾问咨询对公司业绩的因果效应，因为在过去业绩表现更好的公司，在聘请顾问的组里占比过高。要得到顾问效应的无偏估计，你需要观察过往业绩相近的公司。如你所见，这样得出的结果会更合理：
avg_df = (df
.groupby(["consultancy", "profits_prev_6m"])
["profits_next_6m"]
.mean())
avg_df.loc[1] - avg_df.loc[0] ## profits_prev_6m
## 1.0 0.15
## 5.0 0.20
## Name: profits_next_6m, dtype: float64要是你对这些效应取加权平均（权重为每个组的规模），最终会得到平均处理效应的无偏估计。在这里，因为两个组规模相同，这就是简单平均，得出的平均处理效应为\(175000\)。所以，如果你是一名正在决定是否聘请顾问的经理，且看到了上述数据，你可以得出结论：顾问对未来利润的影响约为 \(17.5\) 万美元。当然，要这么做，你得援引条件独立性假设。也就是说，你得假设过往业绩是聘请顾问和未来业绩的唯一共同原因。
你刚刚完整走了一遍这样的示例：把你对因果机制的信念编码成一个图，然后利用这个图找出为了估计平均处理效应需要以哪些变量为条件，即便没有对处理进行随机化。接着，你结合一些数据，看到了实际操作是怎样的 —— 按照调整公式，在假设条件独立性的情况下，对平均处理效应进行了估计。这里用到的工具相当通用，能为你未来遇到的诸多因果问题提供思路。不过，我觉得这还没完。有些图形结构，以及它们所导致的偏差，比其他的要常见得多。仔细研究这些结构是值得的，这样你就能开始体会到在因果推断之旅中前方会遇到的困难。
后门调整并非识别因果效应的唯一可行策略。人们可以利用对因果机制的了解，通过前门来识别因果效应，即便存在无法测量的共同原因：
g = gr.Digraph(graph_attr={"rankdir":"LR"})
g.edge("U", "T")
g.edge("U", "Y")
g.edge("T", "M")
g.edge("M", "Y")
g.render(filename='mediator', format='png', view=True)## 'mediator.png'
采用这种策略时，你必须能够识别处理对中介变量\(（mediator）\)的效应，以及该中介变量对结果的效应。这样一来，处理对结果的因果效应的识别，就变成了这两种效应的结合。不过，在科技行业，很难找到能让这样的图合理成立的应用场景，这就是前门调整不太流行的原因。
混杂偏差
偏差的首要重要来源是混杂。这就是我们到目前为止一直在讨论的偏差，现在给它正式命名。当存在开放的后门路径，使得关联能通过非因果方式传递时，就会发生混杂，通常是因为处理和结果存在共同原因。比如，假设你在人力资源部门工作，想了解新的管理培训项目是否能提高员工的敬业度。但由于培训是自愿参加的，你认为只有原本工作就很出色的经理会参加该项目，而最需要培训的经理反而不会参加。当你衡量参加培训的经理所带领团队的敬业度时，会发现其远高于未参加培训经理所带领团队的敬业度。但很难判断这种差异中有多少是因果效应导致的。因为处理和结果存在共同原因，即便没有因果效应，它们也会呈现出相关变化。
要识别这种因果效应，你需要关闭处理和结果之间所有的后门路径。若能做到，剩下的效应就只有直接效应 \(T \to Y\) 。在我们的例子中，你可以想办法控制经理参加培训前的素质。这样一来，结果的差异就只能归因于培训，因为处理组和对照组中经理培训前的素质会保持一致。
简而言之，要调整混杂偏差，你需要针对处理和结果的共同原因进行调整：
g = gr.Digraph(graph_attr={"rankdir":"LR"})
g.edge("X", "T")
g.edge("X", "Y")
g.edge("T", "Y")
g.edge("Manager Quality", "Training"),## (None,)g.edge("Manager Quality", "Engagement"),## (None,)g.edge("Training", "Engagement")
g.render(filename='confounder', format='jpg', view=True)## 'confounder.jpg'
遗憾的是，并非总能对所有共同原因进行调整。有时，存在未知的原因，或者虽已知但无法测量的原因。经理素质就是这样的例子。尽管付出了诸多努力，我们仍未搞清楚如何衡量管理素质。要是你无法观测到经理素质，就无法以它为条件，培训对敬业度的效应也就无法识别。
替代混杂
在某些情况下，由于存在无法测量的混杂因素，你无法关闭所有后门路径。在下面这个例子中，还是经理素质导致经理选择参加培训，同时影响团队敬业度。所以，处理（培训）和结果（团队敬业度）之间的关系存在混杂。但在这种情况下，你无法以该混杂因素为条件，因为它无法测量。此时，由于混杂偏差，处理对结果的因果效应无法识别。
不过，你有其他可测量的变量，能够作为混杂因素（经理素质）的替代指标。这些变量不在后门路径中，但对其进行控制有助于减少偏差（尽管无法消除偏差 ）。这些变量有时被称为替代混杂因素。
在这个例子中，你无法衡量经理素质，但可以衡量其部分成因，比如经理的任期或教育水平；还能衡量其部分结果，比如团队的人员流动率或绩效。对这些替代变量进行控制虽不足以消除偏差，但确实会有帮助：
g = gr.Digraph()
g.edge("X1", "U")
g.edge("U", "X2")
g.edge("U", "T")
g.edge("T", "Y")
g.edge("U", "Y")
g.edge("Manager Quality", "Team's Attrition")
g.edge("Manager Quality", "Team's Past Performance")
g.edge("Manager's Tenure", "Manager Quality")
g.edge("Manager's Education Level", "Manager Quality")
g.edge("Manager Quality", "Training")
g.edge("Training", "Engagement")
g.edge("Manager Quality", "Engagement")
g.render(filename='surrogate_confounding', format='png', view=True)## 'surrogate_confounding.png'
再谈随机化
在许多重要且高度相关的研究问题中，混杂因素是个大问题，因为你永远无法确定自己是否已经控制了所有混杂因素。但如果你计划主要在行业中运用因果推断，我有个好消息要告诉你。在行业里，你大多关注的是了解公司能够控制的事物的因果效应，比如价格、客户服务和营销预算，这样你就能对它们进行优化。在这些情况下，要知晓混杂因素是什么相当容易，因为企业通常清楚自己在分配处理时用到了哪些信息。不仅如此，即便企业不清楚，对处理变量进行干预几乎也总是可行的。这正是 \(A/B\) 测试的关键所在。当你对处理进行随机化时，你可以从一个存在不可观测混杂因素的图，转变为一个处理的唯一成因是随机性的图：
g = gr.Digraph(graph_attr={"rankdir":"LR"})
g.edge("rnd", "T")
g.edge("T", "Y")
g.edge("U", "Y")
g.render(filename='rnd', format='png', view=True )## 'rnd.png'
因此，除了尝试弄清楚为识别效应需要以哪些变量为条件之外，你还应该问问自己，能采取哪些可能的干预措施，将图形转变为感兴趣的因果量可识别的图形。
当存在不可观测的混杂因素时，并非就毫无办法了。在第四部分\(（Part \space IV）\)，我会介绍可利用数据中的时间结构来应对不可观测混杂因素的方法。第五部分\(（Part \space V）\)会介绍出于同样目的的工具变量的运用。
敏感性分析与部分识别
当你无法测量所有共同原因时，与其直接放弃，将讨论从 “我是否测量了所有混杂因素” 转变为 “未测量的混杂因素需要多强，才能显著改变我的分析结果”，往往会更有成效。这就是敏感性分析的核心思想。若想全面回顾这个主题，我建议你查阅 Cinelli 和 Hazlett 撰写的论文《Making Sense of Sensitivity: Extending Omitted Variable Bias》（理解敏感性：拓展遗漏变量偏差 ）。
此外，即便你关注的因果量无法被点识别，你仍然可以利用可观测数据为其确定边界。这个过程被称为部分识别，目前是一个活跃的研究领域。
选择偏差
要是你觉得混杂偏差已经像因果推断鞋里的一颗狡猾小石子，那等你了解选择偏差后就知道厉害了。混杂偏差是在你未对处理和结果的共同原因进行控制时出现，而选择偏差更多与对共同效应和中介变量进行条件设定相关。
偏差术语
文献中对于偏差的命名并没有统一共识。比如，经济学家往往把各类偏差都称作选择偏差。相反，有些科学家喜欢把我所说的选择偏差进一步细分为碰撞机偏差\(（collider \space bias）\)和中介偏差\(（mediator \space bias）\)。我会采用 Miguel A. Hernán 和 James M. Robins 所著《Causal Inference: What If》（查普曼与霍尔 / CRC 出版社 ）一书中的术语。
以碰撞因子为条件
假设你在一家软件公司工作，想要估计刚刚推出的新功能的影响。为避免出现任何混杂偏差，你对该功能进行了随机推广：10% 的客户被随机选中使用新功能，其余客户则不使用。你想知道这个功能是否让客户更开心、更满意。由于满意度无法直接衡量，你用净推荐值\(（NPS）\)作为其替代指标。为测量 \(NPS\)，你向推广组（处理组）和对照组的客户发送调查问卷，询问他们是否会推荐你的产品。结果出来后，你发现使用新功能且参与 \(NPS\) 调查的客户，其 \(NPS\) 得分高于未使用新功能但也参与调查的客户。你能说这种差异完全是由新功能对 \(NPS\) 的因果效应导致的吗？要回答这个问题，你应该从呈现这种情况的因果图入手：
g = gr.Digraph(graph_attr={"rankdir": "LR"})
g.edge("T", "S")
g.edge("T", "Y")
g.edge("Y", "S")
g.node("S", color="lightgrey", style="filled")
g.edge("RND", "New Feature"),## (None,)g.edge("New Feature", "Customer Satisfaction"),## (None,)g.edge("Customer Satisfaction", "NPS"),## (None,)g.edge("Customer Satisfaction", "Response"),## (None,)g.edge("New Feature", "Response"),## (None,)g.node("Response", "Response", color="lightgrey", style="filled")
g.render(filename='nps', format='png', view=True )## 'nps.png'
长话短说，遗憾的是，答案是否定的。问题在于，你只能测量那些参与 \(NPS\) 调查的客户的 \(NPS\)。你在估计处理组和对照组差异的同时，也以参与 \(NPS\) 调查的客户为条件进行了分析。尽管随机化能让你通过处理组和对照组的结果差异识别平均处理效应，但一旦以共同效应为条件，你就引入了选择偏差。要明白这一点，你可以重新绘制因果图，删除新功能到客户满意度的因果路径，这也会关闭到 \(NPS\) 的直接路径。然后，在以客户回应为条件的情况下，检查 \(NPS\) 是否仍与新功能相关。你会发现确实相关，这意味着关联通过非因果路径在两个变量间传递，而这正是偏差的含义：
nps_model = nx.DiGraph([
("RND", "New Feature"),
# ("New Feature", "Customer Satisfaction"),
("Customer Satisfaction", "NPS"),
("Customer Satisfaction", "Response"),
("New Feature", "Response"),
])
not(nx.d_separated(nps_model, {"NPS"}, {"New Feature"}, {"Response"}))## True另见
选择偏差下的因果识别非常隐蔽。这种删除因果路径并检查处理和结果是否仍相关的方法，并非在所有选择偏差情况中都有效。遗憾的是，在撰写本文时，我不清楚有任何 \(Python\) 库可处理存在选择偏差的图。不过你可以查看 \(DAGitty\)，它可在浏览器中运行，且有用于选择偏差下识别的算法 。
为了培养你对这种偏差的直觉，咱们重新赋予你神一般的能力，假设你能洞察反事实结果的世界。也就是说，你能看到所有客户（即便那些没回复调查的客户）在控制状态下的净推荐值\(NPS\)，即 \(nps\_0\)，以及在处理状态下的净推荐值，即 \(nps_1\) 。咱们还要以一种能知晓真实效应的方式模拟数据。这里，新功能使 \(NPS\) 提高 \(0.4\)（以任何商业标准来看，这数值都很高，但为方便举例，先这么设定 ）。咱们也假设，就像之前的因果图里呈现的那样，新功能和客户满意度都会增加回复 \(NPS\) 调查的可能性。凭借测量反事实的能力，要是你按处理组和对照组汇总数据，会看到这样的情况：
np.random.seed(2)
n = 100000
new_feature = np.random.binomial(1, 0.5, n)
satisfaction_0 = np.random.normal(0, 0.5, n)
satisfaction_1 = satisfaction_0 + 0.4
satisfaction = new_feature*satisfaction_1 + (1-new_feature)*satisfaction_0
nps_0 = np.random.normal(satisfaction_0, 1)
nps_1 = np.random.normal(satisfaction_1, 1)
nps = new_feature*nps_1 + (1-new_feature)*nps_0
responded = (np.random.normal(0 + new_feature + satisfaction, 1) > 1).astype(int)
tr_df = pd.DataFrame(dict(new_feature=new_feature,
responded=responded,
nps_0=nps_0,
nps_1=nps_1,
nps=nps))
tr_df_measurable = pd.DataFrame(dict(new_feature=new_feature,
responded=responded,
nps_0=np.nan,
nps_1=np.nan,
nps=np.where(responded, nps, np.nan)))
tr_df.groupby("new_feature").mean()## responded nps_0 nps_1 nps
## new_feature
## 0 0.183715 -0.005047 0.395015 -0.005047
## 1 0.639342 -0.005239 0.401082 0.401082
首先，注意到使用新功能的客户中，有 \(63\%\) 回复了 \(NPS\) 调查，而对照组中只有 \(18\%\) 回复了调查。接下来，看处理组和对照组的数据行，从 \(NPS_0\) 到 \(NPS_1\) 有 \(0.4\) 的增长。这就说明，新功能对两组客户的效应都是 \(0.4\)。最后，注意到处理组\(（new\_feature=1）\)和对照组\(（new_feature=0）\)之间的 \(NPS\) 差异约为 \(0.4\)。同样，要是你能看到未回复 \(NPS\) 调查客户的 \(NPS\)，直接对比处理组和对照组就能得到真实的平均处理效应。
当然，现实中，你看不到 \(NPS_0\) 和 \(NPS_1\) 这两列数据。你也无法这样查看 \(NPS\) 列，因为你只有回复调查客户的 \(NPS\) 数据（对照组数据行的 \(18\%\) 和处理组数据行的 \(63\%\) 回复了调查 ）：
tr_df_measurable.groupby("new_feature").mean().assign(**{"nps": np.nan})## responded nps_0 nps_1 nps
## new_feature
## 0 0.183715 NaN NaN NaN
## 1 0.639342 NaN NaN NaN
要是你按回复者进一步细分分析，就能看到 \(Response=1\)（回复调查）客户的 \(NPS\)。但注意看，该组中处理组和对照组的差异不再是 \(0.4\)，而只有大约一半\(0.22\)。这是怎么回事呢？这全是选择偏差导致的：
tr_df_measurable.groupby(["responded", "new_feature"]).mean()## nps_0 nps_1 nps
## responded new_feature
## 0 0 NaN NaN NaN
## 1 NaN NaN NaN
## 1 0 NaN NaN 0.314073
## 1 NaN NaN 0.536106
把不可观测的量加回来，你就能明白是怎么回事了（这里聚焦于回复调查的群体 ）：
tr_df.groupby(["responded", "new_feature"]).mean()## nps_0 nps_1 nps
## responded new_feature
## 0 0 -0.076869 0.320616 -0.076869
## 1 -0.234852 0.161725 0.161725
## 1 0 0.314073 0.725585 0.314073
## 1 0.124287 0.536106 0.536106
起初，处理组和对照组在基线满意度 \(Y_0\) 方面是可比的。但一旦以回复调查的客户为条件，处理组的基线满意度就更低，即 \(E[Y_0|T = 0, R = 1] > E[Y_0|T = 1, R = 1]\) 。这意味着，一旦以回复调查的客户为条件，处理组和对照组之间简单的平均值差异并不能识别平均处理效应：
\(\begin{align*} E[Y|T = 1, R = 1] - E[Y|T = 0, R = 1] &= E[Y_1 - Y_0|R = 1]_{\text{ATE}} \\ & \quad + E[Y_0|T = 0, R = 1] - E[Y_0|T = 1, R = 1]_{\text{SelectionBias}} \end{align*}\)
要是结果（客户满意度）会影响回复率，那个偏差项就不会为零。因为满意的客户更有可能回答 NPS 调查，这种情况下就无法识别因果效应。如果处理（新功能）能提高满意度，那么对照组会包含更多基线满意度高于处理组的客户。这是因为处理组会有原本就满意（基线满意度高）的客户，以及原本基线满意度低，但因处理变得更满意并回复调查的客户。
另见
选择偏差是个比我在本章中所能充分阐述的复杂得多的主题。比如，仅仅通过以结果的某个效应为条件，你就可能产生选择偏差，即便该效应并未与处理共享。这种情况被称为 “虚拟对撞机\(（virtual \space collider）\)”。要了解更多相关内容以及其他知识，我强烈建议你查阅 Carlos Cinelli 等人撰写的论文《A Crash Course in Good and Bad Controls》（良好控制与不良控制速成课程 ）。这篇论文涵盖了本章涉及的所有内容，还有更多拓展。而且论文语言通俗易懂，便于阅读。
调整选择偏差
遗憾的是，校正选择偏差绝非易事。在我们一直在讨论的例子中，即便进行了随机对照试验，平均处理效应仍无法识别，原因很简单：一旦以回复调查的客户为条件，你就无法阻断新功能和客户满意度之间非因果的关联传递。要取得进展，你需要做出进一步假设，而这正是图形模型开始发挥优势之处。它能让你非常明确、清晰地阐述这些假设。
比如，你需要假设结果不会导致选择。在我们的例子中，这意味着客户满意度不会使客户更有可能或更不可能回答调查。相反，得有其他可观测变量（一个或多个）同时导致选择和结果。比如，可能唯一导致客户回复调查的因素是他们在应用里花费的时间和新功能。这种情况下，处理组和对照组之间的非因果关联会通过应用使用时长传递：
g = gr.Digraph()
g.edge("U", "X")
g.edge("X", "S")
g.edge("U", "Y")
g.edge("T", "Y")
g.edge("T", "S")
g.node("S", color="lightgrey", style="filled")
g.edge("New Feature", "Customer Satisfaction"),## (None,)g.edge("Unknown Stuff", "Customer Satisfaction"),## (None,)g.edge("Unknown Stuff", "Time in App"),## (None,)g.edge("Time in App", "Response"),## (None,)g.edge("New Feature", "Response"),## (None,)g.node("Response", "Response", color="lightgrey", style="filled")
g.render(filename='adjust_bias', format='png', view=True)## 'adjust_bias.png'
只有专业知识才能判断这个假设的强度。但如果该假设正确，一旦你控制了在应用程序中的使用时长，新功能对满意度的效应就可以识别了。
再次强调，你在这里应用的是调整公式。你只是依据 \(X\) 将数据划分为不同组，让处理组和对照组在这些分组内具有可比性。然后，你可以用每个组的规模作为权重，计算处理组和对照组组内比较的加权平均值。只不过现在，你在做这些的同时，也以选择变量为条件：
\(ATE = \sum_{x} \left\{ \left( E[Y \mid T = 1, R = 1, X] - E[Y \mid T = 0, R = 1, X] \right) P(X \mid R = 1) \right\}\)
一般来说，要调整选择偏差，你必须对导致选择的所有因素进行调整，还必须假设结果和处理都不会直接导致选择，也不会与选择存在隐藏的共同原因。例如，在下面的因果图中，存在选择偏差，因为以 \(S\) 为条件会开启 \(T\) 和 \(Y\) 之间的非因果关联路径：
g = gr.Digraph(graph_attr={"rankdir": "LR"})
g.edge("X1", "U")
g.edge("U", "X2")
g.edge("X5", "S")
g.edge("U", "Y",style="dashed")
g.edge("U", "S",style="dashed")
g.edge("U", "X3")
g.edge("X3", "S")
g.edge("Y", "X4")
g.edge("X4", "S")
g.edge("T", "X5")
g.edge("T", "Y")
g.edge("T", "S", style="dashed")
g.node("S", color="lightgrey", style="filled")
g.render(filename='s', format='png', view=True)## 's.png'
你可以通过调整可测量的、能解释选择的变量 \(X3\)、\(X4\) 和 \(X5\)，关闭其中两条路径。然而，有两条路径你无法关闭（用虚线显示 ）：\(Y \to S \leftarrow T\) 和 \(T \to S \leftarrow U \to Y\) 。这是因为处理会直接导致选择，且结果与选择存在隐藏的共同原因。你可以通过进一步以\(X2\) 和 \(X1\) 为条件，减轻最后这条路径带来的偏差，因为它们解释了\(U\)中的部分变异，但无法完全消除偏差。
这个因果图反映了在选择偏差方面，你会遇到的更合理的情况，就像我们刚才用作例子的回应偏差。在这些情况下，你能做的最佳举措就是以解释选择的变量为条件。这会减少偏差，但无法消除，因为如你所见：\(1）\) 存在导致选择的因素，你不知道或无法测量；\(2）\)结果或处理可能直接导致选择。
实际案例：生存分析中的隐藏偏差
生存分析出现在许多涉及事件持续时间或到某事件发生时间的商业应用中。例如，银行非常想了解贷款规模（贷款金额）如何增加客户拖欠该贷款的可能性。考虑一笔 \(3\) 年期贷款。客户可能在第 \(1\) 年、第 \(2\) 年或第 \(3\) 年拖欠，也可能完全不拖欠。银行的目标是了解贷款金额如何影响 \(P(Default|yr = 1)\)、\(P(Default|yr = 2)\) 和 \(P(Default|yr = 3)\) 。这里，为简单起见，假设银行已将贷款金额随机化。注意，只有在第 \(1\) 年存活（未拖欠 ）的客户才会在第 \(2\) 年被观测到，只有在第 \(1\) 年和第 \(2\) 年都存活的客户才会在第 \(3\) 年被观测到。这种选择导致只有第 \(1\) 年贷款规模的效应可识别。
直观来看，即便贷款金额是随机分配的，也仅在第 \(1\) 年能保持这种随机性。在那之后，如果贷款金额增加了拖欠可能性，拖欠风险较低的客户在高贷款金额区间会过度代表。他们的风险必须足够低，以抵消大额贷款金额导致的风险增加；否则，他们会在第 \(1\) 年拖欠。如果偏差过于极端，甚至可能看起来像是大额贷款导致第 \(1\) 年之后的年份风险降低，这毫无道理。
解决这种选择偏差问题的一个简单办法是聚焦累积结果（生存情况），即 \(Y|time > t\) ，而非年度结果（风险率），即 \(Y|time = t\) 。比如，尽管你无法识别贷款金额对第 \(2\) 年拖欠的影响 \(P(Default|yr = 2)\) ，但你可以轻松识别贷款金额对截至第 \(2\) 年拖欠的影响 \(P(Default|yr \leq 2)\) ：
g = gr.Digraph(graph_attr={"rankdir": "LR"})
g.edge("loan amount", "Default at yr=1")
g.edge("Default at yr=1", "Default at yr=2")
g.edge("Default at yr=2", "Default at yr=3")
g.edge("U", "Default at yr=1")
g.edge("U", "Default at yr=2")
g.edge("U", "Default at yr=3")
g.node("Default at yr=1", color="lightgrey", style="filled")
g.node("Default at yr=2", color="darkgrey", style="filled")
g.render(filename='default', format='png', view=True)## 'default.png'
我也不想让你产生错误的想法，即认为只要控制所有导致选择的因素就是个好主意。在下面的因果图中，以 \(X\) 为条件会开启一条非因果路径 \(Y \to X \leftarrow T\) ：
g = gr.Digraph(graph_attr={"rankdir":"LR"})
g.edge("Y", "X")
g.edge("T", "X")
g.edge("T", "Y")
g## <graphviz.graphs.Digraph object at 0x000001C80EA19F10>
以中介变量为条件
到目前为止，所讨论的选择偏差是由不可避免地对进入某一人群的选择（你被迫以回复调查的人群为条件）导致的，但你也可能无意中造成选择偏差。例如，假设你在人力资源部门工作，想要弄清楚是否存在性别歧视，即资质相当的男性和女性薪酬是否不同。为了进行这样的分析，你可能会考虑控制资历级别；毕竟，你想要比较资质相当的员工，而资历似乎是一个很好的替代指标。换句话说，你认为如果处于相同职位的男性和女性薪酬不同，你就有证据表明公司存在性别薪酬差距。
这种分析存在的问题是，因果图可能大致如下：
g = gr.Digraph(graph_attr={"rankdir": "LR"})
g.edge("T", "M")
g.edge("T", "Y")
g.edge("M", "Y")
g.node("M", color="lightgrey", style="filled")
g.edge("woman", "seniority")
g.edge("woman", "salary")
g.edge("seniority", "salary")
g.node("seniority", color="lightgrey", style="filled")
g.render(filename='seniority', format='png', view=True)## 'seniority.png'
资历级别是处理（女性）与薪酬之间路径的中介变量（mediator）。直观来看，男女薪酬差异有一个直接原因（直接路径：\(woman \to salary\) ），还有一个间接原因，通过资历传递（间接路径：\(woman \to seniority \to salary\) ）。这个因果图说明，女性遭受歧视的一种方式是，不太可能晋升到更高资历级别。男女薪酬差异，一部分是相同资历级别下的薪酬差异，另一部分是资历级别的差异。简单说，\(woman \to seniority \to salary\) 这条路径也是处理与结果之间的因果路径，分析时不应阻断它。如果控制资历级别来比较男女薪酬，你只能识别直接歧视（\(woman \to salary\) ） 。
值得一提的是，以中介变量节点的后代为条件也会引发偏差。这种选择不会完全阻断因果路径，但会部分阻断它：
g = gr.Digraph(graph_attr={"rankdir": "LR"})
g.edge("T", "M")
g.edge("T", "Y")
g.edge("M", "Y")
g.edge("M", "X")
g.node("X", color="lightgrey", style="filled")
g.render(filename='m', format='png', view=True)## 'm.png'
关键要点
在本章中，你主要聚焦于因果推断的识别环节。目标是学习如何使用图形模型，清晰呈现你所做的假设，以及了解这些假设会产生何种关联（因果或非因果关联 ）。为此，你得学习关联在图形中是如何传递的。这份速查表很好地总结了这些结构，所以建议你随时查阅：
接着，你了解到识别就相当于在图形中把因果关联流与非因果关联流区分开来。你可以通过对一些变量进行调整（控制条件），甚至对图形进行干预（比如在随机实验中那样），来阻断非因果的关联路径。像
networkx
这类贝叶斯网络软件在这里特别有用，因为它能在你检查图形中两个节点是否相连时提供帮助。例如，要检查是否存在混杂偏差，你可以简单地移除图形中的因果路径，然后查看处理节点和结果节点在移除该路径后是否仍然相连。如果相连，那就存在需要阻断的后门路径。
最后，你研究了因果问题中两种非常常见的偏差结构。当处理和结果存在共同原因时，就会出现混杂偏差。这个共同原因会形成一个叉形结构，从而在处理和结果之间产生非因果的关联流：
g = gr.Digraph(graph_attr={"rankdir":"LR", "ratio": "0.3"})
g.edge("U", "T")
g.edge("U", "Y")
g.edge("T", "Y")
g.render(filename='confounder', format='png', view=True)## 'confounder.png'
要修正混杂偏差，你需要尝试对共同原因进行调整，直接调整或者借助代理变量调整。这催生了调整公式的思路：
\(ATE = \sum_{x} \left\{ \left( E[Y \mid T = 1, X = x] - E[Y \mid T = 0, X = x] \right) P(X = x) \right\}\)
以及条件独立假设，该假设指出，如果在变量 \(X\) 的分组内，处理分配近似随机，那么你可以通过对 \(X\) 进行条件控制来识别因果量。
或者，如果你能对处理节点进行干预，混杂就会容易处理得多。例如，如果你设计一个随机实验，你会创建一个新的因果图，其中指向处理的所有箭头都被删除，这实际上就消除了混杂偏差。
你还了解了选择偏差，当你以处理和结果之间的共同效应（或共同效应的后代）为条件，或者以中介节点（或中介节点的后代）为条件时，就会出现选择偏差。选择偏差极其危险，因为它不会因实验而消失。更糟糕的是，它可能非常违反直觉，难以察觉：
g = gr.Digraph(graph_attr={"rankdir":"LR"})
g.edge("T", "M")
g.edge("M", "Y")
g.edge("T", "Y")
g.edge("T", "S")
g.edge("Y", "S")
g.node("M", color="lightgrey", style="filled")
g.node("S", color="lightgrey", style="filled")
g.render(filename='m_s', format='png', view=True)## 'm_s.png'
再次强调，了解因果图就像学习一门新语言。你会通过反复观察它、尝试运用它，来掌握大部分相关知识。