ドキュメンテーション

最新のリリースでは、このページがまだ翻訳されていません。 このページの最新版は英語でご覧になれます。

covarianceParameters

クラス: LinearMixedModel

線形混合効果モデルの共分散パラメーターの抽出

構文

  • psi = covarianceParameters(lme)
  • [psi,mse] = covarianceParameters(lme)
  • [psi,mse,stats] = covarianceParameters(lme)
  • [psi,mse,stats] = covarianceParameters(lme,Name,Value)

説明

psi = covarianceParameters(lme) は、変動効果の前の共分散をパラメーター化する、推定された共分散パラメーターを返します。

[psi,mse] = covarianceParameters(lme) は、残差分散の推定も返します。

[psi,mse,stats] = covarianceParameters(lme) は、共分散パラメーターおよび関連する統計を含むセル配列 stats も返します。

[psi,mse,stats] = covarianceParameters(lme,Name,Value) は、1 つ以上の Name,Value のペアの引数によって指定された追加オプションで stats の共分散パラメーターおよび関連する統計を返します。

たとえば、共分散パラメーターの信頼限界に対する信頼度を指定できます。

入力引数

すべて展開する

LinearMixedModel オブジェクトとして返される線形混合効果モデル。

このオブジェクトのプロパティとメソッドについては、LinearMixedModel を参照してください。

名前/値のペアの引数

オプションの Name,Value の引数ペアをコンマ区切りで指定します。ここで、Name は引数名で、Value は対応する値です。Name は単一引用符 (' ') で囲まなければなりません。Name1,Value1,...,NameN,ValueN のように、複数の名前と値のペアの引数を任意の順序で指定できます。

すべて展開する

信頼度。'Alpha' と 0 ~ 1 の範囲にあるスカラー値をコンマ区切りのペアとして指定します。値が α の場合、信頼度は 100*(1-α)% です。

たとえば、99% の信頼区間の場合は、次のように信頼度を指定できます。

例: 'Alpha',0.01

データ型: single | double

出力引数

すべて展開する

変動効果の前の共分散をパラメーター化する共分散パラメーターの推定。長さ R のセル配列として返されます。psi{r} には、グループ化変数 gr (r = 1, 2, ..., R) に関連付けられている変動効果の共分散行列が含まれます。グループ化変数の順序は、モデルを近似するときに入力する順序と同じです。

残差分散の推定。スカラー値として返します。

共分散パラメーターの推定および関連する統計。以下の列で構成されるデータセット配列を含む長さ (R + 1) のセル配列として返します。

Groupグループ化変数名
Name1最初の予測子変数の名前
Name22 番目の予測子変数の名前
Type std (標準偏差)。Name1Name2 が同じ場合

corr (相関)。Name1Name2 が異なる場合

EstimateName1Name2 が同じ場合は、予測子 Name1 または Name2 と関連付けられている変量効果の標準偏差

Name1Name2 が異なる場合は、予測子 Name1 および Name2 と関連付けられている変量効果の間の相関

Lower共分散パラメーターの 95% 信頼区間の下限
Upper共分散パラメーターの 95% 信頼区間の上限

stats{r} は、r 番目のグループ化変数 (r = 1、2、...、R) の共分散パラメーターに対する統計を含むデータセット配列です。stats{R+1} には、残差標準偏差の統計値が格納されます。残差誤差のデータセット配列には、GroupNameEstimateLowerUpper の各フィールドがあります。

すべて展開する

標本データが含まれたフォルダーに移動します。

cd(matlabroot)
cd('help/toolbox/stats/examples')

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

load fertilizer

このデータセット配列には土壌の種類に基づいて土壌が 3 つのブロックに分けられている分割プロット試験のデータが含まれています。土壌の種類は砂質、シルトおよび粘土質です。各ブロックは 5 つのプロットに分割され、5 種類のトマトの苗木 (チェリー、エアルーム、グレープ、枝付き、プラム) がランダムにこれらのプロットに割り当てられます。その後、プロット内のトマトの苗木はサブプロットに分割され、それぞれのサブプロットが 4 つの肥料の中の 1 つにより処置されます。このデータは、シミュレーションされたものです。

