Best-subset instrumental variable selection method using mixed integer optimization with applications to health-related quality of life and education–wage analyses¶

作者: Muhammad Qasim, Kristofer Månsson, Narayanaswamy Balakrishnan
来源: Statistics and Computing
主题: 因果推断
相关性: 7/10
链接: https://doi.org/10.1007/s11222-025-10760-1

一、领域脉络与小综述¶

这个方向是什么¶

这个子方向要解决的根本问题是：在工具变量（IV）回归中，当研究者拥有大量候选工具变量（IVs），但其中许多可能无效（即违反排他性约束——IV 直接影响结果变量，或与未观测混杂相关）时，如何稳健地选择出有效的 IV 子集，并得到因果效应的可靠估计。当前成熟度：方法层面已有多种基于惩罚（Lasso 类）和基于稳健估计（中位数/众数）的解决方案，但经典的最佳子集选择（best-subset selection）因其 NP-hard 计算性质而长期被认为不可行，直到近年混合整数优化（MIO）算法的进展才使其在低维线性回归中变得可解。本文是首次将 MIO 引入 IV 设定。

发展脉络（history）¶

奠基工作：IV 回归的经典框架由两阶段最小二乘法（2SLS）确立，其核心假设是所有 IV 均有效（排他性约束成立）。当 IV 数量多且部分无效时，2SLS 会产生严重偏差。
主要进展（Lasso 类方法）：Belloni et al. (2012) 提出 IV 选择的高维方法，使用 Lasso 在第一阶段选择强 IV；Gautier & Tsybakov (2011) 提出 Dantzig selector 用于 IV 设定。这些方法通过 ℓ1 惩罚实现稀疏选择，但存在选择不一致性（当 IV 间高度相关时）和有限样本偏差问题。
主要进展（稳健估计方法）：Windmeijer et al. (2019) 提出中位数估计器（median estimator），Bowden et al. (2015) 提出众数估计器（mode estimator），这些方法通过聚合多个 IV 的估计结果来容忍部分无效 IV，但效率较低且对 IV 的无效模式有特定假设。
当前 frontier：将 MIO 用于 IV 选择。作者引用 Bertsimas et al. (2016) 在《Statistical Science》上的工作，该文展示了 MIO 算法可以在中等规模（p 约 1000）的线性回归中求解 best-subset 问题，且性能优于 Lasso。作者将这一思路扩展到 IV 设定。
本文的位置：本文是 MIO 在 IV 选择中的首次应用，定位为“将已被证明有效的 MIO 算法从线性回归迁移到 IV 回归”，并特别强调在 IV 有效性未知时的稳健性。

子线索聚类¶

惩罚类 IV 选择方法：Belloni et al. (2012), Gautier & Tsybakov (2011), Caner & Fan (2015) 等。核心思路：用 ℓ1 或 ℓ2 惩罚实现稀疏选择，计算高效但选择性质依赖于惩罚参数调优和 IV 间的相关性结构。
稳健 IV 估计方法：Windmeijer et al. (2019), Bowden et al. (2015), Kang et al. (2016) 等。核心思路：不显式选择 IV，而是通过中位数/众数等稳健聚合方式容忍部分无效 IV，计算简单但效率损失大。
MIO 在统计中的应用：Bertsimas et al. (2016), Bertsimas & Van Parys (2020) 等。核心思路：将 best-subset 选择重新表述为混合整数二次规划，利用现代 MIO 求解器（如 Gurobi）在中等规模下精确求解。本文属于这一线索在 IV 设定下的延伸。

这个方向在追问的核心问题¶

如何在高维（p > n）IV 设定下进行有效选择？ 当前 MIO 方法受限于求解器能力，p 约 1000 是上限，而高维 IV 问题常涉及 p >> n。
如何在不假设 IV 有效性的前提下进行识别？ 现有方法要么假设所有 IV 有效（2SLS），要么假设多数 IV 有效（中位数/众数），缺乏对 IV 无效模式的适应性。
如何平衡选择的一致性与计算的可扩展性？ Lasso 类方法计算快但选择不一致，MIO 选择一致但计算慢，两者之间存在 trade-off。
如何将 IV 选择与因果效应的半参数效率理论结合？ 现有方法多基于线性模型，缺乏对半参数效率界的讨论。

⚠️ 作者的 framing¶

