跳转至

Variable selection for high‐dimensional generalized linear model with block‐missing data

作者: Yifan He, Yang Feng, Xinyuan Song
来源: Scandinavian Journal of Statistics
主题: 高维统计 / 随机矩阵
相关性: 7/10
机构绿灯: Chinese University of Hong Kong（US News 前 50，免分进入精读）
链接: https://doi.org/10.1111/sjos.12632

一、领域脉络与小综述（≥25%）

这个方向是什么

本方向解决的根本问题是：当高维数据的协变量矩阵呈现结构化块缺失（block-wise missing）——即数据缺失以整块变量为单位、而非随机分散在个别观测中——如何同时完成变量选择（识别哪些协变量对响应有影响）与参数估计，并且保证所选模型在统计意义上一致（变量选择一致性、估计一致性）。当前成熟度处于理论初步建立、方法有限、缺失机制假设较强的阶段：已有方法大多只能处理“单块缺失”（整个数据仅一个块缺失）或严重依赖线性模型结构，缺乏在广义线性模型（GLM）框架下对多块缺失的通用解决方案。

发展脉络（history）

将引言中引用的工作串成一条线：

  • 奠基工作：Little & Rubin (2019, 2002) 建立了缺失数据分析的基础框架（缺失机制分类：MCAR、MAR、MNAR），以及基于插补、EM算法、直接似然的通用策略。这是所有缺失数据处理方法的起点。
  • 主要进展（针对高维缺失）：
  • Lounici (2014) 及一系列工作（如Cai & Zhang, 2016; He et al., 2022）将高维协方差矩阵/精度矩阵估计与国际象棋棋盘模式的块缺失结合，发展了基于对缺失模式进行结构建模的估计方法。这些方法依赖于线性/高斯模型，且大多只处理单块缺失。
  • Vo et al. (2020) 提出了一个多项式回归插补方法用于块缺失数据，但作者指出“该方法虽然融入了协变量与响应的关系，但由于使用了多项式函数，只能近似响应与协变量的真实联系，可能引入模型错误设定”。
  • 当前前沿：作者引用了Yan & Bien (2021) 处理多块缺失的线性模型变量选择工作，以及He et al. (2022) 对多块缺失下的稀疏精度矩阵估计工作。但作者指出：“现有方法要么只适用于单块缺失模式，要么过度依赖模型结构，在更一般的广义线性模型框架下仍留有空白。”
  • 本文位置：本文填入这个空白的方法是：利用多块缺失的结构先估计稀疏精度矩阵，然后基于观测块对缺失块进行条件均值插补（即回归插补），最后在插补后的完整数据上运行GLM变量选择，并证明该流程在渐近意义上能保持变量选择一致性和估计一致性。

子线索聚类

  1. 精度矩阵估计 + 块缺失（方法簇）：He et al. (2022), Lounici (2014), Cai & Zhang (2016)，聚焦在利用块缺失结构恢复协方差或精度矩阵的稀疏模式。局限：它们通常得到一个矩阵估计，但不直接处理变量选择；且需要高斯假设。
  2. 变量选择 + 缺失数据（方法簇）：Yan & Bien (2021), Vo et al. (2020)，直接在缺失数据下进行惩罚回归或插补后回归，但仅适用于线性模型。本文试图将其推广到GLM。
  3. 一般插补方法（如多重插补）在高维下的一致性分析：Raghunathan et al. (2001), van Buuren (2018) 提出方法，但作者认为这些方法“未在高维GLM情境下验证变量选择一致性”。

核心追问与瓶颈

  • 核心问题1：块缺失结构下，是否可以仅依靠条件均值插补（而非多重或热卡插补）就达到变量选择一致性？条件需要多强？
  • 核心问题2：块缺失下的精度矩阵估计的一致率（convergence rate）与后续GLM变量选择的一致性之间的传播关系如何？
  • 当前瓶颈：块缺失模式的非随机性（往往是非MAR，取决于其他块的值）使得MAR假设效率受限；且多数理论工作只针对线性模型。

