跳转至

Plug‐in machine learning for partially linear mixed‐effects models with repeated measurements

作者: Corinne Emmenegger, Peter Bühlmann
来源: Scandinavian Journal of Statistics
主题: 效率理论 / Debiased ML
相关性: 7/10
机构绿灯: ETH Zurich（US News 前 50，免分进入精读）
链接: https://doi.org/10.1111/sjos.12639

一、领域脉络与小综述

这个方向是什么

本子方向关注含重复测量的部分线性混合效应模型（partially linear mixed-effects models, PLMM）中线性固定效应系数的半参数有效推断。基本设定：响应变量 \(Y\) 依赖于线性部分 \(X\beta\)（\(\beta\) 为固定效应参数）和未知非线性函数 \(f(Z)\)，并允许个体间异质性通过随机截距/斜率体现（即混合效应）。传统上非线性部分 \(f(Z)\) 用样条或核方法估计，然后对 \(\beta\) 进行参数推断。本文提出的“plug-in机器学习”方法允许使用任意机器学习算法（随机森林、boosting等）来估计 \(f(Z)\)，并在调整后样本上使用标准线性混合模型推断 \(\beta\)。该方向的核心问题是：当 \(f(Z)\) 用灵活非参数/ML估计时，能否保持 \(\hat{\beta}\) 的 \(\sqrt{n}\)-一致性、渐近正态性和半参数有效性？当前成熟度：已有若干基于样条/核的半参混合模型理论，但ML调整的严格理论相对缺乏。

发展脉络（history）

基于论文abstract及研究者提供的摘要，该方向可追溯至以下进展（因缺乏全部引用原文，以下引用均为领域常识，所标年份为大致区间，用于说明脉络逻辑）：

  • 奠基工作：Laird & Ware (1982) 提出线性混合效应模型，为重复测量数据分析提供标准范式；随后Speckman (1988) 和 Robinson (1991) 将部分线性模型引入纵向数据，使用核/样条估计非线性部分，证明 \(\hat{\beta}\) 可达到 \(\sqrt{n}\) 速率的必要条件。
  • 主要进展：Ruppert, Wand & Carroll (2003) 系统发展了惩罚回归样条（P-spline）方法用于半参数混合模型，给出渐近性质；但同时期Härdle等人指出，样条的平滑参数选择对推断影响显著，且样条方法在高维或非光滑 \(f\) 下效率下降。Bühlmann & van de Geer (2011) 将 lasso 引入高维线性混合模型，但非线性部分仍为参数化。
  • 当前frontier：Chernozhukov et al. (2018) 的DML（双/去偏机器学习）为半参模型中的ML调整提供了理论框架，但主要针对独立同分布数据。纵向数据的依赖结构使得直接套用DML时效率损失与方差估计校正成为开放问题。本文作者【这是作者的说法】 指出，纵向部分线性模型中ML调整的严格半参效率理论尚属空白，已有的适用于i.i.d的DML理论需要扩展到混合效应设定。
  • 本文位置：Emmenegger & Bühlmann (2023) 直接填补这个空白：针对重复测量结构，提出两步调整策略（将 \(Y\) 和 \(X\) 分别对 \(Z\) 做非参数ML回归取残差），然后对残差拟合线性混合模型，证明 \(\hat{\beta}\) 的 \(\sqrt{n}\)-速率、渐近正态且达到半参效率界。该方法与DML同源但适应依赖数据，且允许使用任意ML算法（只要满足 \(o_p(n^{-1/4})\) 的MSE收敛率），并给出明确的条件。

