因果面板数据模型的矩阵补全方法
👇 连享会 · 推文导航 | www.lianxh.cn
🍎 Stata：Stata基础 | Stata绘图 | Stata程序 | Stata新命令 📘 论文：数据处理 | 结果输出 | 论文写作 | 数据分享 💹 计量：回归分析 | 交乘项-调节 | IV-GMM | 时间序列 | 面板数据 | 空间计量 | Probit-Logit | 分位数回归 ⛳ 专题：SFA-DEA | 生存分析 | 爬虫 | 机器学习 | 文本分析 🔃 因果：DID | RDD | 因果推断 | 合成控制法 | PSM-Matching 🔨 工具：工具软件 | Markdown | Python-R-Stata 🎧 课程：公开课-直播 | 计量专题 | 关于连享会
连享会 · 五一论文班
作者：王文茂 (中山大学)
邮箱：wangwm8@mail2.sysu.edu.cn
编者按：本文主要摘译自下文，特此致谢！
Source: Athey S, Bayati M, Doudchenko N, et al. Matrix completion methods for causal panel data models[J]. Journal of the American Statistical Association, 2021, 116(536): 1716-1730. -PDF-
目录
1. 引言
2. 基本模型设置
3. 矩阵形式和分类
3.1 数据缺失模式
3.2 瘦矩阵和胖矩阵
3.3 水平回归和垂直回归
3.4 固定效应和因子模型
4. 核范数最小化估计下的矩阵补全
4.1 估计量
4.2 计算估计量
5. 水平和垂直回归的关系
6. 两个图示例证
7. 一般化扩展
8. 参考资料
9. 相关推文
温馨提示： 文中链接在微信中无法生效。请点击底部「阅读原文」。或直接长按/扫描如下二维码，直达原文：
1. 引言
相关性并不意味着因果关系。在进行因果推断时，我们总是希望能够找到一个合理的反事实对照组。已有文献在估计面板数据的平均处理效应 (ATE) 时，为了得到该对照组，采取了非混淆性和合成控制两种方法。其中，
非混淆性可理解为，使用与前一时期观察结果值相似的对照单元观察结果，来推算处理单元丢失的潜在对照结果； 合成控制可理解为，使用对照单元的加权平均结果，来计算处理单元的潜在对照结果，并选择权重，以便对照单元的加权滞后结果与处理单元的滞后结果相匹配。
本文在前人的基础上提出一种新的方法来解决这一问题。借鉴关于因子模型和交互固定效应的计量经济学文献，以及关于矩阵补全的计算机科学和统计学文献，作者提出核范数矩阵补全估计法。在基于真实数据的模拟中，证明该方法的性能优于现有的方法。
本文的贡献主要有以下三方面：
一是推广了矩阵补全文献中的结果。修改了矩阵补全和因子模型文献中的估计量 (提出核范数矩阵补全估计法)，以考虑非正则化的单位和时间固定效应； 二是证明了非混淆性和合成控制，以及核范数最小化方法都可以看作是基于矩阵分解的矩阵补全方法。它们都具有基于 Fröbenius 范数的目标函数。根据这个共同的目标函数，非混淆性和合成控制方法对矩阵分解中的因素施加了不同的限制。核范数最小化方法不施加任何限制，而是使用正则化来表示估计量； 三是将该方法应用于两个真实数据集。人为地指定某些单元和时间段缺失的结果，然后比较不同补全估计量的表现。核范数矩阵补全估计量与另外的方法相比，在一系列情况下表现更为良好。
2. 基本模型设置
在 个时期内观察到 个单位的面板数据中，每个时期每个单位都有两个潜在的结果：(0) 和 。 表示该单元被处理，否则 。对于每个单元和时期，我们观察结果为：
文章集中于估计被处理组的平均效应 (ATE)：
要估计这样的平均处理效应，一种方法是计算缺失的潜在结果。由于已经观察到了 的所有相关值，我们只需要为 的受试者计算其相对应 矩阵的缺失值。这样就可以使 拥有与其相对应的对照项 ，从而能够计算其效应。
如下所示，、 均是 的矩阵。同时，将 定义为对应 中 对应的 中的项集合，将 定义为对应 中 对应的 中的项集合。文章主要就是找到关于在 中缺失值的统计问题 (即公式中的 )，一旦这些值被输入，我们就可以估计平均因果效应 。
其中，
3. 矩阵形式和分类
本节中将讨论矩阵 和 的一些特殊形式并分类。矩阵补全文献主要关注 是完全随机的情况，如上节所示的情形，并且 和 的维度都很大。首先，考虑缺失数据的模式，即 的分布不同于完全随机的分布。其次，考虑矩阵 的不同形状，其中维度 和 的相对大小可能非常不同。第三，考虑计量经济学文献中的一些具体分析，这些分析侧重于缺失数据模式和矩阵形状的特定组合。
3.1 数据缺失模式
数据缺失模式主要分为两种，块状结构和交错处理结构。其中，块状结构如下所示：
块状结构有两种特殊情况。在非混淆性 (Imbens 和 Rubin，2015) 下，估计平均处理效应的文献大多集中在 的情况下，因此仅有的处理单元在最后一个周期，称为单处理周期块结构 (Single-TreatedPeriod-Block Structure)。相反，合成控制文献 (Abadie 等，2010；Abadie，2019) 主要关注具有从时段 开始的多个时段处理的单个处理单元的情况，称为单处理单元块结构，如下所示：
对于仅仅缺少单个 的情况，上面提到的两种方法均适用：
交错处理结构是指不同的单位第一次接触处理的时间可能不同 (Athey 和 Imbens，2018；Shaikh 和 Toulis，2019)，具体如下所示：
3.2 瘦矩阵和胖矩阵
第二种分类涉及矩阵 的形状。由于 的相对大小变化，矩阵的形状就会产生不同。例如，
当 时， 是瘦矩阵：
当 时， 是胖矩阵：
当 时， 是近似正方形矩阵：
3.3 水平回归和垂直回归
依据数据缺失模式和矩阵形状分类的两种特殊组合值得特别关注，因为它们是很多文献的关注点。一是水平回归和非混淆性文献。非混淆性文献 (Rosenbaum 和 Rubin，1983；Rubin，2006；Imbens 和 Wooldridge，2009；Abadie 和 Cattneo，2018) 主要关注具有瘦矩阵 ()、且具有大量处理和控制单元的单处理周期块结构的矩阵。
非混淆性方法的一个简化处理是将上一阶段的结果与滞后的结果进行回归，并使用估计回归来预测缺失的潜在结果，将其称为水平回归。预测结果如下式：
其中，
二是垂直回归和合成控制文献。合成控制文献 (Abadie 等，2010) 主要关注具有相对较胖 () 或近似正方形矩阵 ()、且具有大量预处理周期的单处理单元块结构的矩阵。
Doudchenko 和 Imbens (2016)、Ferman 和 Pinto (2019) 展示了运用合成控制方法进行回归分析，将其称为垂直回归。预测结果如下式：
其中，
3.4 固定效应和因子模型
面板数据中通常会产生个体固定效应和时间固定效应。水平回归中存在时间固定效应，垂直回归中存在个体固定效应。二者结合会产生双向固定效应。常见的双向固定效应模型如下所示：
更一般的因子模型可以写成：
其中， 是 矩阵， 是 矩阵。
早期文献 (Anderson，1958；Goldberger，1972) 主要研究瘦矩阵，渐进估计是基于时期数目固定而单元数目增加。在现在文献中 (Bai，2003；Pesaran，2006；Moon 和 Weidner，2015；Bai 和 Ng，2017)，研究人员考虑了当 和 都增大时更复杂的渐近性，在归一化处理之后，允许一致地估计因子 和载荷 。
在这类文献中，通常假设因子 大小是固定的，Bai 和 Ng (2002)、Moon 和 Weidner (2015) 讨论了估计秩 的方法。Xu (2017) 将这种交互式固定效应方法应用于具有块状分配特点的矩阵补全问题。Gobillon 和 Magnac (2016)、Kim 和 Oka (2014)、以及 Hsiao 等 (2012) 均对此进行了应用。这种分块结构极大地简化了固定秩估计量的计算。但是，在更复杂的数据缺失模式中，这种方法效率不高，在计算上也不简单。
在机器学习和统计领域关于矩阵补全的文献中 (Srebro 等，2005；Candès 和 Recht，2009；Candès 和 Tao，2010；Gross，2011；Rohde 和 Tsybakov，2011)，研究人员从一个不能完全观察到的矩阵 出发，提出了低秩矩阵模型作为矩阵补全的基础。重点不是估计 和 ，而是 的缺失元素。这些估计量依赖于正则化，特别是核范数正则化。
4. 核范数最小化估计下的矩阵补全
在不考虑协变量的情况下，将数据矩阵 模型化为：
其中 为测量误差。
假设 1： 是独立于 ，且 相互独立且服从亚高斯分布。
此时目标就是估计矩阵 ，注意此时固定效应都包含在 里。为了便于描述估计量，定义两个矩阵：
为了估计 ，文中进一步引入矩阵范数和机器学习中的惩罚回归模型，并引入 来避免对个体和时间的固定效应正则化，以此来构造估计量。
4.1 估计量
估计量的一般形式为：
其中，
上式被称为核范数最小化的矩阵补全估计量 (MC-NNM)。式中 为规则参数 (惩罚因子)。
第一项为损失函数，用 Fröbenius 范数形式表示； 第二项为规则项，用矩阵 的核范数来表示。文中指出了其他形式的范数由于不同的原因在这里并不适用，故选用核范数。
核范数 (Nuclear Norm) 是指矩阵奇异值的和。使用核范数的主要优点是可以使用快速凸优化程序来计算所得到的估计量，它可以用来约束低秩的。
秩可以度量相关性，而矩阵的相关性实际上有带有了矩阵的结构信息。如果矩阵之间各行的相关性很强，那么就表示这个矩阵实际可以投影到更低维的线性子空间，也就是用几个向量就可以完全表达了，它就是低秩的。
低秩矩阵每行或每列都可以用其他的行或列线性表出，可见它包含大量的冗余信息。利用这种冗余信息，可以对缺失数据进行恢复，也可以对数据进行特征提取。因为低秩矩阵是非凸的函数，在优化问题里面很难求解，那么就需要寻找它的凸近似来近似它了，而核范数就是它的凸近似。
4.2 计算估计量
首先不考虑固定效应，定义矩阵收缩算子:
以及矩阵 ：
通过交叉验证选择 的最优值。然后可以计算 ，通过平方项中的一阶条件可以得到 和 。
5. 水平和垂直回归的关系
本节讨论了矩阵补全回归与水平回归 (无混淆性)、垂直回归 (合成控制) 和双重差分方法之间的关系，这是本文的第二个贡献。
将矩阵分块后发现，矩阵补全回归、水平回归、垂直回归、合成控制回归、弹性网估计和双重差分回归估计量之间的关系非常密切，它们都可以被看做是基于完全相同的目标函数，只是在正则化和对目标函数的参数的附加限制方面有所不同。定义目标函数：
矩阵补全回归 (MC-NNM)：
水平回归：
垂直回归：
合成控制：
弹性网络垂直回归：
双重差分回归：
作者进一步对此结果进行评述，指出核范数估计受到的限制较少，相比其他的估计适用范围更广。
6. 两个图示例证
本节目的是比较矩阵补全方法 (MC-NNM) 和以前方法的估计精度。在真实数据矩阵 (矩阵中没有缺失值) 中，选择其子集作为假设处理组，并且在随机选择的初始点之后的时间段，估计它们的值。然后，报告了每种方法在伪处理后生成值的平均均方根误差 (RMSE)，以此来评估各种方法表现是否良好和健壮。具体的，作者通过两个实际事例的数据比较了以下 5 种方法：
双重差分回归 (DID) 弹性网络垂直回归 (VT-EN) 弹性网络水平回归 (HR-EN) 原始合成控制方法 (SC-ADH) 核范数最小化矩阵补全法 (MC-NNM)
第一个例子是使用 Abadie 等 (2010) 研究的加州吸烟数据作为控制组 。然后，人为地指定一些要处理的州和时间段，并将其预测值与实际值进行比较。
图 (a) 为随机同时处理方式，图 (b) 为随机交错处理方式。综合两种方式可以看出，核范数最小化矩阵补全法 (MC-NNM) 是相对更好的方法。它的性能改进归功于它使用了额外的观察值 (即处理单元的未经处理时的值)。
第二个例子是使用股票市场 2453 只股票近 10 年的日收益率。为了观察 RMSE 的统计波动，通过 个随机抽样股票的前 个日收益率来创建 50 个子样本。矩阵形状从很瘦到很胖，，，，，。在这里，重点比较了矩阵的形状改变时，弹性网络水平回归 (HR-EN)、弹性网络垂直回归 (VT-EN) 和核范数最小化矩阵补全法 (MC-NNM) 三种估计方法。
图中显示了平均 RMSE。通过下图第一幅图发现，核范数最小化矩阵补全法 (MC-NNM) 无论矩阵胖瘦都可以进行相对较好的估计。下图中第二幅图是用核范数最小化矩阵补全法 (MC-NNM) 来观察矩阵 秩的变化。可以看出，当 相对于 增长时，矩阵的秩变小。这表明，不同股票的日收益率之间的相关性比同一股票不同时间段的日收益率之间的相关性更强。
7. 一般化扩展
这节主要简要讨论了文中的估计量利用方式的扩展。文中提到了 7 种情况，这里简要介绍前 4 种。
一是模型中加入协变量：前文描述的基本模型中并没有考虑协变量的影响。如果考虑到个体内、时间内以及个体与时间之间的相互影响，作者对原有模型修正为：
上式中 表示个体内的协变量， 表示时间内的协变量， 表示个体与时间之间的协变量，、 分别反映其协变量效应，、 分别反映个体固定效应和时间固定效应。此时可将估计量扩展为下式，来估计 ：
二是利用处理组的数据：(1) 当处理效应恒定为一常数时，令 表示该常数，则可将估计量扩展为下式：
(2) 当处理效应是低秩形式时，即它不是常数，而是矩阵 预期是低秩的。此时可以利用估算 的方式去估算 。然后利用 和 来估计处理效应 。
三是自相关误差：原估计量的一个缺点是没有考虑观测的时间序列性质。误差项 在时间上很可能是相关的，使用一阶自回归模型，在不考虑协变量的情况下，目标函数可扩展为：
其中， 可以反映 的自相关程度。
四是加权损失函数：原估计量另一个限制是，对差异结果赋予同样的权重来评估平均因果效应。但实际上各个差异结果的地位并不完全相同，此时就需要对损失函数进行合理加权。在这里，作者采用倾向性得分加权法对原估计量进行了扩展 (矩阵 是表示倾向性得分 的矩阵)。
8. 参考资料
陈强. 高级计量经济学及 Stata 应用. 高等教育出版社, 2014. 第二十八章：处理效应. 陈强. 机器学习及 R 应用. 高等教育出版社, 2020. 第九章：惩罚回归. Hansen, Bruce E..Ecomometrics. 2021. Chapter 29: Machine Learning. 机器学习中的范数规则化之 (一)：L0、L1 与 L2 范数 -Link- 机器学习中的范数规则化之 (二)：核范数与规则项参数选择 -Link- 因果推断学习笔记 -Link-
9. 相关推文
Note：产生如下推文列表的 Stata 命令为：
lianxh 因果推断 合成控制, m
安装最新版lianxh命令：
ssc install lianxh, replace
专题：计量专题 因果推断：哪本教材适合我？ 因果推断新书在线读：Causal Inference-The Mixtape 专题：合成控制法 Stata：纠偏合成控制法介绍-allsynth Synth_Runner命令：合成控制法高效实现 Stata：合成控制法程序分享 Stata：合成控制法-synth-命令无法加载-plugin-的解决办法 合成控制法 (Synthetic Control Method) 及 Stata实现 专题：内生性-因果推断 Stata：因果推断方法综述和Stata操作 fect：基于面板数据的因果推断（上）-T218a fect：基于面板数据的因果推断（下）-T218b 因果推断：未测量混杂因素的敏感性分析-T249 用FE-固定效应模型能做因果推断吗？ 经典文献回顾：政策评价-因果推断的计量方法 因果推断好书：Causal-Inference-Measuring-the-Effect-of-X-on-y
课程推荐：因果推断实用计量方法
主讲老师：丘嘉平教授
🍓 课程主页：https://gitee.com/lianxh/YGqjp
New！ Stata 搜索神器：
lianxh和songblGIF 动图介绍
搜： 推文、数据分享、期刊论文、重现代码 ……
👉 安装：
. ssc install lianxh
. ssc install songbl
👉 使用：
. lianxh DID 倍分法
. songbl all
🍏 关于我们
连享会 ( www.lianxh.cn，推文列表) 由中山大学连玉君老师团队创办，定期分享实证分析经验。 直通车： 👉【**百度一下：**连享会】即可直达连享会主页。亦可进一步添加 「知乎」,「b 站」,「面板数据」,「公开课」 等关键词细化搜索。