マルチレベル分析の基礎〜予測だけでは不満なあなたへ〜

最近はAIを使った予測が流行っていますが、予測よりも現象を解釈したい場合もあるかと思います。

行った施策に対しての効果を分析したり、事象が起こった後に要因を解釈したり...

そんな時に使うのがマルチレベル分析です。現実の仮定を反映しやすく、様々な仮説について解析をすることができます。

今回はそんなマルチレベル分析を導入する必要性や、解析、結果の解釈方法などをまとめていきます!

Rでの実行方法は補足として記事中に記載してあるので、参考にしてみてください!

マルチレベル分析とは?

マルチレベル分析とは、複数の階層性をもつデータを解析する時に用いられる分析手法です。例えばテストの点数と勉強時間との関連を調べる時に対象となる学校を選び、そこから数人の生徒を選んだ場合や、残業時間と従業員の満足度を調べる際に、企業を選び、従業員を選んだ場合、季節ごとに商品の売り上げを調査した場合など、二段階以上のレベルからサンプルを抽出したデータを階層性のあるデータと呼びます。

f:id:h-wadsworth02:20190120161444j:plain

www.medi-08-data-06.work

なぜマルチレベル分析が必要か?

ここでは、テストの点数と勉強時間の関連性を調べることにしましょう。

以下の二つの疑問を考えてみます。

  • 学校間で平均点数に違いがあるか
  • 個人レベルで勉強時間と点数に関連はあるのか

とりあえず、二つ目の疑問から考えてみます。

全部で7校の学校から生徒を数人ずつ選んだと仮定して、テストの点数と勉強時間から回帰直線を求め、散布図に重ねてみます。

f:id:h-wadsworth02:20190211180901j:plain

X軸が勉強時間、Y軸が点数を表しています。ここからは勉強時間が1時間長いと点数が19点高くなることが分かります。

(上記のデータセットは以下のように作成しました)

library(dplyr)
#データセットの作成
set.seed(123)
mu <- list(20,20,30,30,40,50,50)
N_st <- list(10,5,40,25,5,15,10)
ID <- map(N_st,function(n)seq(1:n)) %>% 
  flatten_dbl()
X <- map2(N_st,mu,rnorm,sd=10) %>% 
  flatten_dbl()
sctype <- sample(c(0,1),replace = T,7)
school <- rep(1:7, times=unlist(N_st))

a0 <- 1000
b0 <- 5
a <- a0 + rnorm(length(N_st), mean=0, sd=300) %>% 
  sort(decreasing = T)
b <- b0 + 10*sctype+rnorm(length(N_st), mean=0, sd=10)
sctype <- sctype[order(b)]
b <- b[order(b)]

data_frame(X=X,school = as.factor(school),
           a=a[school],b=b[school],sctype=sctype[school],ID=ID) %>% 
  mutate(Y=rnorm(length(X),a+b*X,100)) %>%
  select(X,ID,school,sctype,Y)-> st_df

Xは勉強時間、IDは学校内での生徒番号、schoolは学校識別ID、sctypeは公立か私立か(後から使います)、Yはテストの点数を表しています。

X ID school sctype Y
14.39524 1 1 0 1453.686
17.69823 2 1 0 1341.678
35.58708 3 1 0 1430.010
20.70508 4 1 0 1476.366
21.29288 5 1 0 1520.955
37.15065 6 1 0 1333.399

今度は、学校ごとに色付けをしてしてみましょう。

f:id:h-wadsworth02:20190211180738j:plain

どうでしょうか?どうやら学校ごとに勉強時間や点数に傾向があり、単純に個人レベルで勉強時間が長いほど点数が伸びるとは言えなさそうです。

さらに、仮に以下のような関係であれば、個人の勉強時間と点数に関連がなくても、学校全体で見ると関連があるように見えてしまいます。

f:id:h-wadsworth02:20190211181454j:plain

これがマルチレベル分析を導入するべき理由です。

マルチレベル分析を利用すると、個人レベルと学校レベルの関連を分離させてデータの解釈をすることができます。

マルチレベル分析の適応

さて、マルチレベル分析を適応すべきかどうかはどのように判断すれば良いでしょうか?もちろんデータが収集された背景を考えて決めても良いですが、何か指標が欲しいところですよね。

ところで、目的変数である点数Yのバラツキには、学校間でのバラツキと個人間でのバラツキが含まれています。

 Yのバラツキ= 学校間のバラツキ+個人間のバラツキ