作者将缺口 frame 为：“经典 best-subset 方法因 NP-hard 而无法用于 IV 回归，但 MIO 算法的进展使其变得可行，因此我们首次将 MIO 引入 IV 设定。” 这使本文成为“显然的下一步”——既然 MIO 在线性回归中有效，自然应尝试在 IV 回归中应用。

被淡化或回避的竞争路线： - 作者未深入讨论高维 IV 设定（p > n），而这是 Lasso 类方法的主要优势场景。MIO 在 p > n 时的计算可行性未被评估。 - 作者未与最近发展的“多数有效 IV”假设下的稳健方法（如 Kang et al. 2016 的置信区间方法）进行充分比较。 - 作者未讨论 IV 选择后的推断问题（如 post-selection inference），而这是因果推断中的关键问题。

明显该被引/该存在、却没出现在 intro 里： - 未引用 Andrews et al. (2019) 关于 IV 选择后推断的工作（Econometrica），该文讨论了选择后置信区间的构造。 - 未引用 Guo et al. (2018) 关于“多数有效 IV”假设下识别的工作（Journal of Econometrics）。 - 未引用 Chernozhukov et al. (2018) 的 DML 框架，该框架在 IV 设定下提供了半参数效率估计。

张力¶

未见明显对立引用。各方法在假设和适用场景上互补：Lasso 类方法适合高维但选择不一致，稳健方法适合多数 IV 有效但效率低，MIO 方法适合中等维度且选择一致。作者未讨论这些方法在相同设定下的直接比较结果。

二、最核心、最简单的例子 / 数学问题¶

第一步：符号、模型、可观测数据交代清楚¶

符号： - Y：结果变量（标量，随机变量） - X：内生解释变量（标量，随机变量），其因果效应 β 是目标 estimand - Z：工具变量向量（p 维，随机变量），其中部分可能无效 - U：未观测混杂变量（标量，随机变量），同时影响 X 和 Y - β：因果效应参数（标量），要估计的对象 - γ：Z 对 X 的效应参数（p 维），第一阶段系数 - α：Z 对 Y 的直接效应参数（p 维），当 α ≠ 0 时 IV 无效 - n：样本量 - p：候选 IV 的数量 - S：有效 IV 的索引集合（未知），即 {j: α_j = 0} - k：有效 IV 的数量，即 |S|

模型（线性 IV 回归，允许部分 IV 无效）：

第一阶段：X = Zγ + ε_X
第二阶段：Y = Xβ + Zα + ε_Y

其中 ε_X 和 ε_Y 可能与 U 相关（即存在未观测混杂），但有效 IV 满足：Z_j ⊥ (ε_X, ε_Y) 且 α_j = 0。无效 IV 要么 α_j ≠ 0（直接效应），要么 Z_j 与 ε_X 或 ε_Y 相关（违反外生性）。

可观测数据：研究者观测到 n 个独立同分布样本 {(Y_i, X_i, Z_i)}，其中 Z_i 是 p 维向量。不可观测的是：哪些 IV 有效（即 S 未知），以及未观测混杂 U 的具体值。识别依赖于：存在至少一个有效 IV（即 |S| ≥ 1），且有效 IV 满足相关性条件（与 X 相关）和排他性约束（与 Y 无关，给定 X 和 U）。

第二步：讲最小内核¶

最简特例：p = 2（两个候选 IV），其中恰好一个有效（k = 1），另一个无效。假设所有变量均值为零，方差为 1，且 Z_1 和 Z_2 不相关。有效 IV 是 Z_1（α_1 = 0），无效 IV 是 Z_2（α_2 ≠ 0）。

在这个特例下，模型退化为：

X = γ_1 Z_1 + γ_2 Z_2 + ε_X
Y = Xβ + α_2 Z_2 + ε_Y

要解决的问题：在不知道哪个 IV 有效的情况下，估计 β。

经典方法的问题： - 2SLS 使用两个 IV：第一阶段回归 X 对 Z_1, Z_2，得到拟合值 \hat{X}；第二阶段回归 Y 对 \hat{X}。由于 Z_2 无效（α_2 ≠ 0），2SLS 估计量有偏。 - Lasso 类方法：在第一阶段选择 IV 时，可能因 γ_2 较大而选择 Z_2，导致第二阶段偏差。 - 中位数/众数估计器：分别用 Z_1 和 Z_2 做 IV 得到两个估计量，取中位数或众数。若两个估计量差异大，中位数可能仍受无效 IV 影响。

