GLMM の基礎

GLMM タブで使われている統計理論の背景です。操作方法は GLMM のページを参照してください。

階層データと独立性の問題

統計モデルの多くは、観測が互いに独立であることを仮定しています。しかし現実のデータには、観測がグループに属する階層構造を持つものが多くあります。

複数の学校に通う生徒のテスト成績（生徒はそれぞれの学校に所属）
複数の病院で治療を受けた患者の回復日数（患者はそれぞれの病院に所属）
同じ被験者に対する反復測定（測定は被験者に所属）

同じグループに属する観測は、グループ固有の要因（学校の教育方針、病院の設備、被験者の体質）を共有するため、互いに似た値を取る傾向があります。この相関を無視して GLM で分析すると、モデルの誤特定（misspecification）になります。データの生成過程にグループ構造があるのに、モデルがそれを表現していない状態です。

この誤特定は複数の問題を引き起こします:

標準誤差の過小推定: グループ内の相関により、実効的なサンプルサイズは見かけの $n$ より小さくなります。独立性を仮定した標準誤差はこれを反映せず、信頼区間が狭くなりすぎます
係数の偏りの可能性: グループレベルの交絡要因がある場合、固定効果の推定が偏ります。たとえば学校の教育方針が勉強時間と成績の両方に影響する場合、学校差を無視した回帰では勉強時間の効果を正しく推定できません。これは省略変数バイアスの問題であり、ランダム切片はグループレベルの交絡を部分的に吸収しますが、完全には解決しません
予測の不確実性の過小評価: グループレベルの変動を無視した予測区間は、実際の変動を過小に見積もります
分散構造の情報の喪失: グループ間とグループ内のばらつきを区別できないため、データの構造を正しく把握できません

混合モデル（mixed model）は、グループ内の相関を変量効果として明示的にモデル化することで、標準誤差の過小推定と分散構造の喪失を直接的に解決します。

固定効果と変量効果

混合モデルの「混合」は、固定効果と変量効果の両方を含むことを指します。

固定効果（fixed effects） は、研究者が推定したい体系的な関係です。たとえば「勉強時間が1時間増えるとテスト成績が何点上がるか」という効果は固定効果です。

変量効果（random effects） は、直接観測できないグループレベルのばらつきをモデル内で表現するための装置です。各学校の平均成績は教育方針や生徒層によって異なります。このばらつきは確率的なメカニズムで生まれたものではなく、具体的な要因に起因しますが、それらの要因をすべて説明変数として測定・投入するのは現実的ではありません。変量効果はこの「観測できないグループ差」を確率分布で近似的に表現します（Hodges, 2014）。

$u_j \sim N(0, \sigma_u^2)$ という仮定は「学校差が正規分布から生成されている」という主張ではありません。グループ間のばらつきの大きさ $\sigma_u^2$ を推定し、それに基づいて縮小推定（後述）を行うための、モデリング上の仮定です。

この解釈では、固定効果/変量効果の区別は「グループがランダムに抽出されたかどうか」ではなく、「グループ間の差をどうモデル化するか」の設計判断です（Gelman, 2005）。グループ差の全体的なばらつきに関心がある場合は変量効果が適切で、特定のグループの効果を直接推定したい場合は固定効果として扱います。

ランダム切片モデル

最も基本的な混合モデルがランダム切片モデルです。各グループに固有の切片（ベースライン）を与えます:

g(\mu_i) = x_i'\beta + u_{j[i]}, \quad u_j \sim N(0, \sigma_u^2)

$g$ はリンク関数、 $x_i'\beta$ は固定効果の線形予測子、 $u_{j[i]}$ は観測 $i$ が属するグループ $j$ のランダム切片です。

直感的には、全グループに共通する回帰直線（ $x_i'\beta$ ）を引き、各グループの直線はそこから $u_j$ だけ上下にずれる、という構造です。正規分布 $N(0, \sigma_u^2)$ の仮定は、このずれの大きさを $\sigma_u^2$ という1つのパラメータで要約するための道具です。各 $u_j$ がどんな値をとるかは、データから推測します（推定と予測を参照）。

Gaussian ファミリーの場合、これは線形混合モデル（LMM）です:

y_i = x_i'\beta + u_{j[i]} + \varepsilon_i, \quad u_j \sim N(0, \sigma_u^2), \quad \varepsilon_i \sim N(0, \sigma_e^2)