マルチレベル分析は、バラツキを学校レベルと個人レベルに分離させて解析するというのが大まかな目的となります。

そこで級内相関というものを使います。級内相関とは同じ群内にある値がどれぐらい似ているかという指標になるのですが、ここでは学校間の点数のバラツキが学校内のバラツキに比べてどれぐらい大きいかを表します。

つまり、個人よりも学校による点数の違いが大きい場合は、マルチレベル分析を適応するべきだということになります。

この級内相関を求める際には、今回は傾きを0とした線形混合モデルの中のランダム化切片モデルというものを使います。難しそうな名前ですが、各学校の点数の平均値(切片)が全体平均を中心とした正規分布に従い、学校内の生徒の点数は学校の平均値を中心とした正規分布に従うことを前提とするモデルのことを指します。

傾きを0としているので、一元配置分散分析と似たようなイメージです。

下の図はX軸を各学校、Y軸を点数としたグラフです。

f:id:h-wadsworth02:20190211201201p:plain:w420

赤色の正規分布の分散が、学校間の分散、黒色が学校内の分散を表しています。

(ちなみに正規分布を仮定しないようなデータは線形混合モデルでは対応できないので、ベイズ推定による階層ベイズモデルを使います。)

www.medi-08-data-06.work

難しい方はざっくりとイメージだけでもつかめればOKです。

式で表すと

 Y_{st} = a_{sc} + \epsilon\\
a_{sc} = a0 + \epsilon_{a}

\epsilon \sim Normal(0\,, \sigma^{2})\\
 \epsilon_{a} \sim Normal (0 \,, \sigma_{sc}^{2})

Y_{st}\,:\,生徒の点数\\
a_{sc}\,:\,学校ごとの平均値\\
\sigma^{2}\,:\,誤差の分散(学校内分散)\\
a0\,:\,切片平均\\
\sigma_{sc}^{2}\,:\,a_{sc}の分散(学校間分散)

となります。

ここから級内相関は、

級内相関 = \dfrac{群間分散}{群間分散+郡内分散}\\
= \dfrac{\sigma_{sc}^{2}}{\sigma_{sc}^{2}+\sigma^{2}}

として求めることができます。

これは、群内のバラツキと比べて、群間のバラツキがどれぐらい大きいか、を表す指標になり、まさにマルチレベル分析を適応すべきかどうかを判断できる指標ですね。

これらのパラメータを求めてみると、


\sigma^{2}\,:\,28116(学校内分散)\\
a0\,:\,  1617.6\\
\sigma_{sc}^{2}\,:\,98742(学校間分散)

となってICCは98742/(98743+28116)=約78%となります。

---補足---

学校内分散や学校間分散の求め方についてご質問頂きましたので、コメントにて回答させて頂きました。

----------

明確な基準はないですが、回帰分析は全ての観測値が独立であることを前提としているので、級内相関が0.1(10%)でもマルチレベル分析を適応すべきではないかとされています。

この級内相関が意味するところは、点数変動の78%は学校間変動によって説明されるという意味になります。

ところで、学校間の違いによる影響が大きいということは、学校間で平均点に大きな違いがあるということです。

つまり、学校間で平均点数に違いがあり、それは切片平均(1617点)を中心に314点(98742の平方根)ぐらいのバラツキがあるということで、第一の疑問に対して答えることができました!

Rで実行するには{lme4}パッケージを使います。{lme4}の使い方についてはこちらも参考にしてみてください。

www.medi-08-data-06.work

library(lme4)

st_df %>% 
  lme4::lmer(data=.,Y~(1|school),REML = F) %>% 
  summary()

>Linear mixed model fit by maximum likelihood  ['lmerMod']
Formula: Y ~ (1 | school)
   Data: .

     AIC      BIC   logLik deviance df.resid 
  1471.5   1479.6   -732.7   1465.5      107 

Scaled residuals: 
     Min       1Q   Median       3Q      Max 
-2.77716 -0.57499 -0.02915  0.56146  3.01817 

Random effects:
 Groups   Name        Variance Std.Dev.
 school   (Intercept) 98742    314.2   
 Residual             28116    167.7   
Number of obs: 110, groups:  school, 7

Fixed effects:
            Estimate Std. Error t value
(Intercept)   1617.6      120.5   13.42