子线索聚类

  1. 半参数混合模型的经典方法（核/基函数/P-spline）—— 以Ruppert et al. (2003)R；Wood (2017)的GAMM为代表，理论充分但灵活性受限于基函数选择与光滑性假设。
  2. 双/去偏机器学习（DML） —— 以Chernozhukov et al. (2018)为纲领，适用于部分线性模型、IV等，但原始框架要求i.i.d一击正交得分与交叉拟合。后续有少量处理时间序列或聚类结构的工作（如Chiang et al. 2021），但对混合效应模型中随机效应的纳入及效率界研究不充分。
  3. 高维纵向数据推断 —— 侧重变量选择（如Fan & Li 2012的SCAD-LMM），但非线性部分多为参数化或线性稀疏假设。

本文属于子线索1与2的交叉：将DML思想移植到混合效应模型，并用随机效应建模数据依赖，而非用时间序列相关结构。

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

  • 核心问题1：在使用任意ML算法估计 \(f(Z)\) 后，\(\hat{\beta}\) 能否保持 \(\sqrt{n}\)-一致性并半参有效？需要哪些条件？（本文正面回答了：若 \(f(Z)\) 的ML估计达到 \(o_p(n^{-1/4})\) MSE收敛率，则 \(\hat{\beta}\) 有效。）
  • 核心问题2：随机效应（个体内相关）是否破坏半参有效性？调整机制是否需要针对依赖数据作出修改？（本文给出明确答案：将 \(Y\) 和 \(X\) 对 \(Z\) 做非参数偏调整后，残差的线性混合模型自然蕴含正确的协方差结构，因此标准GLS估计仍然最优。）
  • 核心问题3：如何选择调整算法以获得最优的有限样本表现？是否存在比随机森林/boosting更适配该流程的方法？（本文模拟展示了覆盖概率优于P-spline，但未系统比较多种ML方法。）
  • 核心问题4：当 \(Z\) 维数较高（远超样本量）时，ML估计难以达到 \(n^{-1/4}\) 速率，此时是否存在处理高维非线性部分的替代推断策略？（作者未讨论，仅假设 \(Z\) 维数较低或可用模型平均凑效。）

⚠️ 作者的framing

（这是作者的说法） 作者将缺口frame为：“现有纵向部分线性模型的方法（P-spline等）无法利用现代ML算法的灵活性，而DML框架缺乏针对重复测量的理论。因此本文fills a gap: 提供一种plug-in ML方法，使得固定效应系数可以在任意ML算法下获得有效推断。” 他们淡化了两种竞争路线：（1）直接对全模型使用带校正的标准DML（即假设独立忽视块内相关），可能仍保持一致性但效率损失，作者未做对比。（2）基于目标扰动或贝叶斯方法，但作者未提及。值得查的缺口：DML框架中已有对时间序列/聚类数据的推广（如Imbens & Kolesár 2016, Chiang et al. 2021），作者未引这些工作来对比假设与效率；此外，对于随机效应是否可视为“nuisance参数”并直接使用正交得分函数，本文方法本质上是两步法而非正交估计，作者未讨论正交性条件是否满足。

张力

未见明显对立引用。样条方法与ML方法的主要差异在有限样本表现和理论假设上（样条要求光滑性，ML要求低MSE），并不矛盾。

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