⚠️ 作者的 framing（必须明确标注成“这是作者的说法”）

  • 作者把缺口 frame 成：“现有方法要么只适用于单块缺失（如Lounici 2014），要么依赖于模型结构（如Yan & Bien 2021为线性模型）；我们提出一种在GLM框架下对多块缺失通用的方法，并通过稀疏精度矩阵估计的插补实现变量选择一致性。”
  • 被淡化的竞争路线：作者未深入讨论多重插补（MI）在高维GLM下的表现，也未提及直接基于缺失似然（如ESM算法）的高维变量选择方法。文中仅一句提到MI“计算代价高且不保证变量选择一致性”，但未引用相关理论。
  • 什么明显该被引却没出现：未见引用关于double robustness（双稳健）插补方法（如Robins et al., 1995），也未引用high-dimensional semiparametric methods for missing data（如van der Laan & Rose, 2011, Targeted Learning），更未引用proximal causal inference中处理隐变量缺失的方法（与块缺失在结构上有类似性）。这些都是值得研究者去查的问题：它们要么提供了更灵活的缺失机制假设，要么在不同结构下解决了高维变量选择。
  • 张力：未见明显对立引用。但作者在实际例子上声称“优于仅适用于单块缺失的方法”——若验证尺度不统一，则此论断可能不够强。

二、最核心、最简单的例子 / 数学问题（≥15%）

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

  • 符号：
  • \(X \in \mathbb{R}^{n \times p}\)：协变量矩阵，\(n\)为样本量，\(p\)为协变量个数（高维情形 \(p \gg n\)）。
  • \(Y \in \mathbb{R}^n\)（或 \(\{0,1\}^n\)）：响应变量向量（连续或二值）。
  • \(\beta \in \mathbb{R}^p\)：回归系数向量，是目标估计量。
  • \(M\)：缺失模式指示矩阵，\(M_{ij}=0\)表示第\(i\)个样本第\(j\)个协变量缺失，\(=1\)表示观测。块缺失：存在若干块集合 \(B_1, B_2, \dots, B_K \subseteq \{1,\dots,p\}\)，使得对每一个块，若一个样本在该块的任意变量缺失，则该块所有变量都缺失（即 \(M_{i,j}=0\) 对所有 \(j \in B_k\) 同时成立）。
  • \(\Omega = \Sigma^{-1} \in \mathbb{R}^{p \times p}\)：协变量 \(X\) 的精度矩阵（假设 \(X\) 服从多元高斯分布）。\(\Omega\) 是稀疏的（大多数元素为0），用于插补缺失块。
  • \(\mathbb{E}[X_i | X_{i,\text{obs}}]\)：给定观测块后缺失块的条件均值，用 \(\Omega\) 的分块结构计算。
  • 模型：
  • 广义线性模型（GLM）：\(\mathbb{E}[Y_i | X_i] = \mu(X_i^\top \beta)\)，其中 \(\mu(\cdot)\) 是一个已知的链接函数（如 logistic 回归的 sigmoidal 函数）。\(Y_i\) 的条件分布在给定 \(X_i\) 后属于指数族（如 Bernoulli, Poisson, Gaussian）。
  • 假设：\(X_i \stackrel{i.i.d}{\sim} \mathcal{N}(0, \Sigma)\)，且 \(\Sigma\) 的精度矩阵 \(\Omega\) 是稀疏的（非零元素的个数远小于 \(p^2\)）。缺失机制为 MAR（缺失依赖观测到的协变量和结果，但不依赖缺失的协变量值本身），更精确地说，块缺失的缺失概率仅依赖于观测块的变量值。
  • 可观测数据：
  • 研究者实际能观测到的是：\(\{(Y_i, X_{i,\text{obs}})\}_{i=1}^n\)，其中 \(X_{i,\text{obs}}\) 是 \(X_i\) 中被观测到的子向量（缺失块整体未观测）。
  • 想要但观测不到的是：缺失块 \(X_{i,\text{mis}}\) 的具体取值。这些值通过条件均值插补（基于 \(\Omega\)）来填充。填充后的完整数据记为 \(\tilde{X}_i\)。

第二步：讲最小内核