REML=Fはパラメーターの推定方法に制限つき最尤法(REsidual Maximum. Likelihood estimation)を使うか、通常の最尤法を使うかを表します。今回は通常の最尤法を使うためFALSEとしました。

結果の見方ですが、Random effectsのschool(intercept)が学校間分散\sigma_{sc}、Residualが全分散から学校間分散を引いた学校内分散\sigma、Fixed effectsの(intecept)が切片平均a0を表しています。

> #ICC
> 98742/(98743+28116)
[1] 0.7783602

試しにランダムなデータで分散分析モデルを試してみると

set.seed(123)
d <- data_frame(x=c(rnorm(70,50,10)),
           school = as.factor(rep(c(1:7),rep(10,7))))

d %>% 
  lmer(data=.,x~(1|school),REML = F) %>% 
  summary()

singular fit
Linear mixed model fit by maximum likelihood . t-tests use
  Satterthwaite's method [lmerModLmerTest]
Formula: x ~ (1 | school)
   Data: .

     AIC      BIC   logLik deviance df.resid 
   512.4    519.2   -253.2    506.4       67 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-2.2648 -0.6313 -0.0134  0.6535  2.3255 

Random effects:
 Groups   Name        Variance Std.Dev.
 school   (Intercept)  0.00    0.000   
 Residual             81.17    9.009   
Number of obs: 70, groups:  school, 7

Fixed effects:
            Estimate Std. Error     df t value Pr(>|t|)    
(Intercept)   50.738      1.077 70.000   47.12   <2e-16 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
convergence code: 0
singular fit

singular fitと出ます。切片の分散(学校間分散)がとても小さいため推定がうまくできていません。つまり、学校間の違いはほとんどないことを表しています。

ランダム化切片モデル

分析を進めていきます。先ほど学校間でテストの平均点に差があることが分かりましたので、仮説を少し変えて以下の疑問に答えていきましょう。

  • 勉強時間の影響を除いた上で、学校間で平均点数に違いがあるか
  • 個人レベルでも勉強時間と点数に関連はあるのか

まずは、勉強時間の影響を除いた上で学校間で平均点に差があるかどうかを明らかにしていきます。

どいうことかというと、テストの点数だけではなく、勉強時間にも学校間で差があることが考えられるため、仮に勉強時間が同じ生徒を比べた時にも、学校間で差があるかどうかを検証します。

イメージは以下のようになります。

f:id:h-wadsworth02:20190211182521j:plain

この切片分散が0でなければ勉強時間の影響を除いた後でも、学校ごとに点数に違いがあることになります。ここでは傾きが学校間で一緒であるという前提のもとに解析するため、共分散分析のイメージに近いですね。

式で書いてみると

 Y_{st} = a_{sc} +b0*X_{st} + \epsilon\\
a_{sc} = a0 + \epsilon_{a}

\epsilon \sim Normal(0\,, \sigma^{2})\\
 \epsilon_{a} \sim Normal (0 \,, \sigma_{sc}^{2})

Y_{st}\,:\,生徒の点数\\
X_{st}\,:\,生徒の勉強時間\\
a_{sc}\,:\,学校ごとの切片\\
b0\,:\, 傾き(勉強時間の点数への影響)\\
\sigma^{2}\,:\,誤差の分散(学校内分散)\\
a0\,:\,切片平均\\
\sigma_{sc}^{2}\,:\,a_{sc}の分散(学校間分散)

これを求めてみると、

b0\,:\, 16\\
\sigma^{2}\,:\,16318\\
a0\,:\,961\\
\sigma_{sc}^{2}\,:33510

となりました。

ここから、学校間分散が0ではないため勉強時間の影響を除いても学校間には183点(33510の平方根)ぐらいの差があると言えそうです。

さらに、級内相関は33510/(33510+16318)=0.67(67%)となり、少し小さくなりました。学校内分散の減少率よりも学校間分散の減少率が大きい、つまり勉強時間は学校レベルの違いに大きく影響したことが分かります。

st_df %>% 
  lme4::lmer(data=.,Y~X+(1|school),REML = F) %>% 
  summary() 


Linear mixed model fit by maximum likelihood [lmerMod]
Formula: Y ~ X + (1 | school)
   Data: .

     AIC      BIC   logLik deviance df.resid 
  1410.0   1420.8   -701.0   1402.0      106 

Scaled residuals: 
     Min       1Q   Median       3Q      Max 
-2.16891 -0.65933 -0.02279  0.63100  2.67048 