$\sigma_u^2$ はグループ間のばらつき、 $\sigma_e^2$ はグループ内（個人レベル）のばらつきです。

パラメータ推定

混合モデルのパラメータ推定は、通常の回帰より複雑です。固定効果 $\beta$ と分散成分（ $\sigma_u^2$ , $\sigma_e^2$ ）を同時に推定する必要があるためです。

REML（制限付き最尤法）

最尤法（ML）は分散成分を過小推定する傾向があります。これは $\beta$ の推定に自由度を消費することを尤度関数が考慮しないためです（標本分散の分母に $n$ ではなく $n-1$ を使うのと同じ理由）。

REML は $\beta$ を射影で除去した残差の尤度を最大化することで、この偏りを補正します。

AIC/BIC の制限

REML の射影は固定効果の構成に依存するため、REML ベースの AIC/BIC は固定効果が異なるモデル間で比較できません。固定効果の構成が異なるモデルを比較するには、ML で再推定した AIC/BIC を使う必要があります。

Laplace 近似と PIRLS

Gaussian ファミリーでは変量効果の積分が解析的に実行できるため REML で十分ですが、Non-Gaussian ファミリー（Binomial, Poisson, Gamma など）では、変量効果を積分除去する操作:

L(\beta, \sigma_u^2) = \int \prod_i f(y_i \mid \beta, u) \cdot \phi(u; 0, \sigma_u^2) \, du

に閉じた形の解がありません。 $\phi$ は正規分布の密度関数です。

Laplace 近似はこの積分を、被積分関数の最大点周辺での2次近似で置き換えます。MIDAS は Bates et al. (2015) が記述した lme4 のアルゴリズムに基づいています。Laplace 近似はグループあたりの観測数が少ない場合や、Binomial でイベントが稀な場合に精度が低下することがあります。

具体的には:

外側ループ: 相対共分散パラメータ $\theta$ を Nelder-Mead 法で最適化
内側ループ: PIRLS（Penalized IRLS）で $(\beta, u)$ を同時推定

PIRLS は GLM の IRLS にランダム効果のペナルティ項を加えたものです。ペナルティ付き逸脱度（Bates et al., 2015 の表記に従い、 $\theta$ は相対共分散パラメータ）:

D_p = \sum_i d(y_i, \mu_i) + \sum_j \frac{u_j^2}{\theta^2}

第1項はデータへの適合度、第2項は $u_j$ が正規分布から大きく外れることへのペナルティです。このペナルティにより、データの少ないグループの推定値は全体平均に引き寄せられます（後述の縮小推定を参照）。

推定と予測

固定効果 $\beta$ は未知の定数なので「推定（estimation）」の対象です。一方、変量効果 $u_j$ はモデル内で確率変数として扱われるため、形式的には「予測（prediction）」の対象になります。

$u_j$ が「予測」と呼ばれるのは、モデルの数学的な構造に由来します。 $u_j$ に確率分布を仮定しているため、データだけでなくその分布からの情報も組み合わせて $u_j$ の値を推測します。 $\beta$ の推定ではデータの有限性だけが不確実性の原因ですが、 $u_j$ の予測ではモデルが仮定した分布（ $N(0, \sigma_u^2)$ ）も情報源として働きます。この2つの情報源の組み合わせが、後述する縮小推定を生み出します。

「予測」という用語は、 $u_j$ が本当にランダムに生成されたことを前提としているわけではありません。グループ差を確率分布で近似するというモデリング上の判断に伴う、形式的な区別です。

変量効果の予測に BLUP（Best Linear Unbiased Predictor）がよく使われます。BLUP は線形予測量のうち $E[\hat u_j - u_j] = 0$ （予測の不偏性）を満たし、平均二乗予測誤差が最小になるもので、Henderson (1950) の混合モデル方程式から導かれます。

Gaussian の場合、 $u_j$ の BLUP は縮小推定（shrinkage estimator）の形をとります:

\hat{u}_j = \frac{n_j}{n_j + \sigma_e^2 / \sigma_u^2} \times \bar{r}_j

$n_j$ はグループ $j$ のサイズ、 $\bar{r}_j$ はグループ $j$ の平均残差（固定効果で説明できない部分）です。

係数 $n_j / (n_j + \sigma_e^2/\sigma_u^2)$ は0から1の値をとり、グループサイズが大きいほど1に近づきます。つまり:

