【GLMM】一般化線形混合モデルについて解説|R
こんにちは、青の統計学です。
GLMMを学ぶ前には、固定効果とランダム効果(変量モデル)、そしてGLM(一般化線形モデル)を理解しておく必要があります。
まだ理解が足りてない方には、まず先に以下のコンテンツをご覧ください。
-GLM-
-固定効果とランダム効果-
【例題つき】固定効果推定と固定効果モデルについて解説|ランダム効果も添えて
一般化線形混合モデル(generalized liner mixed model)
GLMのおさらいをしましょう。
ある分布に従う確率変数Yがあるとして、その分布の平均をμとします。
\(μ\)を変形した形\(f(\mu)\)にすると、目的変数を線形で表すことができます。
これを一般化線形モデル(GLM)と呼びます。
もっと具体的にいうとGLMはリンク関数を用いて、応答変数と線形予測子の間の非線形関係をモデリングします。
これにより、正規分布以外の分布（二項分布、ポアソン分布など）に従う応答変数も扱えるようになります。
ポアソン回帰：目的変数がポアソン分布に従うとして、ポアソン分布の平均λを、\(f(λ)\)として対数の形で表したものをリンク関数として表すことで、右辺を線形として扱う。
ロジスティック回帰：目的変数が二項分布/ベルヌーイ分布に従うとき、\(f(\mu)\)にロジットリンク関数(逆数をとってシグモイド関数にする)であるものを扱う。
$$log\frac{q_i}{1-q_i}=z_i$$
これは、ロジット変換と呼ばれる変換で以下の確率\(q_i\)の対数オッズを取ります。
これが線形になるので解釈がしやすいです。
$$q_i=\frac{1}{1-exp(-z_i)}=logistic(z_i)$$
元々は、確率\(q\)は上のようなロジスティック関数ですね。
さてGLMには、以下のような特徴がありました。
1：確率分布
目的変数が従う確率分布です。ポアソン分布や二項分布などです。
2：リンク関数
下の数式でいうと右辺のことです。
3：線形予測子
下の数式でいうと、左辺のことです。
$$logλ_i=β_i+β_{2}x_i$$
これに加えて、GLMMには4つ目の特徴として、「混合効果」が付け加わります。
混合効果は、固定効果とランダム効果が一緒に同じ統計モデルに加わることです。
混合効果(mixed effect)
固定効果とランダム効果については、以前のコンテンツで解説しました。
今回は、簡単に説明します。
固定効果-fixed effect-
宗教観や住む場所の文化など、定量化しづらいand観測しづらい事象でも、結婚するかどうかなどの目的変数に影響があるものはあります。
特徴として、時間を通じて常に変わらない値です。
どの政党を支持しているか、どの宗教か、どこに住んでいるかなど、ある程度長期に亘り変わらないものを固定効果として採用しています。
$$Y_{it}=\alpha_i+\beta_1X_{1it}+\beta_2X_{2it}+..+\beta_kX_{kit}+u_{it}$$
固定効果は、\(\alpha_i\)です。
時間によって不変なので、\(t(time)\)の添字をつけていません。
ランダム効果-random effect-
ランダム効果は、「個体値」などのサンプル個体ごとのばらつきを表す効果です。
$$Y_{it}=\beta_1X_{1it}+\beta_2X_{2it}+..+\beta_kX_{kit}+\alpha_i+u_{it}$$
そのまま誤差項として扱うことができます。
説明変数とは、相関しないという強めの前提が置かれているため、固定効果推定量よりも信憑性が低いと言われています。
*固定効果とランダム効果の違い
固定効果は、説明変数\(X\)と相関します。それに対して、ランダム効果は説明変数\(X\)と相関しません。
以上の理由で、ランダム効果は統計モデルの誤差項として扱えて、固定効果は誤差項としては扱えないです。
(誤差項\(u_i\)は、説明変数\(X_i\)と相関しないという前提があります)
補足｜説明できる事象の範囲
GLMMを学べば、その先は「ベイズ統計学」と呼ばれる学問に行きます。推測統計学の終点です。
推測統計学の中でもトップクラスの守備範囲を誇ります。
これまでのGLMや線形モデルも含めてどのくらい複雑な現象をカバーできているのかを視覚的にみてみましょう。
まずは最小二乗法を使う線形回帰です。基本的に目的変数(y)が正規分布であることが前提にされています。
- データ\(X\)と\(Y\)が独立かつ同一な分布に従う。
- \(X\)と\(Y\)は線形関係がある。
- 説明変数を条件とすると、誤差項の期待値は0になる。\(E[u_i]=0\)
- \(X_i\)と\(u_i\)の4次モーメントが有限である。
例として、単回帰モデルの最小二乗推定のための仮定を載せました。
最小二乗法では扱えない、目的変数が正規分布以外の世界も扱いたい時にGLMを使います。
主な推定計算方法としては、最尤推定法が挙げられます。忘れてしまった方は、以下のコンテンツをご覧ください。
- 線形予測子
- リンク関数
- 確率分布
この3つの特徴を持つのがGLM(一般化線形モデル)でした。ポアソン分布や二項分布などのカウントデータも扱えるため、汎用性はかなり上がりました。
これが今回扱うGLMM(一般化混合線形モデル)の立ち位置です。
先ほども解説した通り、GLMに固定効果とランダム効果を加えて、観測できない効果をモデルに入れようと試みています。
ロジスティック回帰やポアソン回帰といった GLM では全サンプルの均質性を仮定しているのですが、基本的に現実世界のデータは過分散です。
CODE|R
具体例はそのまま、【二項分布】ロジスティック回帰について(実例つき)の例を扱います。
2055年に亡くなった方の健康寿命を2040年時点から遡ってみました。 中でも、喫煙経験があった25人と25人の計50人を調査し、肺の大きさも測ります。 N:2040年から2055年までの年数。当然全員15年。 y:2040年から測った健康寿命。15以下。 x:腎臓の長さ s:喫煙因子/喫煙してるかしてないかの二択。 N:15 y:10 x:10.43 s:N 例：2050年までは健康に生きれたが2055年に亡くなった。肺の横の長さは10.43cm。喫煙経験はない人。
この時に追加するランダム効果は、被検体の個体差\(r_i\)です。モデルは以下のようになります。
この個体差\(r_i\)は、「正規分布」に従い、かつ個体間の平均が0であると仮定しています。
$$ln\frac{y_i}{1-y_i}=beta_0+beta_1X_i+beta_2s_i+r_i$$
日本語で表すとこんな形です。
【注意】今回は、時間tを考えていないので、上の「固定効果モデル」の話とは異なり、説明変数Xなども含めた部分が固定効果と呼びます。
このため、固定効果推定は行わないことになりますが、固定効果推定について知りたい方は、以下のコンテンツをご覧ください。
【例題つき】固定効果推定と固定効果モデルについて解説|ランダム効果も添えて
第一項から第三項までが、固定効果と呼ばれます。
では、Rでモデルを作ってみましょう。
update.packages()
install.packages("glmmML")
library(glmm)まずはパッケージのインストールを行います。これでGLMMが使えることになりました。
data<-read.csv("health.csv")
data$id
#データの読み込み
今回は、50人にidを振りました。これが、今回のランダム効果の元となります。
glmmML(cbind(y,N-y)~x+s,data=data,family=binomial(),cluster = id)
#GLMMの利用
clusterというのは、「個体ごとに異なる独立なパラメーター」のことです。
今回でいうランダム効果ですね。
global parameter と local parameterという考えがあり、個体ごとに値が変わる\(r\)はlocal parameterと呼び、最尤推定の過程でも一意の値を推定できません。
global parameterというのは、固定効果や各回帰係数で、最尤推定が可能です。
AICを見てみましょう。221.4から67.75になりました。
ランダム効果を加えただけで、以前の【二項分布】ロジスティック回帰について(実例つき)と比べてかなり予測の精度が上がっていることがわかります。
AICに関しては、【良いモデルとは】AIC(赤池情報量基準)について。をご覧ください。
＊ここまで、固定効果とランダム効果の違いを説明してきましたが、定義や認識の仕方が違うのが実情です。
ただ、「GLMMは、GLMに加えて誤差項\(u_i\)を足すことにより、個体差などを表すことができる」ことは覚えていただければと思います。