Random effects:
 Groups   Name        Variance Std.Dev.
 school   (Intercept) 33510    183.1   
 Residual             16318    127.7   
Number of obs: 110, groups:  school, 7

Fixed effects:
            Estimate Std. Error t value
(Intercept)  1171.88      85.43  13.717
X              12.74       1.36   9.365

Correlation of Fixed Effects:
  (Intr)
X -0.558


#説明変数を入れないモデル
Random effects:
 Groups   Name        Variance Std.Dev.
 school   (Intercept) 98742    314.2   
 Residual             28116    167.7   
Number of obs: 110, groups:  school, 7

ここでは、学校間の違いに興味があるためRandom effectsのみ記載しました。

説明変数の中心化

さて、次の疑問である個人レベルでも勉強時間と点数に関連はあるのかに答えるために、中心化の概念を説明しておきます。中心化とは説明変数のバラツキを個人レベルと集団レベルに分離して解析するための手法です。

先ほどのモデルの傾きb0を確認してみると16となっていて、これは勉強時間が1時間長いと点数が16点高くなることを示しています。

しかし、勉強時間には学校レベルと個人レベルの影響が混在しているため、個人の関係を表しているのか、そもそも勉強時間が平均的に長い学校の点数が高いのか分かりません。

そのため、各個人の勉強時間から学校ごとの平均勉強時間を引くことで、学校による勉強時間の違いを除くことができます。 これを集団平均中心化(Centering Within Cluster:CWC)と呼びます。 以下はCWCを行う前とCWCを行った後の勉強と点数の関係を表したグラフです。

f:id:h-wadsworth02:20190211182922j:plain

f:id:h-wadsworth02:20190211182936j:plain

下のグラフでは、勉強時間の平均値が0になっており、学校による勉強時間の違いが排除されていることがわかります。

このCWCを行った後の勉強時間を説明変数とすることで、純粋な個人レベルの勉強時間と点数の関係を推定することができます。

また、全体平均中心化(centerring at the grand mean:CGM)という手法もあります。これは、各個人の勉強時間から全体の勉強時間平均を引く手法です。

さて、説明変数にCWC(集団平均中心化)を行った勉強時間をモデルに加えても傾きが有意であれば、個人レベルで勉強時間と点数に関係があると言えそうです。

式で表すと

 Y_{st} = a_{sc} +b0*(X_{st}-\bar{X_{sc}}) + \epsilon\\
a_{sc} = a0 + \epsilon_{a}

\epsilon \sim Normal(0\,, \sigma^{2})\\
 \epsilon_{a} \sim Normal (0 \,, \sigma_{sc}^{2})

Y_{st}\,:\,生徒の点数\\
X_{st}\,:\,生徒の勉強時間\\
\bar{X_{sc}}\,:\,学校ごとの平均勉強時間\\
a_{sc}\,:\,学校ごとの切片\\
b0\,:\, 傾き(勉強時間の点数への影響)\\
\sigma^{2}\,:\,誤差の分散(学校内分散)\\
a0\,:\,切片平均\\
\sigma_{sc}^{2}\,:\,a_{sc}の分散(学校間分散)

となります。

これを求めてみると、


b0\,:\, 12 (SE=1.39) \\
\sigma^{2}\,:\, 16278\\
a0\,:\,1617\\
\sigma_{sc}^{2}\,:\,99957

となりました。ここから傾き12、標準誤差1.39なので、個人レベルでも勉強時間と点数に関係があることが分かりました!

説明変数を中心化したため、切片平均はそれぞれの学校で勉強時間が学校平均の人の点数の平均を表し、学校間分散はそのバラツキ具合となります。(ややこしい....)

#集団平均中心化のXを作る
st_df %>% 
  group_by(school) %>% 
  mutate(X_sc=mean(X),X_cwc=X-X_sc) -> st_df2

st_df2 %>% 
  lme4::lmer(data=.,Y~X_cwc+(1|school),REML = F)

Linear mixed model fit by maximum likelihood  ['lmerMod']
Formula: Y ~ X_cwc + (1 | school)
   Data: .

     AIC      BIC   logLik deviance df.resid 
  1417.2   1428.0   -704.6   1409.2      106 

Scaled residuals: 
     Min       1Q   Median       3Q      Max 
-2.26259 -0.65575 -0.01902  0.65463  2.67923 

