一般化線形混合効果モデル
一般化線形混合効果モデルとは
GLME (一般化線形混合効果) モデルは、正規分布以外の応答変数分布をもつデータの場合に、1 つ以上のグループ化変数によって変化する可能性のある係数を使用し、応答変数と独立変数との関係を記述します。GLME モデルは、グループに収集および要約されたデータに関する GLM (一般化線形モデル) の拡張と考えられます。あるいは、GLME モデルは、応答変数が正常に分布していないデータに関する LME (線形混合効果モデル) の汎化と考えられます。
混合効果モデルは固定効果および変量効果の項を含みます。固定効果の項は通常、モデルの従来の線形回帰部分です。変量効果の項は、母集団から無作為に抽出された個々の実験単位に関連付けられており、応答に影響を与える可能性のあるグループ間の変動を考慮しています。変量効果には事前分布がありますが、固定効果にはありません。
GLME モデルの式
一般化線形混合効果モデルの標準形式は次のとおりです。
ここで
y は n 行 1 列の応答ベクトルであり、yi はその i 番目の要素です。
b は変量効果ベクトルです。
Distr は、与えられる b に対する y の指定された条件付き分布です。
μ は与えられた b に対する y の条件付き平均であり、μi はその i 番目の要素です。
σ2 は分散パラメーターです。
w は有効な観測値の重みベクトルであり、wi は観測値 i の重みです。
二項分布では、有効な観測値の重みは、
fitglme
において'Weights'
によって名前と値のペアの引数が指定された事前の重みに等しく、'BinomialSize'
名前と値のペアの引数で指定された二項分布のサイズで乗算されています。他のすべての分布の場合、有効な観測値の重みは
fitglme
における'Weights'
の名前と値のペアの引数によって指定された事前の重みと等しくなります。
g(μ) は、平均応答 μ と予測子の線形結合との関係を定義するリンク関数です。
X は n 行 p 列の固定効果計画行列です。
β は p 行 1 列の固定効果ベクトルです。
Z は n 行 q 列の変量効果計画行列です。
b は q 行 1 列の変量効果ベクトルです。
δ はモデル オフセット ベクトルです。
平均応答 μ のモデルは以下です。
ここで、g-1 はリンク関数 g(μ) の逆関数であり、 は次の一般化線形混合効果モデルの固定効果および変量効果の線形予測子です。
GLME モデルは β、θ および σ2 でパラメーター表現されています。
一般化線形混合効果モデルの仮定は次のとおりです。
変量効果ベクトル b は以下の事前分布になります。
ここで、σ2 は分散パラメーターであり、D は、制約のないパラメーター ベクトル θ によりパラメーター表現された、対称な半正定値行列です。
観測値 yi は与えられた b に対して条件付きで独立しています。
モデル近似用データの準備
GLME モデルをデータに当てはめるには、fitglme
を使用します。table
データ型を使用して入力データをフォーマットします。表の各行は 1 つの観測値、各列は 1 つの予測子変数を表します。table
の作成および使用に関する詳細は、table の作成とその table へのデータの代入を参照してください。
入力データは連続したグループ化変数を含むことができます。fitglme
は、以下のデータ型を使用する予測子がカテゴリカルであると仮定します。
logical
categorical
文字ベクトルまたは文字配列
string 配列
文字ベクトルの cell 配列
入力データの表が何らかの NaN
値を含む場合、fitglme
はデータの行全体を近似から除外します。さらにデータ行を除外する場合、モデルを当てはめる際に fitglme
で 'Exclude'
の名前と値のペアの引数を使用できます。
モデルの分布タイプを選択
GLME モデルは、応答データが正規分布に従わないときに使用します。したがって、fitglme
を使用してモデルを当てはめる場合、'Distribution'
名前と値のペアの引数を使用して応答分布タイプを指定しなければなりません。多くの場合、応答データの型はモデルの適切な分布タイプを示します。
応答データの型 | 示される応答分布タイプ |
---|---|
任意の実数 | 'Normal' |
任意の正の数 | 'Gamma' または 'InverseGaussian' |
任意の非負の整数 | 'Poisson' |
ゼロから n までの整数 (n は正の固定値) | 'Binomial' |
モデルのリンク関数の選択
GLME モデルはリンク関数 g を使用し、平均応答と予測子の線形結合との間の関係をマップします。既定では fitglme
は、以下の表に示すように、応答データの指定された分布に基づく、事前定義済みの一般に受け入れられたリンク関数を使用します。ただし、fitglme
の 'Link'
名前と値のペアの引数を使用して、事前定義済みの関数のリストから異なるリンク関数を指定することもでき、また独自の関数を定義することもできます。
値 | 説明 |
---|---|
'comploglog' | g(mu) = log(-log(1-mu)) |
'identity' |
正規分布の正準リンク。 |
'log' |
ポアソン分布の正準リンク。 |
'logit' |
二項分布の正準リンク。 |
'loglog' | g(mu) = log(-log(mu)) |
'probit' | g(mu) = norminv(mu) |
'reciprocal' | g(mu) = mu.^(-1) |
スカラー値 P | g(mu) = mu.^P |
構造体 S | 値が関数ハンドルである 4 つのフィールドを含む構造体は、以下を扱います。
|
モデルをデータに当てはめる場合、fitglme
は既定で正準リンク関数を使用します。
分布 | 既定のリンク関数 |
---|---|
'Normal' | 'identity' |
'Binomial' | 'logit' |
'Poisson' | 'log' |
'Gamma' | -1 |
'InverseGaussian' | -2 |
リンク関数 'comploglog'
、'loglog'
および 'probit'
は主に二項モデルで便利です。
モデル式の指定
fitglme
のモデル仕様では、'y ~ terms'
という形式の文字ベクトルまたは string スカラーであるウィルキンソンの表記法を使用します。ここで y
は応答変数名であり、terms
は以下の表記法で記述されます。
ウィルキンソンの表記法 | 標準表記の因子 |
---|---|
1 | 定数 (切片) 項 |
X^k 、k は正の整数 | X , X2 , ..., Xk |
X1 + X2 | X1 , X2 |
X1*X2 | X1 , X2 , X1.*X2 (element-wise multiplication of X1 and X2) |
X1:X2 | X1.*X2 のみ |
- X2 | X2 は含めない |
X1*X2 + X3 | X1 , X2 , X3 , X1*X2 |
X1 + X2 + X3 + X1:X2 | X1 , X2 , X3 , X1*X2 |
X1*X2*X3 - X1:X2:X3 | X1 , X2 , X3 , X1*X2 , X1*X3 , X2*X3 |
X1*(X2 + X3) | X1 , X2 , X3 , X1*X2 , X1*X3 |
式には既定で定数 (切片) 項が含まれます。モデルから定数項を除外するには、式に –1
を含めます。
一般化線形混合効果モデルについては、式の仕様は 'y ~ fixed + (random1|grouping1) + ... + (randomR|groupingR)'
形式であり、ここで fixed
および random
は固定効果項と変量効果項をそれぞれ含みます。
入力データの表は以下を含むと想定します。
応答変数
y
予測子変数
X1
,X2
, ...,XJ
。ここで J は予測子変数の (連続したグループ化変数を含む) 合計数。グループ化変数
g1
,g2
, ...,gR
。ここで R はグループ化変数の数。
XJ
および gR
のグループ化変数は、categorical 配列、logical 配列、文字配列、string 配列、または文字ベクトルの cell 配列が可能です。
この場合、'y ~ fixed + (random1|g1) + ... + (randomR|gR)'
の形式の式において、項 fixed
は固定効果の計画行列 X
の仕様に対応し、random1
はグループ化変数 g1
に対応する変量効果の計画行列 Z1
の仕様であり、同様に randomR
はグループ化変数 gR
に対応する変量効果の計画行列 ZR
の仕様です。fixed
項および random
項は、以下のようにウィルキンソンの表記法で表現できます。
式 | 説明 |
---|---|
'y ~ X1 + X2' | 切片 X1 および X2 の固定効果。この式は 'y ~ 1 + X1 + X2' と等価です。 |
'y ~ -1 + X1 + X2' | X1 および X2 に対する固定効果を伴う切片はありません。-1 を含めることによって暗黙的な切片の項は抑制されます。 |
'y ~ 1 + (1 | g1)' | グループ化変数 g1 の水準ごとの切片の固定効果と切片の変量効果の和。 |
'y ~ X1 + (1 | g1)' | 固定傾きのランダム切片モデル。 |
'y ~ X1 + (X1 | g1)' | 相関があり得るランダムな切片と傾き。この式は 'y ~ 1 + X1 + (1 + X1|g1)' と等価です。 |
'y ~ X1 + (1 | g1) + (-1 + X1 | g1)' | 切片および傾きに関する独立した変量効果の項。 |
'y ~ 1 + (1 | g1) + (1 | g2) + (1 | g1:g2)' | g1 と g2 に対する独立したメイン効果のあるランダムな切片モデル + 独立した交互作用効果。 |
たとえば、標本データ mfr
は世界で 50 の工場を操業する製造企業からのシミュレートされたデータを含みます。各工場は、完成品の生産のためにバッチ プロセスを実行します。同社は各バッチの欠陥数を減少させるために新たな製造プロセスを開発しました。新しいプロセスの効果をテストするため、同社は実験に参加させる 20 工場を無作為に選びました。10 工場では新プロセスを実施しますが、残りの 10 工場では旧プロセスの実行を続けます。各 20 工場で、同社は 5 つのバッチ (合計 100 バッチ) を実行し、バッチごとに処理時間 (time_dev
)、温度 (temp_dev
)、欠陥数 (defects
) および原材料の供給業者を示すカテゴリカル変数 (supplier
) のデータを記録しました。
新しいプロセス (予測子変数 newprocess
で表現されます) が欠陥数を有意に減少させたかどうかを判断するために、newprocess
、time_dev
、temp_dev
および supplier
を固定効果予測子として使用し、GLME モデルを当てはめます。factory
別にグループ化された変量効果の切片を含めます。これは工場特有の変動に起因する品質に差がある可能性を考慮するためです。応答変数 defects
はポアソン分布になります。
欠陥数はポアソン分布を使用してモデル化できます
これは一般化線形混合効果モデルに対応します
ここで
defectsij は、バッチ j (ここで j = 1, 2, ..., 5) 処理中に工場 i (ここで i = 1, 2, ..., 20) で実行されたバッチにおいて観察された欠陥数です。
μij はバッチ j 処理中の工場 i に対応する欠陥の平均数です。
supplier_Cij および supplier_Bij はそれぞれ、会社
C
またはB
がバッチ j 処理中に工場 i で実行されたバッチ用の加工化学薬品を供給したかどうかを意味するダミー変数です。bi ~ N(0,σb2) は、工場特有の品質変動に相当する、各工場 i の変量効果の切片です。
ウィルキンソンの表記法を使用し、このモデルを以下のように指定します。
'defects ~ 1 + newprocess + time_dev + temp_dev + supplier + (1|factory)'
fitglme
を使用してモデルをデータに当てはめる際、応答変数のポアソン分布を考慮するには、'Distribution'
の名前と値のペアの引数を 'Poisson'
として指定します。既定では、fitglme
はポアソン分布の応答変数に対数リンク関数を使用します。
モデルの表示
近似関数 fitglme
の出力は一般化線形混合効果モデルの情報を提供します。
製造実験データ mfr
を使用し、newprocess
、time_dev
、temp_dev
および supplier
を固定効果予測子とするモデルを当てはめます。応答分布にポアソンを、リンク関数に対数を、近似メソッドにラプラスを指定します。
load mfr glme = fitglme(mfr,... 'defects ~ 1 + newprocess + time_dev + temp_dev + supplier + (1|factory)',... 'Distribution','Poisson','Link','log','FitMethod','Laplace',... 'DummyVarCoding','effects')
glme = Generalized linear mixed-effects model fit by ML Model information: Number of observations 100 Fixed effects coefficients 6 Random effects coefficients 20 Covariance parameters 1 Distribution Poisson Link Log FitMethod Laplace Formula: defects ~ 1 + newprocess + time_dev + temp_dev + supplier + (1 | factory) Model fit statistics: AIC BIC LogLikelihood Deviance 416.35 434.58 -201.17 402.35 Fixed effects coefficients (95% CIs): Name Estimate SE tStat DF pValue '(Intercept)' 1.4689 0.15988 9.1875 94 9.8194e-15 'newprocess' -0.36766 0.17755 -2.0708 94 0.041122 'time_dev' -0.094521 0.82849 -0.11409 94 0.90941 'temp_dev' -0.28317 0.9617 -0.29444 94 0.76907 'supplier_C' -0.071868 0.078024 -0.9211 94 0.35936 'supplier_B' 0.071072 0.07739 0.91836 94 0.36078 Lower Upper 1.1515 1.7864 -0.72019 -0.015134 -1.7395 1.5505 -2.1926 1.6263 -0.22679 0.083051 -0.082588 0.22473 Random effects covariance parameters: Group: factory (20 Levels) Name1 Name2 Type Estimate '(Intercept)' '(Intercept)' 'std' 0.31381 Group: Error Name Estimate 'sqrt(Dispersion)' 1
Model information
表は標本データの観測値の合計 (100)、固定効果および変量効果係数の数 (それぞれ 6 および 20)、共分散パラメーターの数 (1) を表示しています。また、応答変数は Poisson
分布であり、リンク関数は Log
であり、近似メソッドが Laplace
であることもわかります。
Formula
はウィルキンソンの表記法によるモデル仕様を示します。
Model fit statistics
表はモデルの適合度の評価に使用された統計を表します。これには赤池情報量基準 (AIC
)、ベイズ情報量基準 (BIC
) 値、対数尤度 (LogLikelihood
) および逸脱度 (Deviance
) の値が含まれます。
Fixed effects coefficients
表は、fitglme
が 95% の信頼区間を返したことを示します。これには固定効果予測子ごとに 1 行が含まれ、各列にはその予測子に対応する統計が含まれます。列 1 (Name
) には各固定効果係数の名前が含まれ、列 2 (Estimate
) にはその推定値が含まれ、列 3 (SE
) には係数の標準誤差が含まれます。列 4 (tStat
) には係数が 0 に等しいという仮説検定のための t 統計量が含まれます。列 5 (DF
) と列 6 (pValue
) には t 統計量に対応する自由度と p 値がそれぞれ含まれます。最後の 2 列 (Lower
および Upper
) には、各固定効果係数の 95% 信頼区間の下限と上限がそれぞれ表示されます。
Random effects covariance parameters
は各グループ化変数 (ここでは factory
のみ) の表を表示します。これには水準の総数 (20)、共分散パラメーターの型および推定値が含まれます。ここでの std
は、工場の予測子に関連付けられている変量効果の標準偏差が fitglme
から返されることを示します。この推定値は 0.31381 です。また、誤差パラメーターの型 (ここでは分散パラメーターの平方根) およびその推定値 1 を含む表も表示します。
fitglme
により生成される標準表示は変量効果パラメーターの信頼区間を指定しません。covarianceParameters
を使用して、これらの値を計算し表示します。
モデルの操作
fitglme
を用いて GLME モデルを作成した後で、モデルを操作するための追加の関数を使用できます。
係数および信頼区間の検査および検定
固定効果および変量効果係数の推定値、共分散パラメーター、計画行列および関連する統計を抽出するには、次の方法を用います。
fixedEffects
は、推定固定効果係数および関連する統計を近似されたモデルから抽出します。関連する統計には、標準誤差、t 統計量、自由度、各パラメーターが 0 に等しいかどうかを仮説検定するための p 値、および信頼区間が含まれます。randomEffects
は推定変量効果係数および関連する統計を近似された GLME モデルから抽出します。関連する統計には、それぞれの変量効果の推定 EBP (経験的ベイズ予測子)、共分散パラメーターおよび応答が与えられた CMSEP (予測の条件付き平均二乗誤差) の平方根、t 統計量、自由度、および各パラメーターが 0 に等しいかどうかを仮説検定するための p 値、ならびに信頼区間が含まれます。covarianceParameters
は推定共分散パラメーターおよび関連する統計を近似された GLME モデルから抽出します。関連する統計には、共分散パラメーターの推定値および信頼区間が含まれます。designMatrix
は固定効果および変量効果の計画行列またはその指定されたサブセットを、近似された GLME モデルから抽出します。
固定効果および変量効果係数の有意性に関して、カスタマイズされた仮説検定を実行し、カスタムの信頼区間を計算するには、以下を行います。
新たな応答値を生成しモデルを再度当てはめる
近似された GLME モデルに基づき、近似された応答、予測された応答および無作為の応答を含む新たな応答値を生成するには、以下を行います。
残差の検査および可視化
近似された GLME モデルから残差を抽出および可視化するには、以下を行います。
residuals
は、近似されたモデルから生の残差もしくはピアソン残差を抽出します。さらに、条件付き残差または限界残差を計算するかどうかも指定できます。plotResiduals
は、近似されたモデルから生の残差またはピアソン残差により、次のプロットを作成します。残差のヒストグラム
残差と近似値の散布図
残差とラグ付き残差の散布図
参考
fitglme
| GeneralizedLinearMixedModel