fitted

線形混合効果モデルからの近似応答

構文

yfit = fitted(lme)

yfit = fitted(lme,Name,Value)

説明

yfit = fitted(lme) は、線形混合効果モデル lme から近似条件付き応答を返します。

yfit = fitted(lme,Name,Value) は、1 つ以上の Name,Value のペア引数によって指定された追加オプションを使用して、線形混合効果モデル lme から近似応答を返します。

たとえば、近似限界応答を計算するかどうかを指定できます。

例

近似条件付き応答と近似周辺応答の計算

標本データを読み込みます。

load flu

データセット配列 flu には、変数 Date と、インフルエンザ推定罹患率 (Google® 検索から推定される 9 地域の値と疾病対策センター (CDC) による全国の推定値) が格納されている 10 個の変数が含まれています。

線形混合効果モデルを当てはめるには、データが適切な形式のデータセット配列になっていなければなりません。インフルエンザ罹患率を応答として、地域を予測子変数として線形混合効果モデルを当てはめるため、地域に対応する 9 個の列を 1 つの配列にまとめます。新しいデータセット配列 flu2 には、応答変数 FluRate、各推定の元になっている地域を示すノミナル変数 Region およびグループ化変数 Date が含まれなければなりません。

flu2 = stack(flu,2:10,'NewDataVarName','FluRate','IndVarName','Region');
flu2.Date = nominal(flu2.Date);

地域に対する固定効果と、Date で変化するランダム切片で、線形混合効果モデルを当てはめます。

地域はカテゴリカル変数です。モデルを当てはめるときに、名前と値のペア引数 DummyVarCoding を使用してカテゴリカル変数の対比を指定できます。対比を指定しない場合、fitlme は 'reference' 対比を既定で使用します。モデルは切片をもつので、fitlme は最初の地域 NE を参照として受け取り、他の 8 つの地域を表す 8 つのダミー変数を作成します。たとえば、$$I[MidAtl]$$ は地域 MidAtl を表すダミー変数です。詳細については、ダミー変数を参照してください。

対応するモデルは以下のとおりです。

$$\begin{array}{l}{y}_{im}={\beta }_0+{\beta }_{1}I{\left[MidAtl\right]}_i+{\beta }_{2}I{\left[ENCentral\right]}_i+{\beta }_{3}I{\left[WNCentral\right]}_i+{\beta }_{4}I{\left[SAtl\right]}_i\\ +{\beta }_{5}I{\left[ESCentral\right]}_i+{\beta }_{6}I{\left[WSCentral\right]}_i+{\beta }_{7}I{\left[Mtn\right]}_i+{\beta }_{8}I{\left[Pac\right]}_i+b_{0m}+{\epsilon }_{im},m=1,2,...,52,\end{array}$$

ここで、$$y_{im}$$ はグループ化変数 Date の水準 $$m$$ に対応する観測値 $$i$$、$${\beta }_j$$ は固定効果係数 ($$j$$ = 0、1、...、8)、$${\beta }_0$$ は地域 NE の係数です。$$b_{0m}$$ はグループ化変数 Date の水準 $$m$$ に対応する変量効果、$${\epsilon }_{im}$$ は観測値 $$i$$ の観測誤差です。変量効果の事前分布は $$b_{0m}\sim N(0,{\sigma }_b^2)$$、誤差項の分布は $${\epsilon }_{im}\sim N(0,{\sigma }^2)$$ です。

lme = fitlme(flu2,'FluRate ~ 1 + Region + (1|Date)')

lme = 
Linear mixed-effects model fit by ML

Model information:
    Number of observations             468
    Fixed effects coefficients           9
    Random effects coefficients         52
    Covariance parameters                2

Formula:
    FluRate ~ 1 + Region + (1 | Date)

Model fit statistics:
    AIC       BIC       LogLikelihood    Deviance
    318.71    364.35    -148.36          296.71  

Fixed effects coefficients (95% CIs):
    Name                        Estimate    SE          tStat      DF     pValue        Lower        Upper    
    {'(Intercept)'     }          1.2233    0.096678     12.654    459     1.085e-31       1.0334       1.4133
    {'Region_MidAtl'   }        0.010192    0.052221    0.19518    459       0.84534    -0.092429      0.11281
    {'Region_ENCentral'}        0.051923    0.052221     0.9943    459        0.3206    -0.050698      0.15454
    {'Region_WNCentral'}         0.23687    0.052221     4.5359    459    7.3324e-06      0.13424      0.33949
    {'Region_SAtl'     }        0.075481    0.052221     1.4454    459       0.14902     -0.02714       0.1781
    {'Region_ESCentral'}         0.33917    0.052221      6.495    459    2.1623e-10      0.23655      0.44179
    {'Region_WSCentral'}           0.069    0.052221     1.3213    459       0.18705    -0.033621      0.17162
    {'Region_Mtn'      }        0.046673    0.052221    0.89377    459       0.37191    -0.055948      0.14929
    {'Region_Pac'      }        -0.16013    0.052221    -3.0665    459     0.0022936     -0.26276    -0.057514

Random effects covariance parameters (95% CIs):
Group: Date (52 Levels)
    Name1                  Name2                  Type           Estimate    Lower     Upper  
    {'(Intercept)'}        {'(Intercept)'}        {'std'}        0.6443      0.5297    0.78368

Group: Error
    Name               Estimate    Lower      Upper
    {'Res Std'}        0.26627     0.24878    0.285

$$p$$ 値 7.3324e-06 および 2.1623e-10 は、地域 WNCentral および ESCentral のインフルエンザ罹患率の固定効果が地域 NE のインフルエンザ罹患率に対して有意に異なることをそれぞれ示しています。