Random effects:
 Groups   Name        Variance Std.Dev.
 school   (Intercept) 99957    316.2   
 Residual             16278    127.6   
Number of obs: 110, groups:  school, 7

Fixed effects:
            Estimate Std. Error t value
(Intercept) 1617.087    120.504  13.419
X_cwc         12.063      1.394   8.655

Correlation of Fixed Effects:
      (Intr)
X_cwc 0.000 

Fixed effectsの(intercept)とX_CWCが切片と傾きになります。

---2月14日追記---

学校レベルの影響を追加

さて、中心化を行うと学校レベルの影響が排除されることが分かりました。するとこのモデルでは、学校間での比較をすることができなくなります。

そこで、マルチレベル分析において集団平均中心化を用いた際は、その集団平均を集団レベルの変数に加えることが推奨されています。

つまり、平均的に勉強時間が長い学校は点数が高いのかどうか、という学校レベルでの議論できるようになります。

また、集団レベルの変数を投入する際には全体平均中心化(CGM)をすることで、今後の解析が行いやすくなります。CGMを行う詳しい理由については、下記の書籍を参考にしてみてください!

Rで学ぶ マルチレベルモデル[入門編]: 基本モデルの考え方と分析

Rで学ぶ マルチレベルモデル[入門編]: 基本モデルの考え方と分析

以上を式で表すと

 Y_{st} = a_{sc} +b0*(X_{st}-\bar{X_{sc}}) + \epsilon \\
a_{sc} = a0 + ab_{1} * (\bar{X_{sc}} - \bar{X_{st}})+\epsilon_{a}

\epsilon \sim Normal(0\,, \sigma^{2})\\
 \epsilon_{a} \sim Normal (0 \,, \sigma_{sc}^{2})

Y_{st}\,:\,生徒の点数\\
X_{st}\,:\,生徒の勉強時間\\
\bar{X_{sc}}\,:\,学校ごとの平均勉強時間\\
\bar{X}_{st}\,:\,勉強時間の全体平均\\
a_{sc}\,:\,学校ごとの切片\\
b0\,:\, 傾き(勉強時間の点数への影響)\\
\sigma^{2}\,:\,誤差の分散(学校内分散)\\
a0\,:\,切片平均\\
ab_{1}\,:\,学校レベルでの勉強時間の点数への影響\\
\sigma_{sc}^{2}\,:\,a_{sc}の分散(学校間分散)

ここでは個人レベルの切片に関わるパラメーターにはaをつけておきます。 少しややこしですが、これを求めてみると


b0\,:\,12\\
\sigma^{2}\,:\,16297\\
a0\,:\,1595\\
ab_{1}\,:\,25\\
 \sigma_{sc}^{2}\,:\,10723

となります。個人レベルでの傾きは12点と変化していません。しかし、学校レベルでの平均勉強時間の傾きが推定されていて25となっています。つまり、勉強時間が全体平均よりも長い学校にいる生徒は点数が高くなるということが分かります。

さらに学校種(公立、私立)も変数に加えてみましょう。ここでは、学校種によって点数に差があるかどうかを検討できます。

 Y_{st} = a_{sc} +b0*(X_{st}-\bar{X_{sc}}) + \epsilon \\
a_{sc} = a0 + ab_{1} * (\bar{X_{sc}} - \bar{X_{st}})+ab_{2} * SCtype+
\epsilon_{a}

\epsilon \sim Normal(0\,, \sigma^{2})\\
 \epsilon_{a} \sim Normal (0 \,, \sigma_{sc}^{2})

Y_{st}\,:\,生徒の点数\\
X_{st}\,:\,生徒の勉強時間\\
\bar{X_{sc}}\,:\,学校ごとの平均勉強時間\\
\bar{X_{st}}\,:\,勉強時間の全体平均\\
a_{sc}\,:\,学校ごとの切片\\
b0\,:\, 傾き(勉強時間の点数への影響)\\
\sigma^{2}\,:\,誤差の分散(学校内分散)\\
a0\,:\,切片平均\\
ab_{1}\,:\,学校レベルでの勉強時間の点数への影響\\
ab_{2}\,:\,学校種の勉強時間の点数への影響\\
SCtype\,:\,学校種(公立:0,私立:1)\\
\sigma_{sc}^{2}\,:\,a_{sc}の分散(学校間分散)

これを求めると


