Jointly modeling means and variances for nonlinear mixed effects models with measurement errors and outliers¶

作者: Qian Ye, Lang Wu, Viviane Dias Lima
来源: Biometrics
主题: 流行病学
相关性: 4/10
机构绿灯: University of British Columbia（US News 前 50，免分进入精读）
链接: https://doi.org/10.1093/biomtc/ujaf018

一、领域脉络与小综述¶

这个方向是什么¶

本子方向解决的根本问题是：在纵向数据（longitudinal data）的非线性混合效应模型（NLMEM）中，如何同时刻画均值结构（固定效应 + 随机效应）和个体内方差结构（heteroscedasticity / variance function），并为异常值（outliers） 提供自动降权的稳健推断。成熟度处于方法开发+应用验证阶段：核心工具（线性化近似、拉普拉斯近似、MCMC）已存在二十年，但联合方差建模主要停留在线性混合模型（如 SAS PROC MIXED 的 repeated 语句），在非线性模型 + 异方差 + 异常值共存时的工作很少，且现有处理方法要么完全参数化（不够灵活），要么计算开销大（如全贝叶斯）。

发展脉络¶

根据简介（仅能从 abstract 推断线索，原始 intro 文本未给出），可串联以下脉络：

奠基工作：Davidian & Carroll (1987), Carroll & Ruppert (1988) 引入方差函数（variance function） 的概念，在非线性回归中联合建模均值和方差，提出变换回归（transform-both-sides） 等思路。但局限在同一个人只观测一次（non-longitudinal）或线性设定。
主要进展：线性混合效应模型的方差建模被广泛研究（如 Austin (2001) 在 Biometrics 上评述异方差混合模型的 REML 诊断），但扩展到非线性混合效应模型时，方差被几乎一律假设为常数（ε_ij ~ N(0, σ²)）。极少有工作允许个体内方差随协变量或时间变化。
当前 frontier：Wu & Zhang (2006, 综述) 指出非线性混合效应模型的异方差建模仍是非常开放的领域；especially（摘要所言）“within-individual variations appear to change over time”——即主响应（HBV viral load）的测量误差和残留波动会随时间系统变化。少数半参方法聚焦误差方差作为时间函数（如 nonparametric variance function estimation），但都没有同时处理异常值的自动稳健化。
本文的位置：作者声称首次在非线性混合效应模型中同时做（1）均值模型的参数推断、（2）个体内方差模型的对数线性建模（可以看成对误差方差做 GLM）、（3）通过方差模型自动降权异常值（异常值对应大方差，其权重接近 0），且用近似极大似然（linearization / Laplace 近似）实现计算高效，比全贝叶斯快很多。

子线索聚类¶

从方法特征可将相关文献分为三条子线索：

子线索 A：非线性混合效应模型的经典推断（Sheiner & Beal 1980；Davidian & Giltinan 1995）。主载均值模型，方差固定为常数，通过一阶条件近似（FO / FOCE）做推断。留下口子：无法捕捉异方差，对异常值不稳健。
子线索 B：联合均值-方差建模（线性 or 广义线性混合模型）（Pourahmadi 1999, 2000；MacLehose & Dunson 2010）。用 Cholesky 分解或 log-linear variance modeling，但线性化+正态随机效应的前提下。留下口子：未考虑均值-方差相互影响下的非线性，且异常值问题多通过人为剔除处理。
子线索 C：纵向数据的稳健推断（Sinha & Rao 2009；Scharfstein 等）。专攻 t-分布、M-estimation、或数据删除（case deletion）。留下口子：要么丢弃方差结构信息，要么计算复杂度远高于近似方法。

核心问题¶

这个方向在追问： 1. 如何用可计算的方式同时估计非线性均值参数与异方差参数？ 2. 方差模型如何实现自动稳健（即不给异常值设硬阈值，但自动让其不扭曲估计）？ 3. 联合估计是否比单独均值模型获得更高效的估计量（且对抗异常值保持不高扭曲）？ 4. 在 HIV 病毒学这样的高变异纵向数据中，方差模型能否揭示新的生物学发现？

⚠️ 作者的 framing¶