大きなグループ: データが十分にあるので、そのグループ固有の推定値をほぼそのまま使う
小さなグループ: データが少ないので、全体平均（ゼロ）に向かって引き寄せる

これは「情報の借用（borrowing strength）」とも呼ばれます。データの少ないグループは他のグループの情報を借りて推定を安定させます。グループ固有の推定値をそのまま使うと分散が大きくなりますが、全体平均だけを使うとグループの特性を無視してしまいます。BLUP はこのバイアスとバリアンスのトレードオフを最適にバランスさせます。

ICC（級内相関係数）

ICC（Intraclass Correlation Coefficient）は、データの全分散のうちグループ間の違いが占める割合です:

\text{ICC} = \frac{\sigma_u^2}{\sigma_u^2 + \sigma_e^2}

ICC = 0 ならグループ間に差がなく、すべてのばらつきは個人レベルです。ICC = 1 ならグループ内の個人差がなく、すべてのばらつきはグループ間の差です。

ICC が小さければ、グループ構造を無視して GLM で分析してもほぼ同じ結果が得られます。ただし影響の大きさは ICC だけでなくグループサイズにも依存します。設計効果 $\text{DEFF} = 1 + (\bar n - 1) \times \text{ICC}$ （ $\bar n$ はグループの平均サイズ）が標準誤差の膨張倍率の目安になります。ICC が大きい場合は混合モデルを使う必要があります。

Non-Gaussian の場合

Non-Gaussian ファミリーでは残差分散 $\sigma_e^2$ が直接定義されないため、潜在変数の分散を代用します。たとえば Binomial + Logit リンクでは、潜在連続変数がロジスティック分布に従うと考え、その分散 $\pi^2/3$ を使います。

ファミリー + リンク	潜在尺度の残差分散
Binomial + Logit	$\pi^2/3 \approx 3.29$
Poisson + Log	$1$ （対数正規-Poisson 混合モデルにおける近似）
Gamma + Log	$\phi$ （推定された分散パラメータ）

MIDAS の実装

MIDAS の GLMM は分布ファミリーに応じて異なるアルゴリズムを使います。

Gaussian（LMM）: Profile REML。分散成分の比 $\theta^2 = \sigma_u^2/\sigma_e^2$ を Nelder-Mead 法で数値最適化し、各 $\theta$ に対して $\beta$ と $\sigma_e^2$ を閉じた形で求めます。

Non-Gaussian（Binomial, Poisson, Gamma 等）: Bates et al. (2015) に基づく Laplace 近似。外側ループで相対共分散パラメータ $\theta$ を Nelder-Mead 法で最適化し、内側ループで PIRLS により $(\beta, u)$ を同時推定します。分散パラメータ $\phi$ は Poisson/Binomial では1に固定、Gamma では profiled deviance（ $D_p/n$ ）から推定します。

固定効果の初期値には GLM（ランダム効果なし）の推定結果を使います。説明変数は内部的にスケーリング（標準化）され、推定後に元のスケールに逆変換されます。

ランダム傾きと交差ランダム効果

ランダム切片モデルはグループごとのベースラインの違いだけをモデル化します。説明変数の効果もグループごとに異なる場合は、ランダム傾き（random slope）モデルが必要です:

g(\mu_i) = x_i'\beta + u_{0j} + u_{1j} x_{1i}

$u_{0j}$ はランダム切片、 $u_{1j}$ はランダム傾きで、両者は多変量正規分布に従います。

また、観測が複数のグルーピング変数に属する場合（たとえば生徒が学校と地域の両方に属する）、交差ランダム効果（crossed random effects）が使えます。

MIDAS は現在ランダム切片モデルのみ対応しています。

参考文献

Henderson, C. R. (1950). Estimation of genetic parameters (Abstract). Annals of Mathematical Statistics, 21, 309-310.
Gelman, A. (2005). Analysis of variance—why it is more important than ever. The Annals of Statistics, 33(1), 1-53.
Hodges, J. S. (2014). Richly Parameterized Linear Models: Additive, Time Series, and Spatial Models Using Random Effects. Chapman & Hall/CRC.
Bates, D., Mächler, M., Bolker, B., & Walker, S. (2015). Fitting linear mixed-effects models using lme4. Journal of Statistical Software, 67(1), 1-48.