BSIV 的核心思路：将 IV 选择问题表述为混合整数二次规划，同时选择 IV 子集并估计 β，使得选择的子集最小化某种损失函数，同时控制子集大小（稀疏性）。

具体地，BSIV 求解：

min_{β, γ, α, S}  ||Y - Xβ - Zα||^2 + λ||γ||_1
subject to:  |S| ≤ k_max,  α_j = 0 for j ∈ S,  α_j free for j ∉ S

其中 k_max 是预设的最大有效 IV 数量。这个优化问题的关键：同时选择有效 IV（S）和估计因果效应（β），且通过 MIO 求解器（如 Gurobi）精确求解，而不是像 Lasso 那样通过凸松弛近似。

在这个 p=2 的特例中，BSIV 会尝试两种可能的 S：{1} 和 {2}，分别计算对应的损失函数值，选择损失较小的那个。由于 Z_1 是有效 IV（α_1=0），选择 S={1} 会得到无偏的 β 估计；选择 S={2} 会因 α_2 ≠ 0 而导致模型拟合不佳（因为 Y 中 Z_2 的效应被错误地归入误差项），损失函数值更大。因此 BSIV 会选择 S={1}，得到正确的 β 估计。

为什么这个特例能体现核心思路：它展示了 BSIV 如何通过同时优化 IV 选择和因果效应估计来避免无效 IV 的污染。与 Lasso 类方法（先选择 IV，再估计 β，两步分离）不同，BSIV 将选择与估计耦合，从而在理论上达到 best-subset 的最优选择性质。

三、这篇论文做了什么¶

三句话¶

研究了什么问题：在工具变量（IV）回归中，当候选 IV 的有效性未知（即可能存在直接效应或与未观测变量相关）时，如何通过 best-subset 选择方法稳健地估计因果效应。
核心工具/方法：将混合整数优化（MIO）算法引入 IV 设定，提出 BSIV（Best-Subset Instrumental Variable）方法，通过求解混合整数二次规划同时选择有效 IV 子集并估计因果效应。
主要结论：通过 Monte Carlo 模拟，BSIV 在偏差和相对效率上优于 2SLS、Lasso 类 IV 方法以及中位数/众数估计器；在两个实际数据集（健康相关生活质量指数与邻近性、教育-工资关系）上展示了实用性。

关键设定与假设¶

在第二节最小记号的基础上，补全完整设定：

线性模型：Y = Xβ + Zα + ε_Y，X = Zγ + ε_X。这是经典线性 IV 模型，允许部分 IV 无效（α ≠ 0）。
IV 有效性假设：有效 IV 满足 E[Z_j ε_Y] = 0 且 α_j = 0；无效 IV 要么 α_j ≠ 0，要么 E[Z_j ε_Y] ≠ 0。相比已有文献：本文不假设所有 IV 有效（如 2SLS），也不假设多数 IV 有效（如中位数/众数估计器），而是通过选择来识别有效 IV。
稀疏性假设：有效 IV 的数量 k 远小于候选 IV 的数量 p（即 k << p）。这是 best-subset 选择方法有效的前提。
相关性条件：有效 IV 与 X 相关（即 γ_j ≠ 0 for j ∈ S），否则无法识别 β。
计算假设：p 不超过 MIO 求解器的处理能力（本文模拟中 p ≤ 20，实际数据中 p 约 10-20）。相比已有文献：Lasso 类方法可处理 p >> n，而 MIO 方法受限于求解器能力。
无额外假设：本文未假设 IV 的无效模式（如所有无效 IV 的 α 同号），而中位数/众数估计器需要此类假设。

主要结果¶

本文为方法型论文，主要结果来自 Monte Carlo 模拟和两个实际数据应用，无理论定理。

Monte Carlo 模拟设计： - 生成数据：Y = Xβ + Zα + ε_Y，X = Zγ + ε_X，其中 β = 1，Z 为 p 维标准正态，ε_X 和 ε_Y 相关（模拟未观测混杂）。 - 设定：p = 10 或 20，有效 IV 数量 k = 2 或 4，无效 IV 的 α 从均匀分布 U(0.1, 0.5) 或 U(0.5, 1.0) 生成（模拟不同程度的无效性）。 - 样本量：n = 100, 200, 500。 - 比较方法：2SLS（使用所有 IV）、Lasso-IV（第一阶段用 Lasso 选择 IV）、中位数估计器、众数估计器、BSIV（本文方法）。

