covarianceParameters
線形混合効果モデルの共分散パラメーターの抽出
構文
説明
例
切片に対する 2 つの変量効果
標本データを読み込みます。
load('fertilizer.mat');
このデータセット配列には土壌の種類に基づいて土壌が 3 つのブロックに分けられている分割プロット試験のデータが含まれています。土壌の種類は砂質、シルトおよび粘土質です。各ブロックは 5 つのプロットに分割され、5 種類のトマトの苗木 (チェリー、エアルーム、グレープ、枝付き、プラム) がランダムにこれらのプロットに割り当てられます。その後、プロット内のトマトの苗木はサブプロットに分割され、それぞれのサブプロットが 4 つの肥料の中の 1 つにより処置されます。このデータは、シミュレーションされたものです。
実用目的でこのデータを ds
という名前のデータセット配列に保存し、Tomato
、Soil
および Fertilizer
をカテゴリカル変数として定義します。
ds = fertilizer; ds.Tomato = nominal(ds.Tomato); ds.Soil = nominal(ds.Soil); ds.Fertilizer = nominal(ds.Fertilizer);
線形混合効果モデルを当てはめます。Fertilizer
は固定効果変数であり、平均収穫量はブロック (土壌の種類) とブロック内のプロット (土壌の種類の中のトマトの種類) によって独立して変化します。このモデルは以下に対応します。
ここで = 1、2、...、60 は観測値、 = 2、...、5 はトマトの種類に対応し、 = 1、2、3 はブロック (土壌) に対応します。 は 番目の土壌の種類を表し、 は 番目の土壌の種類で入れ子にされている 番目のトマトの種類を表します。 はトマトの種類 の水準を表すダミー変数です。
変量効果と観測誤差の事前分布は、、 および です。
lme = fitlme(ds,'Yield ~ Fertilizer + (1|Soil) + (1|Soil:Tomato)');
変量効果項の共分散パラメーターの推定値 ( と の推定値) を計算します。
psi = covarianceParameters(lme)
psi=2×1 cell array
{[2.7730e-17]}
{[ 352.8481]}
残差分散 () を計算します。
[~,mse] = covarianceParameters(lme)
mse = 151.9007
相関する可能性がある変量効果の項
標本データを読み込みます。
load('weight.mat');
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 を基準として使用し、必要なダミー変数 を作成します。このモデルは以下に対応します。
ここで、 は観測番号 に、 は被験者番号 に対応します。 は固定効果係数、、および 、 は変量効果です。 は初期重量を意味し、 はプログラムの種類を表すダミー変数です。たとえば、 はプログラム B を表すダミー変数です。
変量効果と観測値の誤差の事前分布は次のとおりです。
および
lme = fitlme(ds,'y ~ InitialWeight + Program + (Week|Subject)');
変量効果の共分散パラメーターの推定を計算します。
[psi,mse,stats] = covarianceParameters(lme)
psi = 1x1 cell array
{2x2 double}
mse = 0.0105
stats=2×1 cell array
{3x7 classreg.regr.lmeutils.titleddataset}
{1x5 classreg.regr.lmeutils.titleddataset}
mse
は推定された残差分散です。これは、 の推定値です。
変量効果項 (、 および ) に対する共分散パラメーターの推定値を確認するには、psi
のインデックスを指定します。
psi{1}
ans = 2×2
0.0572 0.0490
0.0490 0.0624
切片に対する変量効果項の分散の推定値 は 0.0572 です。週に対する変量効果項の分散の推定値 は 0.0624 です。切片と週に対する変量効果項の共分散の推定値 は 0.0490 です。
stats
は 2 行 1 列の cell 配列です。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
)、変量効果変数 (Name1
、Name2
)、共分散パラメーターの種類 (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) です。
2 つのグループ化変数
標本データを読み込みます。
load carbig
ガロンあたりの走行マイル数 (MPG) の線形混合効果モデルを当てはめます。加速度および重量は固定効果で、モデル年度によってグループ化される切片と加速度は相関された変量効果の可能性があり、重量は独立変量効果であり、自動車の原点でグループ化されます。このモデルは以下に対応します。
ここで、 は変数 Model_Year
の水準を表し、 は変数 Origin
の水準を表します。 は、 番目の観測値に対応する 番目のモデル年および 番目の生産国のガロンあたりの走行マイル数です。変量効果項と観測値の誤差の事前分布は次のとおりです。
ここで、変量効果項 は 1 番目のグループ化変数の水準 における 1 番目の変量効果を表します。変量効果項 は、1 番目のグループ化変数の 番目の水準 () における切片 (0) に対する 1 番目の変量効果項 (1) に対応します。同様に、 は 1 番目の変量効果項 (1) における 1 番目の予測子 (1) の水準 です。
同様に、 は 2 番目のグループ化変数の水準 における 2 番目の変量効果項を表します。
は切片に対する変量効果項の分散、 は予測子 Acceleration に対する変量効果項の分散、 は切片と予測子 Acceleration に対する変量効果項の共分散です。 は 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=2×1 cell array
{2x2 double }
{[6.6770e-08]}
mse = 9.0753
stats=3×1 cell array
{3x7 classreg.regr.lmeutils.titleddataset}
{1x7 classreg.regr.lmeutils.titleddataset}
{1x5 classreg.regr.lmeutils.titleddataset}
残差分散 mse
は 9.0755 です。psi
は 2 行 1 列の cell 配列であり、stats
は 3 行 1 列の cell 配列です。内容を確認するには、これらの cell 配列のインデックスを指定しなければなりません。
まず、psi
の最初のセルのインデックスを指定します。
psi{1}
ans = 2×2
8.2649 -0.8698
-0.8698 0.1157
psi
の 1 番目のセルには、切片に対する相関性がある変量効果の共分散パラメーター として 8.5160、加速度に対する として 0.1087 が格納されています。切片と加速度に対する変量効果項の共分散の推定値 は -0.8387 です。
次に、psi
の 2 番目のセルのインデックスを指定します。
psi{2}
ans = 6.6770e-08
psi
の 2 番目のセルには、重量に対する変量効果項の分散の推定値 が格納されています。
stats
の最初のセルのインデックスを指定します。
stats{1}
ans = COVARIANCE TYPE: FULLCHOLESKY Group Name1 Name2 Type Estimate Lower Upper Model_Year {'Intercept' } {'Intercept' } {'std' } 2.8749 1.0486 7.882 Model_Year {'Acceleration'} {'Intercept' } {'corr'} -0.88949 -0.98666 -0.3251 Model_Year {'Acceleration'} {'Acceleration'} {'std' } 0.34015 0.19355 0.59777
この表では、切片と加速度に関する変量効果の項の標準偏差が推定して示されています。標準偏差の推定は、psi
の最初のセルの対角要素の平方根であることに注意してください。具体的には、2.9182 = sqrt(8.5160) および 0.32968 = sqrt(0.1087) です。相関は、切片と加速度の共分散および切片と加速度の標準偏差の関数です。切片と加速度の共分散は、psi の 1 番目のセルの非対角値 -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.0002584 9.0882e-05 0.00073469
stats
の 2 番目のセルには、標準偏差の推定と、Weight
の変量効果の項の標準偏差の 95% 信頼限界が含まれます。グループ化変数は Origin
です。
stats
の 3 番目のセルのインデックスを指定します。
stats{3}
ans = Group Name Estimate Lower Upper Error {'Res Std'} 3.0125 2.8025 3.2383
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.8749 0.76378 10.821 Model_Year {'Acceleration'} {'Intercept' } {'corr'} -0.88949 -0.99322 0.0026856 Model_Year {'Acceleration'} {'Acceleration'} {'std' } 0.34015 0.16213 0.71364
stats{2}
ans = COVARIANCE TYPE: FULLCHOLESKY Group Name1 Name2 Type Estimate Lower Upper Origin {'Weight'} {'Weight'} {'std'} 0.0002584 6.5446e-05 0.0010202
stats{3}
ans = Group Name Estimate Lower Upper Error {'Res Std'} 3.0125 2.7395 3.3127
入力引数
lme
— 線形混合効果モデル
LinearMixedModel
オブジェクト
線形混合効果モデル。fitlme
または fitlmematrix
を使用して構築した LinearMixedModel
オブジェクトとして指定します。
名前と値の引数
オプションの引数のペアを Name1=Value1,...,NameN=ValueN
として指定します。ここで Name
は引数名、Value
は対応する値です。名前と値の引数は他の引数の後ろにする必要がありますが、ペアの順序は関係ありません。
R2021a より前では、名前と値をそれぞれコンマを使って区切り、Name
を引用符で囲みます。
例: [psi,mse,stats] = covarianceParameters(lme,'Alpha',0.01);
Alpha
— 有意水準
0.05 (既定値) | 0 ~ 1 の範囲のスカラー値
有意水準。'Alpha'
と 0 ~ 1 の範囲にあるスカラー値から構成されるコンマ区切りのペアとして指定します。値が α の場合、信頼水準は 100 × (1 – α)% です。
たとえば、99% の信頼区間の場合は、次のように信頼水準を指定できます。
例: 'Alpha',0.01
データ型: single
| double
出力引数
psi
— 共分散パラメーターの推定
cell 配列
変動効果の前の共分散をパラメーター表現する共分散パラメーターの推定。長さ R の cell 配列として返されます。psi{r}
には、グループ化変数 gr (r = 1, 2, ..., R) に関連付けられている変動効果の共分散行列が含まれます。グループ化変数の順序は、モデルを当てはめるときに入力する順序と同じです。
mse
— 残差分散の推定
スカラー値
残差分散の推定。スカラー値として返します。
stats
— 共分散パラメーターの推定および関連する統計
cell 配列
共分散パラメーターの推定および関連する統計。以下の列で構成されるデータセット配列を含む長さ (R + 1) の cell 配列として返します。
Group | グループ化変数名 |
Name1 | 最初の予測子変数の名前 |
Name2 | 2 番目の予測子変数の名前 |
Type |
|
Estimate |
|
Lower | 共分散パラメーターの 95% 信頼区間の下限 |
Upper | 共分散パラメーターの 95% 信頼区間の上限 |
stats{r}
は、r 番目のグループ化変数 (r = 1、2、...、R) の共分散パラメーターに対する統計を含むデータセット配列です。stats{R+1}
には、残差標準偏差の統計値が格納されます。残差誤差のデータセット配列には、Group
、Name
、Estimate
、Lower
、Upper
の各フィールドがあります。
バージョン履歴
R2013b で導入
MATLAB コマンド
次の MATLAB コマンドに対応するリンクがクリックされました。
コマンドを MATLAB コマンド ウィンドウに入力して実行してください。Web ブラウザーは MATLAB コマンドをサポートしていません。
Select a Web Site
Choose a web site to get translated content where available and see local events and offers. Based on your location, we recommend that you select: .
You can also select a web site from the following list:
How to Get Best Site Performance
Select the China site (in Chinese or English) for best site performance. Other MathWorks country sites are not optimized for visits from your location.
Americas
- América Latina (Español)
- Canada (English)
- United States (English)
Europe
- Belgium (English)
- Denmark (English)
- Deutschland (Deutsch)
- España (Español)
- Finland (English)
- France (Français)
- Ireland (English)
- Italia (Italiano)
- Luxembourg (English)
- Netherlands (English)
- Norway (English)
- Österreich (Deutsch)
- Portugal (English)
- Sweden (English)
- Switzerland
- United Kingdom (English)