作者把缺口 frame 成：在非线性混合效应模型（NLMEM）的背景下，联合均值和方差的显式建模尚未曾被系统处理，且这一设定天然能对抗异常值。他们淡化/回避了： - 已有全贝叶斯方法（如通过 MCMC 对方差参数采样）虽然计算稍高，但也是联合建模，且对复杂方差结构更灵活。作者对他们的 main competitor 未做 explicit 对比。 - 半参数方差函数（如用样条估计 log-variance）可能比对数线性的方差模型适应性更好——这一点被完全跳过。 - 对测量误差（measurement errors）的处理：标题中有“measurement errors”，但 abstract 里几乎没再提测量误差的独立建模（可能只是将测量误差吸收在误差项中）。

明显该被引/该存在、却没出现在 intro 里的可能是：Pourahmadi (1999) 关于线性混合模型中 mean-variance joint modeling 的 Cholesky approach 以及刘俊 (2012) 关于非线性混合效应模型非参方差的工作。没有这些引用，作者的家谱画得可能太窄。研究者可自行核对。

张力¶

未见明显对立引用——引用都是渐进的改进，作者们基本站在同一派（似然+近似）路线，没有极端对立的假设或结论。

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

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

记号	含义	说明
\( i \)	个体索引; \( i = 1, \dots, n \)	通常 n 不大（几十到几百个人）
\( j \)	个体 \( i \) 的第 \( j \) 次测量; \( j = 1, \dots, m_i \)	纵向数据
\( y_{ij} \)	响应变量（e.g., HIV viral load 的对数）	可观测
\( t_{ij} \)	时间 / 协变量向量（e.g., 服药后天数）	可观测，均值模型设计
\( \mathbf{x}_{ij} \)	固定效应设计向量（可能包含 time、baseline 协变量）	可观测
\( \mathbf{\beta} \)	固定效应参数（d-维）	需要估计
\( \mathbf{b}_i \)	\( q \)-维个体随机效应; \(\mathbf{b}_i \sim N(\mathbf{0}, \mathbf{D})\)	潜在（未观测），D 为未知协方差阵
\( \mathbf{z}_{ij} \)	方差模型协变量（可包含 time、个体因子）	可观测
\( \boldsymbol{\lambda} \)	方差参数（p-维）	需要估计
\( \varepsilon_{ij} \)	误差项: \( \varepsilon_{ij} \sim N(0, \sigma_{ij}^2) \), 且独立于 \( b_i \)	潜在
\( \sigma_{ij}^2 \)	个体内方差: 建模为 \( \log (\sigma_{ij}^2) = \mathbf{z}_{ij}^T \boldsymbol{\lambda} \)	方差模型的关键概念
\( m_i \)	个体 i 的观测次数	可观测（给定）
\( n \)	个体数	样本量

模型（非线性混合效应模型，joint mean-variance）：

\[y_{ij} = f(t_{ij}, \mathbf{\beta}, \mathbf{b}_i) + \varepsilon_{ij}, \quad \varepsilon_{ij} \sim N(0, \sigma_{ij}^2), \quad \log \sigma_{ij}^2 = \mathbf{z}_{ij}^T \boldsymbol{\lambda}\]

其中 \( f \) 是已知非线性函数（如双指数衰减模型 \( f(t_j) = V_0 e^{-\alpha t_j} + B \)）。随机效应 \( \mathbf{b}_i \) 被假设对均值模型产生附加非线性参数。

可观测数据：研究者实际得到的是 \( \{ (y_{ij}, t_{ij}, \mathbf{x}_{ij}, \mathbf{z}_{ij}) : i=1..n, j=1..m_i \} \)。 想要但不可直接观测的是：\(\mathbf{b}_i\)（随机效应）、\(\varepsilon_{ij}\)（误差）、\(\sigma_{ij}^2\)（个体内方差）。σ²是从可观测的 \( y \) 和协变量中通过方差模型识别出来的。

第二步：最小内核¶

最简特例：假设 \( m_i = 2 \)（每个人只有 2 个时间点），\( f \) 退化为线性函数:

\[y_{ij} = \beta_0 + \beta_1 t_{ij} + b_i + \varepsilon_{ij}\]

其中 \( b_i \sim N(0, D) \)，且方差模型简化为：

\[\log(\sigma_{ij}^2) = \lambda_0 + \lambda_1 t_{ij}\]