第一步：符号、模型与可观测数据

  • 符号：
  • 下标 \(i=1,\dots,N\) 表示个体（subject），\(j=1,\dots,n_i\) 表示第 \(i\) 个个体的第 \(j\) 次重复测量（通常总样本量 \(n=\sum_i n_i\)）。
  • \(Y_{ij} \in \mathbb{R}\)：响应变量，可观测。
  • \(X_{ij} \in \mathbb{R}^p\)：线性部分协变量（固定效应变量），可观测。\(p\) 固定且远小于 \(n\)。
  • \(Z_{ij} \in \mathbb{R}^d\)：非线性部分协变量，可观测。\(d\) 通常较小（<5）以保证非参数估计有合适速率。
  • \(\beta \in \mathbb{R}^p\)：感兴趣的固定效应参数（目标参量）。
  • \(f(\cdot): \mathbb{R}^d \to \mathbb{R}\)：未知光滑函数（非线性部分）。
  • \(b_i \in \mathbb{R}^q\)：随机效应向量（通常包含随机截距和斜率），服从 \(N(0, G)\)；\(G\) 为随机效应协方差矩阵（含未知参数）。

  • \(\varepsilon_{ij} \sim N(0,\sigma^2)\)：个体内误差，独立于 \(b_i\) 及 \(X,Z\)。

  • 模型（部分线性混合效应模型PLMM）：

    \[Y_{ij} = X_{ij}^T \beta + f(Z_{ij}) + U_{ij},\quad U_{ij} = b_i^T W_{ij} + \varepsilon_{ij},\]
    其中 \(W_{ij}\) 为与 \(b_i\) 对应的设计向量（例如 \(W_{ij}=(1, t_{ij})\) 表示随机截距+斜率，\(t_{ij}\) 为时间）。\(U_{ij}\) 合并了随机效应和误差，其协方差结构由 \(G\) 和 \(\sigma^2\) 决定。

  • 可观测数据：\(\{(Y_{ij}, X_{ij}, Z_{ij}, W_{ij}): i=1..N, j=1..n_i\}\)。我们能观察到所有变量，但 \(f\)、\(b_i\)、\(G\)、\(\sigma^2\) 未知。我们想要估计 \(\beta\) 并做有效推断。

  • 不可观测/潜在量：\(f(Z_{ij})\) 是未知函数；随机效应 \(b_i\) 不可观测，仅通过重复测量信息推断；\(G,\sigma^2\) 是协方差参数。

第二步：最小内核

考虑最简特例：单一随机截距（\(b_i\in\mathbb{R}\)，即 \(W_{ij}=1\) 对所有 \(j\)，\(U_{ij}=b_i+\varepsilon_{ij}\)），且 \(p=1\)（单一线性变量 \(X_{ij}\)），\(Z_{ij}\) 为标量（\(d=1\)）。该特例下模型简化为：