変量効果項の標準偏差 $${\sigma }_b$$ の信頼限界 (0.5297, 0.78368) には、変量効果項が有意であることを示す 0 が含まれていません。compare メソッドを使用して、変量効果の項の有意性をテストすることもできます。

任意の観測におけるモデルからの近似条件付き応答には、固定効果および変量効果からの寄与が含まれます。たとえば、2005 年 10 月 9 日の週における地域 WNCentral の罹患率に対する最良線形不偏予測量 (BLUP) の推定値は次のようになります。

$$\begin{array}{l}{\underset{}{\overset{ˆ}{y}}}_{WNCentral,10/9/2005}={\underset{}{\overset{ˆ}{\beta }}}_0+{\underset{}{\overset{ˆ}{\beta }}}_{3}I\left[WNCentral\right]+{\underset{}{\overset{ˆ}{b}}}_{10/9/2005}\\ =1.2233+0.23687-0.1718\\ =1.28837.\end{array}$$

固定効果と変量効果の両方による推定に対する寄与が含まれているので、これは近似された条件付き応答です。この値を以下のように計算することもできます。

beta = fixedEffects(lme);
[~,~,STATS] = randomEffects(lme); % Compute the random-effects statistics (STATS)
STATS.Level = nominal(STATS.Level);
y_hat = beta(1) + beta(4) + STATS.Estimate(STATS.Level=='10/9/2005')

y_hat = 1.2884

上記の計算で、beta(1) は $${\beta }_0$$ の推定値に、beta(4) は $${\beta }_3$$ の推定値に対応します。fitted メソッドを使用して、近似値を簡単に表示できます。

F = fitted(lme);
F(flu2.Date == '10/9/2005' & flu2.Region == 'WNCentral')

ans = 1.2884

2005 年 10 月 9 日の週における地域 WNCentral の推定限界応答は次のようになります。

$$\begin{array}{l}{\underset{}{\overset{ˆ}{y}}}_{WNCentral,10/9/2005}^{(marginal)}={\underset{}{\overset{ˆ}{\beta }}}_0+{\underset{}{\overset{ˆ}{\beta }}}_{3}I\left[WNCentral\right]\\ =1.2233+0.23687\\ =1.46017.\end{array}$$

近似限界応答を計算します。

F = fitted(lme,'Conditional',false);
F(flu2.Date == '10/9/2005' & flu2.Region == 'WNCentral')

ans = 1.4602

残差と近似値のプロット

標本データを読み込みます。

load('weight.mat');

weight には長期間の調査によるデータが含まれています。そこには 20 人の被験者が 4 つの運動プログラムにランダムに割り当てられ、体重の減少が 6 回の 2 週間の期間にわたって記録されています。このデータは、シミュレーションされたものです。

データをテーブルに保存します。Subject および Program をカテゴリカル変数として定義します。

tbl = table(InitialWeight,Program,Subject,Week,y);
tbl.Subject = nominal(tbl.Subject);
tbl.Program = nominal(tbl.Program);

線形混合効果モデルを当てはめます。初期体重、プログラムの種類、週、週とプログラムの種類の間の交互作用は固定効果です。切片と週は被験者ごとに異なります。

lme = fitlme(tbl,'y ~ InitialWeight + Program*Week + (Week|Subject)');

近似値および生の残差を計算します。

F = fitted(lme);
R = residuals(lme);

残差と近似値の対比をプロットします。

plot(F,R,'bx')
xlabel('Fitted Values')
ylabel('Residuals')

次に、プログラム別にグループ化された、残差と近似値の対比をプロットします。

figure()
gscatter(F,R,Program)

入力引数

lme — 線形混合効果モデル

LinearMixedModel オブジェクト

線形混合効果モデル。fitlme または fitlmematrix を使用して構築した LinearMixedModel オブジェクトとして指定します。

名前と値の引数

オプションの引数のペアを Name1=Value1,...,NameN=ValueN として指定します。ここで Name は引数名、Value は対応する値です。名前と値の引数は他の引数の後ろにする必要がありますが、ペアの順序は関係ありません。

R2021a より前では、名前と値をそれぞれコンマを使って区切り、Name を引用符で囲みます。

例: yfit = fitted(lme,'Conditional',false)

Conditional — 条件付き応答のインジケーター

true (既定値) | false

条件付き応答のインジケーター。'Conditional' と以下のいずれかで構成されるコンマ区切りのペアとして指定します。

true	固定効果と変量効果の両方からの寄与 (条件付き)
false	固定効果のみからの寄与 (限界)

例: 'Conditional',false

データ型: logical

出力引数

yfit — 近似応答値

n 行 1 列のベクトル

n 行 1 列のベクトルとして返される、近似応答値。ここで、n は観測値の数です。

詳細

近似条件付き応答と近似限界応答

条件付き応答には、固定効果と変量効果の両方からの寄与が含まれますが、限界応答には固定効果からの寄与のみが含まれます。

線形混合効果モデル lme には、n 行 p 列の固定効果の計画行列 X と、n 行 q 列の変量効果の計画行列 Z があるとします。また、推定した p 行 1 列の固定効果ベクトルが $$\widehat{\beta }$$、変量効果について推定した q 行 1 列の最良線形不偏予測量 (BLUP) ベクトルが $$\widehat{b}$$ であるとします。近似条件付き応答は以下のようになります。

$${\widehat{y}}_{Cond}=X\widehat{\beta }+Z\widehat{b},$$

また、近似限界応答は以下のようになります。

$${\widehat{y}}_{Mar}=X\widehat{\beta },$$

バージョン履歴

R2013b で導入

参考

LinearMixedModel | residuals | response