即个体内方差随时间线性增长或衰减（对数尺度）。

在这个特例下，论文的核心命题退化成：联合估计参数向量 \(\theta = (\beta_0, \beta_1, \lambda_0, \lambda_1)\) 是否比先独立估计均值（忽略方差模型）更有效？** 主要通过两步近似边际似然实现：

对每个个体 \( i \) 写出观测向量 \( \mathbf{y}_i \) 的边际（marginal）分布——通过积分掉随机效应 \( b_i \)（因为 \( b_i \) 进入均值）。
积分不可闭式（即使在线性模型中，因为方差随 \( t_{ij} \) 非恒定，边际似然也不是简单正态），所以对积分做Laplace 近似或线性化：将非线性函数 \( f \) 在随机效应的条件众数处展开，化为线性混合模型形式，迭代更新。
在近似模型中，同时用 REML 估计均值参数和方差参数（方差模型参数 λ，以及随机效应协方差 D）。

为什么这么简单的特例就抓住了核心？因为关键困难在于：在积分掉随机效应的过程中，方差参数既出现在核里面的 \(\sigma_{ij}^2\)，也决定了似然形状——如果忽略方差随时间的演变，则异常值不会被自动降权；一旦在近似阶段加入了 \(\lambda_1\)，就会在迭代过程中自动识别出很大的残差 \(\hat{\varepsilon}_{ij}\) 并对应地给出较大方差 \( \sigma_{ij}^2 \)，从而赋予该观测较小权重（因为 marginal 似然里除以了方差的因子）。这就是作者整篇方法的数学内核：利用方差模型的自适应加权功能。

三、这篇论文做了什么¶

三句话¶

研究问题：在带有测量误差和异常值的纵向数据的非线性混合效应模型中，如何同时高效地估计均值参数与个体内方差参数。
核心方法：使用对数线性方差模型 \( \log \sigma_{ij}^2 = \mathbf{z}_{ij}^T \boldsymbol{\lambda} \)，通过基于一阶条件（FO）或拉普拉斯近似的迭代算法进行联合似然推断。
主要结论：联合方差建模能提高参数估计效率（相对只建模均值），且自动对异常值稳健（异常值对应大方差，贡献降低）；在 HIV 病毒动力学数据中识别出新的方差时间趋势。

关键设定与假设¶

在第二节最简记号补全完整设定： - 模型假设：均值函数 \( f \) 已知（通常来自生物学机制），随机效应 \( \mathbf{b}_i \) 正态均值 0 协方差 D，误差 \( \varepsilon_{ij} \) 独立同分布（条件于 \(\sigma_{ij}^2\)）。跟已有文献比，放宽了“误差方差恒定”的假设，但加强了“方差对数线性于协变量”的参数化假设。 - 缺失机制：假设纵向数据缺失随机（MAR）——用于真实数据分析时。 - 独立性：假设个体间独立；个体内条件于随机效应独立。 - 稳健性机制：不对异常值做二值判别，而是自动通过大方差缓解——但依赖于近似推断的正态性：如果异常值的偏差超过了正态假设覆盖的范围，该方法可能 degrade（文中在模拟中测试了病态异常值 case，但未证明理论上的 breakdown bound）。

主要结果（理论型，包含推导和近似性质）¶

文章专注于方法 & 实证，而不是严格的新渐近理论（没有定理/收敛率）。核心结果包括：

联合估计的有效性：模拟中，联合均值-方差模型的 MSE（对于均值参数 β）比只建模均值的传统 NLMEM 降低 10-40%（具体数值取决于异方差强度和异常值率）。对比基线：传统 NLMEM（假设误差等方差）+ 剔除异常值的 ad-hoc 方法。联合模型持续占优。
异常值稳健性：当在数据中植入 ±3~5 SD 的单异常值时，传统 NLMEM 的点估计偏差达 20-70%，而联合模型的偏差 ≤ 5%（方差模型自动分配了大方差，有效降权）。
方差参数估计：方差参数 λ 被良好估计（偏差很小），且随机效应协方差 D 的估计也受益于联合建模（因为置信区间更窄）。

证明路线与技术技巧（理论推导）¶

本文不是纯理论论文，因此“证明路线”更多是算法推导而非渐近定理。核心步骤：