最简特例： - 设定：\(p=2\)，协变量只有两块：\(X_{1}\)（完全观测）与 \(X_{2}\)（某些样本的整个块缺失）。响应 \(Y\) 服从 logistic 回归 \(\mathbb{P}(Y=1|X) = 1/(1+e^{-X^\top\beta})\)。符号：\(\Omega = \begin{pmatrix} 1 & \rho \\ \rho & 1 \end{pmatrix}\)（标准化后）。缺失机制：\(X_2\) 的缺失仅依赖于 \(X_1\)（MAR）。可观测数据：对每个样本 \(i\)，要么观测到完整的 \((X_{i1}, X_{i2})\)，要么只观测到 \((X_{i1}, Y_i)\)。 - 本文的核心思路退化为： 1. 估计 \(\Omega\)：利用观测到 \(X_2\) 的样本来估计 \(\rho\)（比如通过样本相关系数），剩下的样本无法直接估计。但块缺失结构允许我们使用所有样本中 \(X_1\) 的完整数据与部分样本中 \(X_2\) 的数据来联合估计 \(\rho\)。这相当于：用所有 \(n\) 个样本的 \(X_1\) 方差估计（为1，已知），以及 \(\approx n_2\) 个有完整数据的样本的协方差估计，通过最大化似然或求解Yule-Walker方程来估计 \(\rho\)。若 \(n_2\) 足够大且缺失不严重，则估计的一致成立。 2. 回归插补：对缺失 \(X_2\) 的样本，用 \(\hat{\mathbb{E}}[X_2|X_1] = \hat{\rho} X_1\) 填充 \(X_2\)。这实质上就是将缺失块替换为基于观测块的线性预测。 3. 变量选择：在填充后的完整数据 \(\{Y_i, (X_{i1}, \hat{X}_{i2})\}\) 上，运行带有 \(\ell_1\) 惩罚的logistic回归（LASSO），并证明变量选择一致性（即当\(n\to\infty\)，正确地将系数非零的变量保留在模型中，零系数变量剔除）。

为什么这个例子抓住了核心：在这种两变量、单一缺失块的特例下，所有一般性都退化为简单设定，但本文的三个关键步骤（精度矩阵估计、条件均值插补、惩罚GLM）依然完整呈现，而且插补的准确性直接取决于 \(\hat{\rho}\)的收敛速度。论文的数学困难在于：当 \(p \gg n\)、有多个块时，\(\Omega\) 的估计（通过graphical lasso）不再简单，但其本质与这个特例是一致的——利用块缺失结构来联合估计所有块间的协方差。

三、这篇论文做了什么（≥45%）

三句话概括

  • 研究问题：对多块缺失数据（block-wise missing data），在广义线性模型框架下，如何同时实现变量选择与参数估计的一致性？
  • 核心方法：提出一个两阶段方法——先利用块缺失结构估计稀疏精度矩阵（graphical lasso + 对缺失模式修正的似然）进行条件均值插补，再在插补后的完整数据上运行LASSO惩罚GLM进行变量选择。关键在于插补本身不依赖响应变量（\(Y\)），从而避免了因插补引入内生性偏差。
  • 主要结论：在MAR缺失机制、稀疏精度矩阵假设、以及一定的正则性条件下，所提方法达到了变量选择一致性和参数估计一致性；且当块缺失比例不过高时，其渐近率与完全观测数据的LASSO GLM相当。

关键设定与假设

（补充第二节未完整列出的重要假设，在完整记号下）

  • 假设1（缺失机制）：\(P(\text{块缺失} \mid X, Y) = P(\text{块缺失} \mid X_{\text{obs}})\)，即MAR（随机缺失），缺失仅依赖于观测到的协变量。相比已有文献：这是弱于MNAR但强于MCAR的常见假设。而现有工作（如Yan & Bien 2021）要求MCAR或非常近似。
  • 假设2（精度矩阵稀疏性）：\(\Omega\) 的非零对角线与上三角元素个数为 \(s\)，且有 \(s = O(n / (\log p)^{3/2})\) 或类似条件。这是graphical lasso一致性的标准条件。
  • 假设3（GLM链接函数的光滑性与对数似然的强凸性）：链接函数二阶可微、Hessian矩阵在真实 \(\beta_0\) 附近有正定界、且为强凸。这是LASSO GLM一致性所需条件。
  • 假设4（缺失率有界）：在任一块的缺失率（缺失样本占总样本比例）有界于一个小于1的常数。这是为了确保精度矩阵的样本协方差矩阵估计可逆（经过缺失调整后）。
  • 假设5（分块结构）：缺失模式遵循预先定义好的块划分（如来自不同研究的变量列表），且块之间不重叠。这是块缺失的本质结构。

主要结果

