线性混合效应模型

线性混合效应模型用于涉及依赖数据的回归分析。当处理纵向研究和其他对每个受试者进行多次观测的研究设计时，会出现此类数据。一些特定的线性混合效应模型包括

  • 随机截距模型，其中一个组中所有响应都相加一个特定于该组的值。

  • 随机斜率模型，其中一个组中的响应遵循一个（条件）均值轨迹，该轨迹在观测协变量中是线性的，斜率（以及可能的截距）根据组而变化。

  • 方差成分模型，其中一个或多个分类协变量的水平与从分布中抽取的样本相关联。这些随机项根据其协变量值相加，从而确定每个观测值的条件均值。

statsmodels 中 LME 的实现主要基于组，这意味着不同组中响应的随机效应必须是独立实现的。在我们的混合模型实现中，有两种类型的随机效应：（i）随机系数（可能是向量），这些系数具有未知的协方差矩阵，以及（ii）随机系数，这些系数是从共同单变量分布中独立抽取的样本。对于 (i) 和 (ii)，随机效应通过它们与组特定设计矩阵的矩阵/向量积影响组的条件均值。

一个简单的随机系数示例，如 (i) 中所示，是

\[Y_{ij} = \beta_0 + \beta_1X_{ij} + \gamma_{0i} + \gamma_{1i}X_{ij} + \epsilon_{ij}\]

这里，\(Y_{ij}\) 是受试者 \(i\) 的第 \(j^\rm{th}\) 个测量响应，而 \(X_{ij}\) 是此响应的协变量。所有受试者共用“固定效应参数” \(\beta_0\) 和 \(\beta_1\)，而误差 \(\epsilon_{ij}\) 与其他所有内容无关，并且独立同分布（均值为零）。“随机效应参数” \(\gamma_{0i}\) 和 \(\gamma_{1i}\) 遵循均值为零的二元分布，由三个参数描述：\({\rm var}(\gamma_{0i})\)、\({\rm var}(\gamma_{1i})\) 和 \({\rm cov}(\gamma_{0i}, \gamma_{1i})\)。还存在一个参数用于 \({\rm var}(\epsilon_{ij})\)。

一个简单的方差成分示例，如 (ii) 中所示，是

\[Y_{ijk} = \beta_0 + \eta_{1i} + \eta_{2j} + \epsilon_{ijk}\]

这里，\(Y_{ijk}\) 是在条件 \(i, j\) 下的第 \(k^\rm{th}\) 个测量响应。唯一的“均值结构参数”是 \(\beta_0\)。\(\eta_{1i}\) 独立同分布，均值为零，方差为 \(\tau_1^2\)，而 \(\eta_{2j}\) 独立同分布，均值为零，方差为 \(\tau_2^2\)。

statsmodels MixedLM 处理大多数非交叉随机效应模型，以及一些交叉模型。要在模型中包含交叉随机效应，需要将整个数据集视为一个组。然后可以使用模型的方差成分参数定义具有各种交叉和非交叉随机效应组合的模型。

statsmodels LME 框架目前支持通过 Wald 检验和系数置信区间、概要似然分析、似然比检验和 AIC 进行估计后推断。

示例

In [1]: import statsmodels.api as sm

In [2]: import statsmodels.formula.api as smf

In [3]: data = sm.datasets.get_rdataset("dietox", "geepack").data

In [4]: md = smf.mixedlm("Weight ~ Time", data, groups=data["Pig"])

In [5]: mdf = md.fit()

In [6]: print(mdf.summary())
         Mixed Linear Model Regression Results
========================================================
Model:            MixedLM Dependent Variable: Weight    
No. Observations: 861     Method:             REML      
No. Groups:       72      Scale:              11.3669   
Min. group size:  11      Log-Likelihood:     -2404.7753
Max. group size:  12      Converged:          Yes       
Mean group size:  12.0                                  
--------------------------------------------------------
             Coef.  Std.Err.    z    P>|z| [0.025 0.975]
--------------------------------------------------------
Intercept    15.724    0.788  19.952 0.000 14.179 17.268
Time          6.943    0.033 207.939 0.000  6.877  7.008
Group Var    40.394    2.149                            
========================================================

可以在此处找到详细示例

Wiki 上有一些笔记本示例：MixedLM 的 Wiki 笔记本

技术文档

数据被划分为不相交的组。组 \(i\) 的概率模型为

\[Y = X\beta + Z\gamma + Q_1\eta_1 + \cdots + Q_k\eta_k + \epsilon\]

其中

  • \(n_i\) 是组 \(i\) 中的观测值数量

  • \(Y\) 是一个 \(n_i\) 维响应向量

  • \(X\) 是一个 \(n_i * k_{fe}\) 维固定效应系数矩阵

  • \(\beta\) 是一个 \(k_{fe}\) 维固定效应斜率向量

  • \(Z\) 是一个 \(n_i * k_{re}\) 维随机效应系数矩阵

  • \(\gamma\) 是一个 \(k_{re}\) 维随机向量，均值为 0，协方差矩阵为 \(\Psi\)；请注意，每个组都会获得 gamma 的独立实现。

  • \(Q_j\) 是一个 \(n_i \times q_j\) 维设计矩阵，用于第 \(j^\rm{th}\) 个方差成分。

  • \(\eta_j\) 是一个 \(q_j\) 维随机向量，包含方差为 \(\tau_j^2\) 的独立同分布值。

  • \(\epsilon\) 是一个 \(n_i\) 维向量，包含均值为 0，方差为 \(\sigma^2\) 的独立同分布正态误差；\(\epsilon\) 值在组内和组间都是独立的

\(Y, X, \{Q_j\}\) 和 \(Z\) 必须完全被观测到。\(\beta\)、\(\Psi\) 和 \(\sigma^2\) 使用 ML 或 REML 估计进行估计，而 \(\gamma\)、\(\{\eta_j\}\) 和 \(\epsilon\) 是随机的，因此定义了概率模型。

边缘均值结构为 \(E[Y|X,Z] = X*\beta\)。如果仅关注边缘均值结构，GEE 是混合模型的良好替代方法。

符号

  • \(cov_{re}\) 是随机效应协方差矩阵（上面称为 \(\Psi\)），而 \(scale\) 是（标量）误差方差。每个方差成分也只有一个估计的方差参数 \(\tau_j^2\)。对于一个组，给定 exog 的 endog 的边缘协方差矩阵为 \(scale*I + Z * cov_{re} * Z\)，其中 \(Z\) 是一个组中随机效应的设计矩阵。

参考文献

实现细节的主要参考资料是

  • MJ Lindstrom, DM Bates (1988). 重复测量数据的线性混合效应模型的牛顿-拉夫森和 EM 算法。美国统计协会杂志。第 83 卷，第 404 期，第 1014-1022 页。

另请参阅这篇较新的文档

所有似然、梯度和 Hessian 计算都严格遵循 Lindstrom 和 Bates 的方法。

以下两篇文档更多地是从用户的角度编写的

模块参考

模型类是

MixedLM(endog, exog, groups[, exog_re, ...])

线性混合效应模型

结果类是

MixedLMResults(model, params, cov_params)

用于包含线性混合效应模型拟合结果的类。
最后更新时间：2024 年 10 月 3 日