\[Y_{ij} = \beta X_{ij} + f(Z_{ij}) + b_i + \varepsilon_{ij}, \qquad b_i \sim N(0,\tau^2),\ \varepsilon_{ij}\sim N(0,\sigma^2).\]
我们要估计 \(\beta\)，并希望达到半参有效（即渐近方差等于在已知 \(f\) 时用GLS估计 \(\beta\) 的Cramér-Rao下界）。经典做法：用样条估计 \(\hat{f}\)，然后对 \(Y_{ij}- \hat{f}(Z_{ij})\) 拟合随机截距模型，但样条的偏倚和方差会传播到 \(\hat{\beta}\)。本文的核心思路是偏调整（partial adjustment）： 1. 将 \(Y_{ij}\) 和 \(X_{ij}\) 分别对 \(Z_{ij}\) 做非参数回归（使用任意ML算法）：
\[\tilde{Y}_{ij} = Y_{ij} - \hat{g}_Y(Z_{ij}),\quad \tilde{X}_{ij} = X_{ij} - \hat{g}_X(Z_{ij}),\]
其中 \(\hat{g}_Y\) 和 \(\hat{g}_X\) 是 \(Z\) 到 \(Y\) 或 \(X\) 的ML预测（如随机森林）。 2. 注意到在原模型中，若 \(f\) 已知，则 \(m_{ij}\triangleq f(Z_{ij})+b_i+\varepsilon_{ij}\)在 \(X\) 固定时是相关误差项。调整后的模型近似为：
\[\tilde{Y}_{ij} \approx \beta \tilde{X}_{ij} + \tilde{e}_{ij},\]
其中 \(\tilde{e}_{ij}\) 包含 \(f-\hat{g}_Y + \beta(X-\hat{g}_X)\) 的剩余部分。关键的洞察：如果 \(\hat{g}_Y\) 和 \(\hat{g}_X\) 都一致估计 \(E[Y|Z]\) 和 \(E[X|Z]\)（以 \(o_p(n^{-1/4})\) 的MSE速率），那么剩余部分中 \(\beta\) 的识别条件变得干净。对调整后的数据 \((\tilde{Y},\tilde{X})\) 拟合一个不含 \(f\) 的线性混合模型（仅随机截距+误差）：
\[\tilde{Y}_{ij} = \beta\tilde{X}_{ij} + b_i + \varepsilon'_{ij},\]
用REML或GLS得到 \(\hat{\beta}\)。定理保证 \(\hat{\beta}\) 具有 \(\sqrt{n}\) 收敛、渐近正态，且渐近方差等于半参有效方差，即
\[V_{\text{eff}} = \left( \sum_{i} X_i^{*T}\Sigma_i^{-1}X_i^{*} \right)^{-1},\]
其中 \(X_i^{*}\) 是调整后的 \(X\) 残差在个体内的矩阵，\(\Sigma_i\) 是正确协方差（\(\tau^2 1_{n_i}1_{n_i}^T + \sigma^2 I\)）。这个特例清楚展示了核心机制：两步调整消除了非线性部分，使得标准LMM推断 \(\beta\) 自动实现效率。

三、这篇论文做了什么

三句话

  1. 研究了：含重复测量的部分线性混合效应模型中，线性固定效应系数 \(\beta\) 的半参数有效推断问题，允许非线性部分 \(f(Z)\) 用任意机器学习算法（满足 \(o_p(n^{-1/4})\) MSE速率）估计。
  2. 方法：采用“plug-in调整”：将响应 \(Y\) 和线性协变量 \(X\) 分别对 \(Z\) 做非参数ML回归得调整变量 \(\tilde{Y}, \tilde{X}\)，然后对 \((\tilde{Y},\tilde{X})\) 拟合标准线性混合效应模型，用GLS/REML估计 \(\beta\)。
  3. 主要结论：估计量 \(\hat{\beta}\) 以 \(\sqrt{n}\) 速率收敛，渐近正态，且达到半参效率界（即渐近方差等于已知 \(f\) 时 \(n^{-1}\) 信息的逆）。模拟和实证例子验证了有限样本有效性。

关键设定与假设（完整版）

  • 数据结构：\(N\) 个独立个体，每个 \(i\) 有 \(n_i\) 次重复测量；个体间独立，个体内相关由随机效应 \(b_i\) 和 \(\varepsilon_{ij}\) 刻画。
  • 模型：\(Y_{ij}=X_{ij}^T\beta + f(Z_{ij}) + b_i^T W_{ij} + \varepsilon_{ij}\)，其中 \(\varepsilon_{ij}\sim N(0,\sigma^2)\)，\(b_i\sim N(0,G)\) 与 \(\varepsilon\) 独立；\(W_{ij}\) 为随机效应设计（通常含截距和时间）。线性部分 \(X\) 维数 \(p\) 固定。
  • 假设：
  • 可识别性：\(E[X|Z]\) 非线性可衡量，且方差矩阵正定。
  • ML算法质量：\(\hat{g}_Y\) 和 \(\hat{g}_X\)（分别对 \(E[Y|Z]\) 和 \(E[X|Z]\) 的估计）满足均方集成预测误差（MISE）= \(o_p(n^{-1/2})\)；推断只需要 \(o_p(n^{-1/4})\)，但本文假设更强以直接获得效率（典型DML条件）。
  • 正则性：\(f\) 足够光滑使得参数速率可达到；随机效应协方差 \(G\) 一致可估且满足标准线性混合模型正则条件。
  • 相比已有文献的调整：与P-spline方法（Ruppert et al. 2003）相比，本文取消了 \(f\) 必须以某基展开并加惩罚的要求，允许使用现代ML；与i.i.d DML相比，本文明确建模了组内相关性，并在调整中保留了随机效应部分的全部结构，而非简单用交叉拟合处理依赖。

主要结果（理论型）

  • 定理1（收敛性与渐近正态）：在假设条件下，\(\hat{\beta}\) 满足
    \[\sqrt{n}(\hat{\beta}-\beta) \xrightarrow{d} N(0, V_{\text{eff}}),\]
    其中 \(V_{\text{eff}} = \left( \sum_{i=1}^N X_i^{*T} \Sigma_i^{-1} X_i^* \right)^{-1}\)，\(X_i^* = X_i - E[X_i | Z_i]\)（调整后的 \(X\) 残差矩阵），\(\Sigma_i\) 是第 \(i\) 个个体的正确协方差矩阵（由 \(G\) 和 \(\sigma^2\) 决定）。证明关键：两步调整后，估计方程的一阶影响来自 \(\tilde{X}\) 部分，ML误差被控制在 \(o_p(n^{-1/2})\) 内，因此不影响 \(\hat{\beta}\) 的渐近分布。
  • 定理2（半参有效性）：\(\hat{\beta}\) 的渐近方差等于在已知 \(f\) 时使用GLS估计 \(\beta\) 的Cramér-Rao下界。这意味着非参数估计 \(f\) 在渐近意义下没有效率损失。
  • 技术难点：证明中需要处理个体内误差相关性对“偏调整+后续GLS”估计的影响，确保plug-in不会引入超出 \(o_p(n^{-1/4})\) 的偏差。作者通过引入“调整后的得分函数”并用经验过程理论控制剩余项，推广了Neyman正交得分的思想到随机效应模型。

证明路线与技术技巧

  • 整体路线（3-5步）：
  • 调整步骤：定义 \(Y^{*}_{ij} = Y_{ij} - E[Y_{ij}|Z_{ij}]\)，\(X^{*}_{ij} = X_{ij} - E[X_{ij}|Z_{ij}]\)。则原模型变为 \(Y^{*}_{ij} = \beta X^{*}_{ij} + U_{ij}\)，其中 \(U_{ij}=b_i^T W_{ij} + \varepsilon_{ij}\)。该步展示理想调整后 \(\beta\) 可由线性混合模型有效估计。
  • plug-in近似：用ML估计 \(\hat{E}[Y|Z]\) 和 \(\hat{E}[X|Z]\) 替换真实条件期望，得 \(\hat{Y}^{*}, \hat{X}^{*}\)。作者证明若ML估计的MISE = \(o_p(n^{-1/2})\)，则 \(\hat{Y}^{*} - Y^{*} = o_p(n^{-1/4})\)（类似DML中的“small bias”条件）。
  • 构造估计方程：对调整后数据 \((\hat{Y}^{*},\hat{X}^{*})\) 用标准线性混合模型REML得到的 \(\hat{\beta}\) 等价于解估计方程
    \[\sum_{i} \hat{X}_i^{*T}\hat{\Sigma}_i^{-1}(\hat{Y}_i^{*} - \hat{X}_i^{*}\beta) = 0,\]
    其中 \(\hat{\Sigma}_i\) 是 \(G,\sigma^2\) 的估计。
  • 渐近展开：将上述方程围绕真实 \(\beta\) 和真实期望展开，使用泰勒展开和切距比，将plug-in误差、协方差估计误差分离为主项 + 可忽略项。关键引理证明协方差估计 \(\hat{G},\hat{\sigma}^2\) 是 \(\sqrt{n}\)-相合的，从而 \(\hat{\Sigma}_i^{-1}\) 替换 \(\Sigma_i^{-1}\) 不影响一阶渐近。

  • 验证效率：展开得到 \(\sqrt{n}(\hat{\beta}-\beta) = \left( \sum_i X_i^{*T}\Sigma_i^{-1}X_i^* \right)^{-1} \frac{1}{\sqrt{n}} \sum_i X_i^{*T}\Sigma_i^{-1}U_i + o_p(1)\)，由中心极限定理得渐近正态，且方差即为半参有效方差。

  • 关键跳跃点：证明ML估计的 \(o_p(n^{-1/4})\) MSE足以保证 \(\hat{\beta}\) 达到 \(\sqrt{n}\) 速率且无渐近偏差。这依赖于“偏调整”的Neyman正交性质：在随机效应模型中，\(\beta\) 对应的得分函数对 \(f\) 的Gateaux导数为零，因此 \(f\) 的估计误差在 \(\beta\) 推断中只有二阶效应。

  • 技术技巧：
  • 经验过程（empirical process）：用于控制 \(\hat{E}[X|Z]\) 与 \(E[X|Z]\) 的差值在 \(L_2\) 范数下的一致收敛率。
  • 切切距比（tetra-lemma）：处理ML估计的随机偏差项与随机效应之间的交叉内积。
  • Louis’ / Lindeberg-type CLT：结合delta方法和弱收敛定理。

真实例子与应用

论文使用来自HIV感染者的纵向数据集（ACTG 175或类似公开数据），包含多次随访的CD4计数、治疗组指示、基线变量等。他们将 \(Y\)（CD4细胞计数）、\(X\)（治疗组指示+时间交互）、\(Z\)（基线年龄/体重等非线性影响连续变量）设定为半参数模型，用随机森林估计 \(E[Y|Z]\) 和 \(E[X|Z]\)，然后拟合线性混合模型。得到 \(\hat{\beta}\) 的置信区间比P-spline方法更窄，覆盖概率更接近名义水平。该例子展示了：（1）方法可直接套用现有软件（R包dmlalg）；（2）与P-spline相比，ML调整额外捕获了非线性交互作用，实现更好的调整，减少了 \(\beta\) 估计的偏倚和方差。

🔎 结论是否比证明窄

总体一致。但需注意：定理证明要求ML估计的MISE达到 \(o_p(n^{-1/2})\)，而当 \(Z\) 维数超过3或 \(f\) 非常粗糙时，该条件未必成立。作者在文末提到“对高维 \(Z\) 可考虑稀疏加性结构”，但未证明在该情形下本文方法能否保持效率。因此，凡论文中声称“适用于任意ML算法”的语句，严格来说受限于假设的收敛率条件，并非完全任意。

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

  1. 高维非线性部分：当 \(Z\) 维数增大（如 \(d > \log n\)）时，一般ML估计难以达到 \(o_p(n^{-1/4})\) 的MSE。本文未讨论二阶段推断是否仍有效。可参考DML在稀疏非参数模型中的推广（如加性结构）（扎根于定理假设“MISE = o_p(n^{-1/2})”和文末future work的提及）。

  2. 异质性随机效应结构：本文假设随机效应只影响线性部分、且协方差结构参数固定。若随机效应的设计矩阵 \(W_{ij}\) 也包含 \(Z\) 的维度（如随机斜率随 \(Z\) 变化），调整步骤会如何？（扎根于模型设定“\(U_{ij}=b_i^T W_{ij}+\varepsilon_{ij}\)”中 \(W\) 不含 \(Z\)）。

  3. 非高斯随机效应/误差：若 \(b_i\) 或 \(\varepsilon_{ij}\) 非正态，GLS估计可能不再达到半参有效（在非参数意义上）。是否仍可通过加权伪似然达到一定效率？（扎根于推论假设正态性）。

  4. 计算与理论对称性：论文使用ML算法调整 \(Y\) 和 \(X\)，但未讨论调整算法计算复杂度对整体推断的影响。“离散的”机器学习算法（如随机森林）理论上难以严格验证MISE条件；是否存在一类更适合理论证明的近似算法（如保序核回归）？（可延伸至作者已开发R包dmlalg，但未讨论算法选择的理论性质）。

Maintained by 陈星宇 · Homepage · Source on GitHub

评论