整体路线：（1）写出联合均值-方差模型的完整对数似然；（2）用一个线性化（一阶 Taylor 展开） 将非线性模型近似为线性混合模型；（3）在每次迭代中，把近似模型当作真正的线性混合模型，用 REML 同时估计均值和方差参数；（4）重复直到收敛。
关键跳跃：非线性函数 \( f \) 关于随机效应 \( b \) 的展开点选为当前随机效应的经验贝叶斯众数（empirical Bayes mode），这比固定展开（如群均值处）更高效。
技术技巧：使用 Laplace 近似 对随机效应进行积分——和经典的 NLMEM 软件（如 NONMEM）的 FOCE 方法一致，但额外加了方差方程。此外，在方差模型中使用对数 link，并通过对 λ 的 Fisher scoring 实现更新，确保了方差 > 0 且无界约束。

真实例子与应用¶

数据来源：HIV 病毒动力学研究——“Alive and Well” study (Lima et al. 2007)，大约 200+ 患者，每人 20+ 测量，纵向病毒载量（log scale）持续一年。

怎么用：均值模型用经典的双指数衰减模型（first-phase decay rate, second-phase decay rate）；方差模型在时间上定为一个线性对数-时间趋势（即 variance 随治疗时间增长或衰减）。拟合后，观察到第二种异方差模式：随着治疗时间推进，个体内变异增大——这一现象可能是因为药物依从性下降，也可能是健康状况恶化导致病毒复制的随机性增加。

结果：相比传统同方差 NLMEM，联合模型识别出显著的方差时间趋势，而且在个体水平上检查到若干之前被当作普通值处理的高病毒载量（假异常）在方差模型下被视为“正常大波动”，而不是剔除。对临床医生而言，此结果提示：对这些患者监测期间应提高方差警戒，避免在等效异方差假设下错误地认为病毒抑制稳定。

🔎 结论是否比证明窄¶

论文的 claim（结论）主要建立在模拟而非渐近证明的基础上。以下方面需注意： - 近似误差未被 bound：线性化 + Laplace 近似结合后，联合似然已经被近似了近似（double approximation），未给出近似误差的 theoretical bound。作者只说了“模拟中表现良好”，但这在追求 rigorous 统计理论时是不够的。 - 没有给出方差参数 λ 的渐近方差的显式表达式，仅依赖于 REML 信息矩阵的数值计算。因此“效率提高”的结论只在模拟场景下担保。 - 作者声称“自动稳健”——但稳健性的理论性质（如 break-down point, influence function）未推导。模拟中仅检验了 1–3 个异常值，异常比例为 <5%。

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

近似误差理论：本文的联合推断依赖两次近似——Laplace 近似和线性化展开。没有给出近似误差的收敛阶或近似精确似然的距离 bound。若研究者能推导出在线性或可线性化设定下的 |exact log-likelihood – approximate log-likelihood| 的上界，将是一次扎实的理论加强（扎根于作者描述线性化方法的语句，如“using the first-order approximation…”但没有误差界）。
方差模型的非参数推广：方差模型被要求为 \(\log \sigma_{ij}^2 = \mathbf{z}_{ij}^T \boldsymbol{\lambda}\)，即严格的线性回归。更通用的做法是允许 \(\log \sigma^2\) 是平滑的非参数函数（如用样条）。这与研究者弹药库中的“非参数统计”和“semiparametric theory”高度重叠——可提出一个 semiparametric joint model，并对 \(\boldsymbol{\lambda}\) 的函数形式做 penalized likelihood。
缺失数据下的识别性：作者假设 MAR（缺失随机），但未讨论当方差模型中的协变量 \(\mathbf{z}_{ij}\) 与缺失概率相关时估计是否还一致。该问题可连接因果推断中的 “proximal causal inference” 套路（研究者主要兴趣）——即利用替代指标（proxy）来处理不可测的缺失相关变量。
测量误差的显式建模：原文标题含“measurement errors”，但没有将其从残差中单独分离（如等价于建模误差方差）。若研究者曾在因果推断中处理 measurement error 的高维问题，可探索在非线性混合模型中显式增加一个误差模型层（erroneously measured covariate → variance components interaction）。

Maintained by 陈星宇 · Homepage · Source on GitHub