挑2个最关键的理论结果：

  • 定理1（插补后协方差一致性）：设 \(\hat{\Omega}\) 是通过修正的graphical lasso（考虑了块缺失）得到的精度矩阵估计，则存在常数 \(C>0\) 使得以高概率有 \(\|\hat{\Omega} - \Omega\|_2 \leq C \sqrt{\frac{s \log p}{n}}\)。直觉：插补使用了正确的结构（插补本身不引入新不确定性），所以精度矩阵的收敛率与完全观测数据下的graphical lasso几乎一致（只有常数因子增加，因为缺失减少了有效样本量）。必要条件：缺失块的数量 \(K\) 固定且小（不是无穷多块），否则有效样本量太少。
  • 定理2（变量选择一致性）：在插补后的数据上，LASSO惩罚GLM（选择合适的正则化参数 \(\lambda \asymp \sqrt{\frac{\log p}{n}}\)）所筛选的模型以概率趋于1等于真实稀疏模型，且估计误差 \(\|\hat{\beta} - \beta_0\|_2 = O_p\left(\sqrt{\frac{s \log p}{n}}\right)\)。直觉：由于插补是一致的（均方误差可忽略），插补引入的偏差不会破坏变量选择的相合性。技术难点：需要控制插补误差与LASSO自噪声的交互，并通过restricted eigenvalue condition（RE条件）来证明插补后的设计矩阵仍满足RE条件（这需要额外论证，因为插补引入了依赖关系）。

证明路线与技术技巧（理论型）

整体路线（3-5步）： 1. 精度矩阵的修正估计：利用块缺失模式构造一个可观测结构下的伪样本协方差矩阵——只使用有完整观测的交叉块元素，缺失块对应位置填入0作为占位。作者证明了此伪协方差矩阵从谱范数收敛到真实\(\Sigma\)。然后在其上运行graphical lasso得到\(\hat{\Omega}\)。（技巧：用missing data adapated covariance estimator + graphical lasso） 2. 向后回归插补的一致性：对每个缺失样本，利用\(\hat{\Omega}\)的分块矩阵公式计算条件均值 \(\hat{\mathbb{E}}[X_{\text{mis}}|X_{\text{obs}}] = - \hat{\Omega}_{22}^{-1} \hat{\Omega}_{21} X_{\text{obs}}\)（其中\(\Omega\)的分块与缺失块对应）。作者证明这个插补的L2误差以高概率有界于 \(O(\sqrt{\frac{s \log p}{n}})\)。 3. 插补后设计矩阵的RE条件：证明在插补后，观测到的设计矩阵 \(\tilde{X}\)（含插补值）仍以高概率满足LASSO的restricted eigenvalue condition，尽管插补引入了样本间的依赖（对共享缺失模式的样本，插补值使用了相同的\(\hat{\Omega}\)，因此有相关性）。作者通过coupling论证：将插补视为基于一个“几乎独立”的精度矩阵估计，进而应用经典的随机矩阵理论。 4. 变量选择一致性：利用标准LASSO GLM证明路线（如Negahban et al., 2012, van de Geer, 2008），将插补误差视为附加的approximation error，并证明它不影响一致性的主导项。关键跳跃：作者必须证明\(\|\tilde{X}\beta_0 - X\beta_0\|_2\)（插补引入的预测偏差）足够小，以至于其L2范数仍能被正则化参数所包容。 5. 最终输出：结合以上步骤，通过theorem of oracle inequality for M-estimators with imperfect covariates得出结论。

技术技巧点名： - missing-data-adapted covariance estimation：使用仅基于观测元素构造的伪协方差矩阵，然后在其上运用graphical lasso。这是He et al. (2022)引入的技术。 - 分块矩阵逆的条件均值公式：利用精度矩阵的分块直接计算条件均值，无需显式计算协方差逆，从而避免了因缺失导致的协方差矩阵求逆困难。 - Coupling + 随机矩阵理论：处理因共享插补模型带来的样本依赖，通过将插补步骤与独立估计耦合，分离出独立邻域。使用了Marchenko-Pastur law和Bai-Yin定理的变体来证明特征值边界。 - RE条件的稳健性论证：证明插补仅改变了设计矩阵的谱范数界（以概率接近常数），不改变其受限本征值性质。这用到了Weyl不等式和Slepian引理的某种形式（对随机矩阵的扰动界）。

真实例子与应用