b0\,:\,12\\
\sigma^{2}\,:\,16315 \\
a0\,:\,1595\\
ab_{1}\,:\,22\\
ab_{2}\,:\,83(SE=201)\\
 \sigma_{sc}^{2}\,:\,10247

となりました。学校種の勉強時間の影響は、公立(0)に比べて私立(1)になると83点高くなります。

と言いたいところですが標準誤差が非常に大きく傾きが有意ではないため、学校種によって点数に違いがあるとは言えなさそうです。

これは学校種をモデルに加えても、学校間分散や学校内分散に大きな変化がみられなかったことからも推測できます。

このように、同じモデルの中で学校レベルと個人レベルの二つのレベルでデータを解釈することができます!

#勉強時間の学校平均に全体平均中心化を行う
st_df2 %>% 
  ungroup(school) %>% 
  mutate(X_sc_cgm = X_sc - mean(X))->st_df3

#実行する
st_df3 %>% 
  lme4::lmer(data=.,Y~X_cwc+X_sc_cgm+(1|school),REML = F) %>% 
  summary()

Linear mixed model fit by maximum likelihood [
lmerMod]
Formula: Y ~ X_cwc + X_sc_cgm + (1 | school)
   Data: .

     AIC      BIC   logLik deviance df.resid 
  1404.5   1418.0   -697.3   1394.5      105 

Scaled residuals: 
     Min       1Q   Median       3Q      Max 
-2.33690 -0.65150 -0.03606  0.65453  2.82370 

Random effects:
 Groups   Name        Variance Std.Dev.
 school   (Intercept) 10723    103.5   
 Residual             16297    127.7   
Number of obs: 110, groups:  school, 7

Fixed effects:
            Estimate Std. Error t value
(Intercept) 1595.514     42.135  37.867
X_cwc         12.063      1.395   8.650
X_sc_cgm      25.636      3.611   7.099

Correlation of Fixed Effects:
         (Intr) X_cwc 
X_cwc     0.000       
X_sc_cgm -0.085  0.000


#学校種追加

st_df3 %>% 
  lme4::lmer(data=.,Y~X_cwc+X_sc_cgm+sctype+(1|school),REML = F) %>% 
  summary()

Linear mixed model fit by maximum likelihood [
lmerMod]
Formula: 
Y ~ X_cwc + X_sc_cgm + sctype + (1 | school)
   Data: .

     AIC      BIC   logLik deviance df.resid 
  1406.4   1422.6   -697.2   1394.4      104 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-2.3264 -0.6497 -0.0284  0.6601  2.8163 

Random effects:
 Groups   Name        Variance Std.Dev.
 school   (Intercept) 10247    101.2   
 Residual             16315    127.7   
Number of obs: 110, groups:  school, 7

Fixed effects:
            Estimate Std. Error t value
(Intercept) 1564.128     86.921  17.995
X_cwc         12.063      1.395   8.645
X_sc_cgm      22.429      8.560   2.620
sctype        83.016    201.773   0.411

Correlation of Fixed Effects:
         (Intr) X_cwc  X_sc_c
X_cwc     0.000              
X_sc_cgm  0.784  0.000       
sctype   -0.880  0.000 -0.910

まとめ

今回はマルチレベル分析の必要性とその基本であるランダム化切片モデルをまとめました。

マルチレベル分析は、説明変数の中心化や結果の解釈などが少し複雑ですが、適切な解析をするには必須の分析法になると思います。

今回は扱わなかった、傾きにもランダム性を仮定したランダム切片傾きモデルや、複数のモデルから適切なモデルを選択する方法などもあるので、次回まとめていきます!

www.medi-08-data-06.work

※本記事は筆者が個人的に学んだことをまとめた記事なります。数学の記法や詳細な理論、筆者の勘違い等で誤りがあった際はご指摘頂けると幸いです。

参考

こちらの内容を参考にしました。マルチレベル分析の基礎から具体的な使い方がまとめてあり、とても分かりやすい内容になってます。実践的に使いたい方にはおすすめです!

Rで学ぶ マルチレベルモデル[入門編]: 基本モデルの考え方と分析

Rで学ぶ マルチレベルモデル[入門編]: 基本モデルの考え方と分析

こちらも同じくマルチレベル分析の方法について書かれていますが、今回扱った線形混合モデルに加え構造方程式モデルなども解説されていて、心理学、社会学をテーマにした内容になっています。

個人と集団のマルチレベル分析

個人と集団のマルチレベル分析