線性混合效應模型

線性混合效應模型用於分析涉及相依資料的迴歸。當處理縱向研究或其他對每個受試者進行多次觀察的研究設計時,就會出現這種資料。一些特定的線性混合效應模型包括:

  • 隨機截距模型:其中群組中的所有反應都以特定於該群組的值進行加法偏移。

  • 隨機斜率模型:其中群組中的反應遵循一個(條件)平均軌跡,該軌跡在觀察到的共變數中是線性的,斜率(以及可能的截距)會因群組而異。

  • 變異數組件模型:其中一個或多個類別共變數的水平與分佈中的抽取相關。這些隨機項會根據其共變數值,以加法方式決定每個觀察的條件平均值。

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\) 是平均值為 0 且共變異數矩陣為 \(\Psi\)\(k_{re}\) 維度隨機向量;請注意,每個群組都會獲得其自己的獨立實現的 gamma。

  • \(Q_j\) 是第 \(j^\rm{th}\) 個變異數組件的 \(n_i \times q_j\) 維度設計矩陣。

  • \(\eta_j\) 是包含獨立且具有相同分佈的值的 \(q_j\) 維度隨機向量,變異數為 \(\tau_j^2\)

  • \(\epsilon\) 是平均值為 0 且變異數為 \(\sigma^2\) 的 i.i.d 常態誤差的 \(n_i\) 維度向量;\(\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)。Newton Raphson and EM algorithms for linear mixed effects models for repeated measures data。美國統計協會雜誌。第 83 卷,第 404 期,第 1014-1022 頁。

另請參閱此較新的文件

所有概似、梯度和 Hessian 計算都密切遵循 Lindstrom 和 Bates。

以下兩個文件從使用者的角度撰寫

模組參考

模型類別是

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

線性混合效應模型

結果類別是

MixedLMResults(model, params, cov_params)

包含擬合線性混合效應模型結果的類別。


上次更新:2024 年 10 月 03 日