论文使用了阿尔茨海默病神经影像学倡议（ADNI）数据。

  • 数据/场景：ADNI包含来自不同成像中心的结构MRI和PET数据。由于不同中心可能只报告某些成像特征，导致数据具有天然的块缺失结构——例如，中心A提供了海马体体积、中心B提供了皮层厚度，部分病人只在中心A成像，部分在两个中心都成像。本文以此作为多块缺失的应用场景。
  • 方法应用：作者将成像特征（约80个协变量）按模态（MRI、PET）分为两个块，构建块缺失模式。响应变量为阿尔茨海默病的状态（认知正常vs.轻度认知障碍，或认知正常vs.痴呆）。然后应用所提插补+LASSO logistic regression筛选重要成像特征。使用BIC类型准则选择正则化参数 \(\lambda\)。将方法回归插补与三种基线比较：① 仅使用完全观测样本（complete-case analysis）的LASSO；② 多重插补（MI）后的LASSO；③ 均值插补后的LASSO。
  • 结果：
  • 回归插补的变量选择获得了在验证集上更稳定的AUC（准确率），且最终模型选择了更少的变量（更稀疏），但解释力不减。
  • 基线方法（complete-case）由于样本量大幅减少，标准误很大，所选模型不稳定。
  • 均值插补表现最差，因为它完全破坏了条件独立性假设，导致变量选择一致性被破坏。
  • 例子想说明什么：验证理论预期的几种现象——① 回归插补较完备的基于模型的插补保留了变量选择一致性；② 在多块缺失场景下，它相较于仅使用完全观测数据的方法能有效利用更多样本，扩大有效样本量；③ ADNI数据块缺失结构为本文的方法提供了一个真实、清晰的应用案例。

🔎 结论是否比证明窄

  • 论文证明的是在高斯协变量且缺失机制为MAR下的变量选择一致性。但结论文中声称“对多种缺失机制（如MCAR、MAR）具有稳健性”，却未证明MNAR下的表现（尽管模拟中用了MNAR生成的情境）。因此结论比证明窄——MNAR下的稳健性属于conjecture而非定理。
  • 模拟中设置的缺失比例最高为40%，但理论证明并未明确限制缺失率上界；当缺失率接近100%（几乎整块缺失）时，插补实际上退化为猜测，结论不再成立。作者在定理中隐含了“缺失率有界于一个严格小于1的常数”，但没有定量给出可接受的上界。
  • 此外，论文的精度矩阵估计部分的证明假设所有变量是联合高斯的，而结论文中并未明确限制协变量的分布，只是说“基于块缺失的结构估计”。这意味着该结论在高斯假设之外（如\(t\)分布、混合分布）是否成立仍是开放问题。

四、开放问题（≤10%，点到为止，扎根具体语句）

以下是本文留下的、可被研究者亲自验证的开放问题：

  1. 非高斯协变量的泛化：证明中协变量 \(X\) 需要满足联合高斯假设（否则精度矩阵的条件均值插补公式不再成立）。未来的工作可以放宽为非参数的高斯copula或椭圆形分布。扎根点：定理1的假设（“\(X\) 服从多元正态分布”），以及条件均值公式 \(\mathbb{E}[X_{\text{mis}}|X_{\text{obs}}] = -\Omega_{22}^{-1} \Omega_{21} X_{\text{obs}}\) 的事实。研究者可以将非高斯扩展与自己的逆问题与非参数工具结合开发。

  2. 块缺失数目与稀疏性的交互：证明假设块数目 \(K\) 是固定且小的（定理的证明需要块间独立的结构）。若块数目随 \(p\) 增长（如 \(K \to \infty\)，每块包含少量变量），则估计收敛性会退化。扎根点：定理1的证明需要“缺失模式是预定义的、且块不重叠”，未讨论块数增长的情形。若研究者对 信息-计算权衡感兴趣，可以考察在大量块下的计算复杂度与统计效率的折中。

  3. 有限样本误差界（non-asymptotic error bounds）：论文给出的是渐近一致率（\(n \to \infty\)），但对于有限样本（如ADNI中常见的\(n \approx 400, p \approx 80\)），没有给出精确的误差界或置信区间。扎根点：定理2的尾部概率陈述：“以概率趋于1……”——但未给出有限样本的常数。这对实际应用有意义：插补后的LASSO是否能像完全数据那样有可计算的 \(\lambda\) 选择规则？研究者可利用RE条件的高概率版本，结合self-normalized bounds推导出有限样本界。

  4. 处理多分类响应或生存数据：论文仅限于二值与连续GLM，但ADNI中可能有多分类认知状态或复发事件。作者在讨论中提到“我们的方法也可推广到其他GLM”，但未给出理论或模拟。扎根点：结论部分的 “our method can be extended to other GLM without essential change”——但这并不显然，因为多分类logit模型的Ｈessian矩阵的结构复杂，插补后的变量选择一致性需要重新验证。研究者可以针对多分类或生存分析场景做专门的推广证明。

Maintained by 陈星宇 · Homepage · Source on GitHub

评论