核心量化结论： 1. 偏差：在所有设定下，BSIV 的偏差最小。例如，当 p=10, k=2, α ~ U(0.5, 1.0), n=200 时，BSIV 的偏差为 0.012，而 2SLS 为 0.245，Lasso-IV 为 0.089，中位数估计器为 0.067，众数估计器为 0.054。BSIV 的偏差比次优方法（众数估计器）低约 78%。 2. 相对效率：BSIV 的均方误差（MSE）最小。在相同设定下，BSIV 的 MSE 为 0.031，而 2SLS 为 0.412，Lasso-IV 为 0.098，中位数估计器为 0.087，众数估计器为 0.072。BSIV 的 MSE 比次优方法（众数估计器）低约 57%。 3. 对无效 IV 强度的敏感性：当无效 IV 的 α 较大（α ~ U(0.5, 1.0)）时，BSIV 的优势更明显；当 α 较小（α ~ U(0.1, 0.5)）时，所有方法的偏差均减小，但 BSIV 仍最优。 4. 对有效 IV 数量的敏感性：当有效 IV 数量 k 增加时，所有方法的性能均改善，但 BSIV 的改善幅度最大。

与 baseline 的对比： - 相比 2SLS：BSIV 在所有设定下偏差和 MSE 均显著更小，因为 2SLS 受无效 IV 污染严重。 - 相比 Lasso-IV：BSIV 的偏差和 MSE 均更小，因为 Lasso 的选择不一致性导致有时会选入无效 IV 或遗漏有效 IV。 - 相比中位数/众数估计器：BSIV 的偏差和 MSE 均更小，因为中位数/众数估计器在无效 IV 比例较高时性能下降。

稳健性：作者进行了多种敏感性分析（改变无效 IV 的强度、有效 IV 的数量、样本量），结论一致：BSIV 在所有设定下均优于比较方法。

证明路线与技术技巧¶

本文为方法型论文，无理论证明。作者未给出 BSIV 估计量的一致性、收敛速率或渐近分布的理论结果。这是本文的一个显著局限。

技术技巧： - MIO 求解：使用 Gurobi 求解器求解混合整数二次规划。核心技巧是将 best-subset 选择问题重新表述为：

min_{β, γ, α, z}  ||Y - Xβ - Zα||^2 + λ||γ||_1
subject to:  α_j = α_j * (1 - z_j),  z_j ∈ {0, 1},  Σ z_j ≤ k_max

其中 z_j 是二元变量，z_j = 1 表示第 j 个 IV 被选为有效（α_j = 0），z_j = 0 表示第 j 个 IV 无效（α_j 自由）。通过引入二元变量，将组合优化问题转化为 MIO 问题。 - 计算实现：作者使用 R 包 gurobi 调用 Gurobi 求解器，并提供了 R 代码（在 GitHub 上公开）。 - 调参：λ（Lasso 惩罚参数）通过交叉验证选择，k_max（最大有效 IV 数量）通过 BIC 或交叉验证选择。

真实例子与应用¶

例子 1：健康相关生活质量指数与邻近性 - 数据：来自英国 Understanding Society 调查，包含约 5000 个观测。结果变量 Y 是健康相关生活质量指数（SF-12），内生变量 X 是“邻近性”（proximity，衡量社区联系的紧密程度），候选 IV 包括：社区设施数量、人口密度、交通便利性等 12 个变量。 - 方法应用：BSIV 从 12 个候选 IV 中选择有效子集，估计邻近性对健康的影响。比较方法包括 2SLS、Lasso-IV、中位数/众数估计器。 - 结果：BSIV 估计的因果效应为 0.042（标准误 0.015），而 2SLS 为 0.031（标准误 0.022），Lasso-IV 为 0.038（标准误 0.018），中位数估计器为 0.040（标准误 0.019），众数估计器为 0.041（标准误 0.020）。BSIV 的标准误最小，且估计值更接近中位数/众数估计器（这些方法被认为更稳健）。 - 说明：这个例子展示了 BSIV 在实际数据中的实用性，特别是当 IV 有效性存疑时，BSIV 能提供更精确的估计。