実用目的でこのデータを ds という名前のデータセット配列に保存し、TomatoSoil および Fertilizer をカテゴリカル変数として定義します。

ds = fertilizer;
ds.Tomato = nominal(ds.Tomato);
ds.Soil = nominal(ds.Soil);
ds.Fertilizer = nominal(ds.Fertilizer);

線形混合効果モデルを近似します。Fertilizer は固定効果変数であり、平均収穫量はブロック (土壌の種類) とブロック内のプロット (土壌の種類の中のトマトの種類) によって独立して変化します。このモデルは以下の式に対応します。このモデルは次の式に対応します。

yijk=β0+j=25β2jI[T]ij+b0kSk+b0jk(S*T)jk+εijk,

ここで、i = 1, 2, ..., 60 は観測値に、j = 2, ..., 5 はトマトの種類に、k = 1, 2, 3 はブロック (土壌) に対応します。Sk は k 番目の土壌の種類を、(S*T)jk は k 番目の土壌の種類で入れ子にされた j 番目のトマトの種類を表します。I[T]ij は、トマトの種類のレベル j を表すダミー変数です。

変量効果と観測値の誤差の事前分布は次のとおりです。b0k~N(0,σ2S)、b0jk~N(0,σ2S*T) および εijk ~ N(0,σ2)。

lme = fitlme(ds,'Yield ~ Fertilizer + (1|Soil) + (1|Soil:Tomato)');

変量効果項の共分散パラメーターの推定 (σ2S および σ2S*T の推定) を計算します。

psi = covarianceParameters(lme)
psi = 

    [4.4026e-17]
    [  352.8481]

残差分散 (σ2) を計算します。

[~,mse] = covarianceParameters(lme)
mse =

  151.9007

標本データが含まれたフォルダーに移動します。

cd(matlabroot)
cd('help/toolbox/stats/examples')

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

load weight

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

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

ds = dataset(InitialWeight,Program,Subject,Week,y);
ds.Subject = nominal(ds.Subject);
ds.Program = nominal(ds.Program);

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

'reference' ダミー変数コーディングの場合、fitlme はプログラム A を基準として使用して、必須のダミー変数 I[.] を作成します。このモデルは以下の式に対応します。

yim=β0+β1IWi+β2Weeki+β3I[PB]i+β4I[PC]i+β5I[PD]i+b+0mb1mWeekim+εim,

ここで、i は観測番号に対応し、i = 1, 2, ...,120 です。m は被験者番号に対応し、m = 1, 2, ..., 20 です。βj は固定効果係数で、j = 0, 1, ..., 8、b0m および b1m は変量効果です。IW は初期体重を表し、I[.] はプログラムの種類を表すダミー変数です。たとえば、I[PB]i はプログラム B を表すダミー変数です。

変量効果と観測値の誤差の事前分布は次のとおりです。(b0mb1m)~N(0,(σ02σ0,1σ0,1σ12)) および εim ~ N(0,σ2)。

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

変量効果の共分散パラメーターの推定を計算します。

[psi,mse,stats] = covarianceParameters(lme)
psi = 

    [2x2 double]


mse =

    0.0105


stats = 

    [3x7 classreg.regr.lmeutils.titleddataset]
    [1x5 dataset                             

mse は推定された残差分散です。これは σ2 に対する推定です。

変量効果の項 (σ20、σ21 および σ20,1) に対する共分散パラメーターの推定を確認するには、psi にインデックスを指定します。

psi{1}
ans =

    0.0572    0.0490
    0.0490    0.0624           

変量効果の項の切片 σ20 に推定される分散は 0.0572 です。週に対する変量効果の項 σ21 に推定される分散は 0.0624 です。切片と週に対する変量効果の項 σ0,1 に推定される共分散は 0.0490 です。

stats は 2 行 1 列のセル配列です。stats の最初のセルには、変量効果の標準偏差の信頼区間と、切片と週に関する変量効果間の相関が含まれます。これらを表示するには、stats のインデックスを指定します。

stats{1}
ans = 


    Covariance Type: FullCholesky

    Group      Name1                Name2                Type          Estimate    Lower      Upper  
    Subject    '(Intercept)'        '(Intercept)'        'std'         0.23927     0.14364    0.39854
    Subject    'Week'               '(Intercept)'        'corr'        0.81971     0.38662    0.95658
    Subject    'Week'               'Week'               'std'          0.2497     0.18303    0.34067       

表示では、グループ パラメーター (Group)、変量効果変数 (Name1Name2)、共分散パラメーターの種類 (Type)、各パラメーターの推定 (Estimate)、パラメーターの 95% 信頼区間 (Lower, Upper) が示されています。この表の推定は、次のようにして psi の推定に関連付けられます。

切片の変量効果の項の標準偏差は、0.23927 = sqrt(0.0527) です。同様に、週の変量効果の項の標準偏差は、0.2497 = sqrt(0.0624) です。最終的に、切片と週の変量効果の項の間にある相関は、0.81971 = 0.0490/(0.23927*0.2497) です。

この表示では、モデルを近似するときに使用された共分散パターンも示されていることに注意してください。この例では、共分散パターンは FullCholesky です。変量効果の項の共分散パターンを変更するには、モデルを近似するときに 'CovariancePattern' 名前と値のペアの引数を使用しなければなりません。

stats の 2 番目のセルには、残差標準偏差に関する同様の統計が含まれます。2 番目のセルの内容を表示します。

stats{2}
ans = 

    Group    Name             Estimate    Lower       Upper  
    Error    'Res Std'        0.10261     0.087882    0.11981

残差標準偏差の推定は、mse の平方根 0.10261 = sqrt(0.0105) です。

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

load carbig

ガロンあたりの走行マイル数 (MPG) の線形混合効果モデルを近似します。加速度および重量は固定効果で、モデル年度によってグループ化される切片と加速度は相関された変量効果の可能性があり、重量は独立変量効果であり、自動車の原点でグループ化されます。このモデルは次の式に対応します。

MPGimk=β0+β1Acci+β2Weighti+b10m+b11mAcci+b21kWeighti+εimk,m=1,2,...,13,k=1,2,...,8,

ここで、m は変数 Model_Year のレベルを、k は変数 Origin のレベルを表します。MPGimk は i 番目の観測値、m 番目のモデル年度および k 番目の生産国 (i 番目の観測値に対応) におけるガロンあたりの走行マイル数です。変量効果項と観測値の誤差の事前分布は次のとおりです。

b1m=(b10mb11m)~N(0,(σ102σ10,11σ10,11σ112)),b2k~N(0,σ22),εimk~N(0,σ2).

ここで変量効果項 b1m は、1 番目のグループ化変数のレベル m における 1 番目の変量効果を表します。変量効果の項 b10m は、最初のグループ化変数の m 番目のレベル (m) における、切片 (0) に対する、最初の変量効果の項 (1) に対応します。同じく、b11m は、最初の変量効果の項 (1) での最初の予測子 (1) のレベル m です。

同様に、b2k は、2 番目のグループ化変数のレベル k における 2 番目の変量効果の項を表します。

σ210 は切片に関する変量効果項の分散、σ211 は予測子の加速度に関する変量効果項の分散、σ10,11 は切片と予測子の加速度に関する変量効果項の共分散です。σ22 は 2 番目の変量効果項の分散で、σ2 は残差分散です。

最初に、線形混合効果モデルを近似するための計画行列を準備します。

X = [ones(406,1) Acceleration Weight];
Z = {[ones(406,1) Acceleration],[Weight]};
Model_Year = nominal(Model_Year);
Origin = nominal(Origin);
G = {Model_Year,Origin};

計画行列を使用してモデルを近似します。

lme = fitlmematrix(X,MPG,Z,G,'FixedEffectPredictors',....
{'Intercept','Acceleration','Weight'},'RandomEffectPredictors',...
{{'Intercept','Acceleration'},{'Weight'}},'RandomEffectGroups',{'Model_Year','Origin'});

変量効果の共分散パラメーターの推定を計算します。

[psi,mse,stats] = covarianceParameters(lme)
psi = 

    [2x2 double]
    [6.7989e-08]


mse =

    9.0755


stats = 

    [3x7 classreg.regr.lmeutils.titleddataset]
    [1x7 classreg.regr.lmeutils.titleddataset]
    [1x5 dataset                             ]

残差分散 mse は 9.0755 です。psi は 2 行 1 列のセル配列であり、stats は 3 行 1 列のセル配列です。内容を確認するには、これらのセル配列のインデックスを指定しなければなりません。

まず、psi の最初のセルのインデックスを指定します。

psi{1}
ans =

    8.5160   -0.8387
   -0.8387    0.1087

psi の 1 番目のセルには、切片 σ210 (8.5160) および加速度 σ211 (0.1087) の相関された変量効果の共分散パラメーターが含まれます。切片と加速度の変量効果項の共分散に対する推定 σ10,11 は -0.8387 です。

次に、psi の 2 番目のセルのインデックスを指定します。

psi{2}
ans =

   6.7989e-08

psi の 2 番目のセルには、重量 σ22 に関する変量効果の項での分散推定が含まれます。

stats の最初のセルのインデックスを指定します。

stats{1}
ans = 


    Covariance Type: FullCholesky

    Group         Name1                 Name2                 Type          Estimate    Lower       Upper   
    Model_Year    'Intercept'           'Intercept'           'std'           2.9182      1.1552      7.3716
    Model_Year    'Acceleration'        'Intercept'           'corr'        -0.87172    -0.98267    -0.30082
    Model_Year    'Acceleration'        'Acceleration'        'std'          0.32968     0.18863     0.57619

この表では、切片と加速度に関する変量効果の項の標準偏差が推定して示されています。標準偏差の推定は、psi の最初のセルの対角要素の平方根であることに注意してください。具体的には、2.9182 = sqrt(8.5160) および 0.32968 = sqrt(0.1087) です。相関は、切片と加速度の共分散および切片と加速度の標準偏差の関数です。切片と加速度の共分散は、psi の最初のセルの非対角値 –0.8387 です。したがって、相関は –.8387/(0.32968*2.92182) = –0.87 です。

切片と加速度のグループ化変数は Model_Year です。

stats の 2 番目のセルのインデックスを指定します。

stats{2}
ans = 


    Covariance Type: FullCholesky

    Group     Name1           Name2           Type         Estimate      Lower         Upper     
    Origin    'Weight'        'Weight'        'std'        0.00026075    9.2158e-05    0.00073775

stats の 2 番目のセルには、標準偏差の推定と、Weight の変量効果の項の標準偏差の 95% 信頼限界が含まれます。グループ化変数は Origin です。

stats の 3 番目のセルのインデックスを指定します。

stats{3}
ans = 

    Group    Name             Estimate    Lower     Upper
    Error    'Res Std'        3.0126      2.8028    3.238

stats の 3 番目のセルには、残差標準偏差の推定および 95% の信頼限界が含まれます。残差標準偏差の推定は、mse の平方根 sqrt(9.0755) = 3.0126 です。

共分散パラメーターに対する 99% 信頼区間を作成します。

[~,~,stats] = covarianceParameters(lme,'Alpha',0.01);
stats{1}
ans = 


    Covariance Type: FullCholesky

    Group         Name1                 Name2                 Type          Estimate    Lower       Upper   
    Model_Year    'Intercept'           'Intercept'           'std'           2.9182     0.86341      9.8633
    Model_Year    'Acceleration'        'Intercept'           'corr'        -0.87172    -0.99089    0.013164
    Model_Year    'Acceleration'        'Acceleration'        'std'          0.32968     0.15828     0.68669
stats{2}
ans = 


    Covariance Type: FullCholesky

    Group     Name1           Name2           Type         Estimate      Lower         Upper    
    Origin    'Weight'        'Weight'        'std'        0.00026075    6.6466e-05    0.0010229
stats{3}
ans = 

    Group    Name             Estimate    Lower    Upper 
    Error    'Res Std'        3.0126      2.74     3.3123
この情報は役に立ちましたか?