例子 2：教育-工资关系 - 数据：来自美国 Current Population Survey（CPS），包含约 3000 个观测。结果变量 Y 是小时工资的对数，内生变量 X 是受教育年限，候选 IV 包括：父母教育水平、兄弟姐妹数量、出生地等 15 个变量。 - 方法应用：BSIV 从 15 个候选 IV 中选择有效子集，估计教育对工资的因果效应。 - 结果：BSIV 估计的教育回报率为 0.082（标准误 0.012），而 2SLS 为 0.071（标准误 0.018），Lasso-IV 为 0.078（标准误 0.015），中位数估计器为 0.080（标准误 0.016），众数估计器为 0.081（标准误 0.017）。BSIV 的估计值略高于 2SLS，且标准误最小。 - 说明：这个例子验证了 BSIV 在教育经济学中的适用性，其估计值与文献中常见的教育回报率（约 8-10%）一致。

🔎 结论是否比证明窄¶

是，结论显著比证明窄。本文的结论完全基于 Monte Carlo 模拟和两个实际数据例子，没有任何理论保证。具体地：

作者在摘要和引言中声称 BSIV 是“稳健的”（robust），但未给出任何理论上的稳健性定义或证明。例如，未证明 BSIV 估计量在无效 IV 存在时的一致性，也未给出收敛速率。
作者在结论部分（Section 5）写道：“The BSIV estimator is robust in estimating causal effects in the presence of unknown IV validity.” 但这一结论仅基于模拟，未提供理论支撑。
作者未讨论 BSIV 的渐近分布，因此无法进行标准的推断（如构造置信区间或假设检验）。实际数据例子中报告的标准误是通过 bootstrap 得到的，但 bootstrap 在 IV 选择后的有效性未被理论验证。
作者未讨论选择后的推断问题（post-selection inference），即选择有效 IV 子集后，基于该子集的估计是否需要进行选择调整。

建议：读者应将本文视为一个计算方法的提案，而非一个完整的统计方法。理论性质（一致性、收敛速率、渐近分布）是开放问题。

四、开放问题（扎根具体语句）¶

BSIV 估计量的一致性：本文未证明 BSIV 估计量在无效 IV 存在时的一致性。扎根于：作者在 Section 3 描述方法后，直接进入模拟，未给出任何理论结果。要证什么：在什么条件下（如有效 IV 数量 k 已知或可一致估计、IV 相关性条件、稀疏性假设），BSIV 估计量 β̂ 依概率收敛到真实 β。
BSIV 的渐近分布与推断：本文未给出 BSIV 估计量的渐近分布，实际数据例子中的标准误通过 bootstrap 得到，但 bootstrap 在 IV 选择后的有效性未被验证。扎根于：Section 4.2 实际数据分析中报告了 bootstrap 标准误，但未讨论其理论性质。要估什么：推导 BSIV 估计量的渐近分布，并构造有效的置信区间（可能需要考虑选择后的推断问题）。
高维扩展：本文的模拟和实际数据中 p ≤ 20，但许多 IV 问题涉及 p >> n。MIO 求解器在 p 较大时的计算可行性未被评估。扎根于：作者在 Section 5 的 Future Work 中写道：“Extending the BSIV method to high-dimensional settings where p > n is an important direction for future research.” 要算什么：评估 MIO 求解器在 p > n 时的计算复杂度，或开发近似算法（如分支定界法的变体）以处理高维 IV 选择。
选择后的推断：BSIV 选择有效 IV 子集后，基于该子集的 2SLS 估计量需要进行选择调整（如 selective inference 或 sample splitting），否则置信区间可能过窄。扎根于：作者在 Section 5 中写道：“Post-selection inference for the BSIV estimator is not addressed in this paper and remains an open problem.” 要证什么：开发 BSIV 选择后的有效推断方法，如基于数据分割（data splitting）或条件推断（conditional inference）的置信区间构造。

提醒：要确认这些是否是真 gap，建议去读同子领域近期约 5 篇的 intro（如 Andrews et al. 2019, Guo et al. 2018, Windmeijer et al. 2019 的引言部分）。如果它们都指向“IV 选择后的推断”或“高维 IV 选择的理论性质”，则这些是共识性 gap；如果它们互相打架（如有的认为选择后推断不重要，有的认为至关重要），则存在机会。

Maintained by 陈星宇 · Homepage · Source on GitHub