ドキュメンテーション

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

線形回帰

データの準備

回帰の近似を開始するには、データを近似関数に望ましい形式にします。すべての回帰手法は、配列 X の入力データと独立したベクトル y の応答データか、テーブルまたはデータセット配列 tbl 内の入力データと tbl の列としての応答データで始まります。入力データの各行が、1 つの観測値を表します。各列が 1 つの予測子 (変数) を表します。

テーブルまたはデータセット配列 tbl では、次のように 'ResponseVar' の名前と値のペアで応答変数を示します。

mdl = fitlm(tbl,'ResponseVar','BloodPressure');
% or
mdl = fitglm(tbl,'ResponseVar','BloodPressure');

応答変数は既定で最後の列です。

数値の "カテゴリカル" 予測子を使用できます。カテゴリカル予測子は可能性のある固定セットからの値をとります。

  • 数値配列 X では、'Categorical' 名前と値のペアでカテゴリカル予測子を示します。たとえば、6 つの予測子から 23 がカテゴリカルであることを示すには、次のようにします。

    mdl = fitlm(X,y,'Categorical',[2,3]);
    % or
    mdl = fitglm(X,y,'Categorical',[2,3]);
    % or equivalently
    mdl = fitlm(X,y,'Categorical',logical([0 1 1 0 0 0]));
  • テーブルまたはデータセット配列 tbl では、近似関数はこれらのデータ型がカテゴリカルであることを想定しています。

    • 論理値

    • カテゴリカル (ノミナルまたは順序)

    • 文字列または文字配列

    数値予測値がカテゴリカルであることを示すには、'Categorical' 名前と値のペアを使用します。

欠損数値データは NaN で表されています。他のデータ型用の欠損データを表すには、グループ化変数の欠損値を参照してください。

入力と応答データのデータセット配列

Excel® スプレッドシートから数値配列を作成するには、次のようにします。

ds = dataset('XLSFile','hospital.xls',...
    'ReadObsNames',true);

ワークスペース変数からデータセット配列を作成するには、次のようにします。

load carsmall
ds = dataset(MPG,Weight);
ds.Year = ordinal(Model_Year);

入力および応答データのテーブル

Excel スプレッドシートからテーブルを作成するには、次のようにします。

tbl = readtable('hospital.xls',...
    'ReadRowNames',true);

ワークスペース変数からテーブルを作成するには、次のようにします。

load carsmall
tbl = table(MPG,Weight);
tbl.Year = ordinal(Model_Year);

入力データの数値配列、応答の数値ベクトル

たとえば、ワークスペース変数から数値配列を作成するには、次のようにします。

load carsmall
X = [Weight Horsepower Cylinders Model_Year];
y = MPG;

Excel スプレッドシートから数値配列を作成するには、次のようにします。

[X Xnames] = xlsread('hospital.xls');
y = X(:,4); % response y is systolic pressure
X(:,4) = []; % remove y from the X matrix

sex などの数値以外のエントリは X には表示されません。

近似手法の選択

モデルをデータにあてはめるには、以下の 3 つの方法があります。

最小二乗近似

fitlm を使用して、モデルのデータへの最小二乗近似を作成します。この方法は、モデルの型を十分に把握し、パラメーターを見つけることが主な目的である場合に最適です。いくつかのモデルを試す場合にも便利です。最小二乗近似を使う場合、外れ値を破棄するためにデータを手動で調べなければなりませんが、このタスクに役立つ方法もあります (残差 — 学習データのモデル品質を参照)。

ロバスト近似

外れ値にほとんど影響されないモデルを作成するには、fitlmRobustOpts 名前と値のペアと併用します。ロバスト近似は、外れ値を手動で破棄する手間が省けます。ただし、ロバスト近似では step が使えません。そのため、ロバスト近似を使用する場合、適正なモデルでステップワイズを検索できません。

ステップワイズ近似

モデルを検索し、パラメーターをそのモデルに適合させるには、stepwiselm を使用します。stepwiselm は、定数などの 1 つのモデルから開始し、一度に 1 つずつ項を増減する貪欲法で、それ以上改良できなくなるまで、最適な項を毎回選択します。ステップワイズ近似を使用して、有効な項のみを含む適正なモデルを見つけます。

この結果は、開始モデルに依存します。通常、定数モデルで開始すると、小さいモデルになります。さらに多くの項で開始すると、複雑なモデルになりますが、二乗平均誤差はより小さくなります。大小のステップワイズ モデルの比較を参照してください。

ロバスト オプションはステップワイズ近似では使用できません。ステップワイズ近似の実行後に、モデルで外れ値を探してください (残差 — 学習データのモデル品質を参照)。

単一モデルまたはモデル範囲の選択

線形回帰のモデルを指定する方法はいくつかあります。最も便利だと思う方法を使用してください。

fitlm で指定するのは近似のモデルです。モデル仕様を指定しない場合の既定値は 'linear' です。

stepwiselm では、指定するモデルの仕様は、ステップワイズ手順で改善を試みる最初のモデルです。モデルの仕様を指定しない場合の既定の最初のモデルは 'constant' で、既定の上限モデルは 'interactions' です。上限のモデルを変更するには、Upper 名前と値のペアを使用します。

簡単な文字列

文字列モデル タイプ
'constant'モデルは定数 (切片) 項だけを含みます。
'linear'モデルは各予測子に対して切片と線形項を含みます。
'interactions'モデルは、切片、線形項、異なる予測子のペアのすべての積 (二乗項なし) を含みます。
'purequadratic'モデルは、切片、線形項、二乗項を含みます。
'quadratic'モデルは、切片、線形項、交互作用、二乗項を含みます。
'polyijk'モデルは多項式であり、最初の予測子は次数 i まで、2 番目の予測子は次数 j まで、3 番目以降もすべて同様に 0 から 9 までの数値を使用します。たとえば、'poly2111' には、1 つの定数のほかにすべての線形項と積項があり、また、予測子 1 の二乗の項を含んでいます。

たとえば、行列予測子をもつ fitlm を使って交互作用モデルを指定するには、次のようにします。

mdl = fitlm(X,y,'interactions');

stepwiselm と予測子のテーブルまたはデータセット配列 tbl を使ってモデルを指定するには、1 つの定数から始め、線形モデル上限があると想定します。tbl 内の応答変数が 3 番目の列にあるとします。

mdl2 = stepwiselm(tbl,'constant',...
    'Upper','linear','ResponseVar',3);

項の行列

項行列は、モデル内の項を指定する t 行 (p + 1) 列の行列です。t は項の数、p は予測子変数の数、+ 1 は応答変数を示します。

T(i,j) の値は、項 i の変数 j の指数です。ABC の 3 つの予測子変数があると仮定します。

[0 0 0 0] % Constant term or intercept
[0 1 0 0] % B; equivalently, A^0 * B^1 * C^0
[1 0 1 0] % A*C
[2 0 0 0] % A^2
[0 1 2 0] % B*(C^2)
各項の最後の 0 は、応答変数を表します。一般に、以下のようになります。

  • テーブルまたはデータセット配列に変数がある場合、応答変数の位置に応じて 0 で応答変数を表さなければなりません。次の例はこれを説明します。

    標本データを読み込み、データセット配列を定義します。

    load hospital
    ds = dataset(hospital.Sex,hospital.BloodPressure(:,1),hospital.Age,...
    hospital.Smoker,'VarNames',{'Sex','BloodPressure','Age','Smoker'});

    線形モデル 'BloodPressure ~ 1 + Sex + Age + Smoker' を項行列で表します。応答変数はデータセット配列の 2 列目にあるので、応答変数を示す 0 の列が項行列の 2 列目になければなりません。

    T = [0 0 0 0;1 0 0 0;0 0 1 0;0 0 0 1]
    
    T =
    
         0     0     0     0
         1     0     0     0
         0     0     1     0
         0     0     0     1

    データセット配列を再定義します。

    ds = dataset(hospital.BloodPressure(:,1),hospital.Sex,hospital.Age,...
    hospital.Smoker,'VarNames',{'BloodPressure','Sex','Age','Smoker'});
    

    これで、応答変数はデータセット配列の最初の項になりました。項行列を使用して、同じ線形モデル 'BloodPressure ~ 1 + Sex + Age + Smoker' を指定します。

    T = [0 0 0 0;0 1 0 0;0 0 1 0;0 0 0 1]
    T =
    
         0     0     0     0
         0     1     0     0
         0     0     1     0
         0     0     0     1
  • 行列および列ベクトルに予測子と応答変数がある場合、各項の最後に応答変数を示す 0 がなければなりません。次の例はこれを説明します。

    標本データを読み込み、予測子の行列を定義します。

    load carsmall
    X = [Acceleration,Weight];
    

    項行列を使用してモデル 'MPG ~ Acceleration + Weight + Acceleration:Weight + Weight^2' を指定し、このモデルをデータにあてはめます。このモデルには、変数 Acceleration および Weight の主効果および 2 因子交互作用項、変数 Weight の 2 次の項が含まれます。

    T = [0 0 0;1 0 0;0 1 0;1 1 0;0 2 0]
    
    T =
    
         0     0     0
         1     0     0
         0     1     0
         1     1     0
         0     2     0
    

    線形モデルを近似します。

    mdl = fitlm(X,MPG,T)
    mdl = 
    
    Linear regression model:
        y ~ 1 + x1*x2 + x2^2
    
    Estimated Coefficients:
                       Estimate       SE            tStat      pValue    
        (Intercept)         48.906        12.589     3.8847    0.00019665
        x1                 0.54418       0.57125    0.95261       0.34337
        x2               -0.012781     0.0060312    -2.1192      0.036857
        x1:x2          -0.00010892    0.00017925    -0.6076         0.545
        x2^2            9.7518e-07    7.5389e-07     1.2935       0.19917
    
    Number of observations: 94, Error degrees of freedom: 89
    Root Mean Squared Error: 4.1
    R-squared: 0.751,  Adjusted R-Squared 0.739
    F-statistic vs. constant model: 67, p-value = 4.99e-26

    変数 Weight に対応する切片および x2 の項のみが 5% の有意水準を示しています。

    ここで、定数モデルを開始モデルとし、交互作用項をもつ線形モデルを上位モデルとして、ステップワイズ回帰を実行します。

    T = [0 0 0;1 0 0;0 1 0;1 1 0];
    mdl = stepwiselm(X,MPG,[0 0 0],'upper',T)
    1. Adding x2, FStat = 259.3087, pValue = 1.643351e-28
    
    mdl = 
    
    Linear regression model:
        y ~ 1 + x2
    
    Estimated Coefficients:
                       Estimate      SE           tStat      pValue    
        (Intercept)        49.238       1.6411     30.002    2.7015e-49
        x2             -0.0086119    0.0005348    -16.103    1.6434e-28
    
    Number of observations: 94, Error degrees of freedom: 92
    Root Mean Squared Error: 4.13
    R-squared: 0.738,  Adjusted R-Squared 0.735
    F-statistic vs. constant model: 259, p-value = 1.64e-28

    ステップワイズ回帰の結果は、前の手順の fitlm の結果と整合性があります。

モデル仕様の式は次のような形式の文字列です

'Y ~ terms'

  • Y は応答名です。

  • terms 次が含まれます

    • 変数名

    • + 次の変数を含みます

    • - 次の変数を除外します

    • : 項の積である交互作用を定義します

    • * 交互作用とすべての次数の低い項を定義します

    • * で繰り返されるとおり、^ は予測子をべき乗にし、^ は低い次数の項も含みます

    • () 項をグループ化します

    ヒント:   式には既定で定数 (切片) 項が含まれます。モデルから定数項を除外するには、式に -1 を含めます。

次に例を示します。

'Y ~ A + B + C' は切片のある 3 変数線形モデルです。
'Y ~ A + B + C - 1' は切片のない 3 変数線形モデルです。
'Y ~ A + B + C + B^2' は切片と B^2 項のある 3 変数モデルです。
'Y ~ A + B^2 + C' は、B^2B 項を含むため、前の例と同じです。
'Y ~ A + B + C + A:B'A*B 項を含みます。
'Y ~ A*B + C'A*B = A + B + A:B であるため前の例と同じです。
'Y ~ A*B*C - A:B:C' には 3 因子交互作用を除き、ABC 間のすべての交互作用があります。
'Y ~ A*(B + C + D)' にはすべての線形項があり、さらに他の変数のそれぞれとの A の積があります。

たとえば、行列予測子をもつ fitlm を使って交互作用モデルを指定するには、次のようにします。

mdl = fitlm(X,y,'y ~ x1*x2*x3 - x1:x2:x3');

stepwiselm と予測子のテーブルまたはデータセット配列 tbl を使ってモデルを指定するには、1 つの定数から始め、線形モデル上限があると想定します。tbl の応答変数が 'y'、予測子変数が 'x1''x2' および'x3' という名前であるとします。

mdl2 = stepwiselm(tbl,'y ~ 1','Upper','y ~ x1 + x2 + x3');

モデルのデータへのあてはめ

最も一般的な近似のオプション引数は次のとおりです。

  • fitlm のロバスト回帰で 'RobustOpts' 名前と値のペアを 'on' にします。

  • stepwiselm で適切な上限モデルを指定します。たとえば、'Upper''linear' に設定します。

  • 'CategoricalVars' 名前と値のペアを使用してどの変数がカテゴリカルであるか指定します。予測子 1 および 6 がカテゴリカル変数であることを指定するには、[1 6] などの列番号をもつベクトルを入力します。または、変数がカテゴリカルであることを示す 1 のエントリをもつ、データ列と同じ長さの論理ベクトルを入力します。7 つの予測子があり、予測子 16 がカテゴリカルである場合、logical([1,0,0,0,0,1,0]) と指定します。

  • テーブルまたはデータセット配列では 'ResponseVar' 名前と値のペアで応答変数を示します。既定の設定は、配列の最後の列です。

次に、例を示します。

mdl = fitlm(X,y,'linear',...
    'RobustOpts','on','CategoricalVars',3);
mdl2 = stepwiselm(tbl,'constant',...
    'ResponseVar','MPG','Upper','quadratic');

品質の調査と近似モデルの調整

モデルで近似した後に、結果を確認します。

モデル表示

線形回帰モデルは、その名前を入力するとき、または disp(mdl). を入力するときにいくつかの診断を表示します。この表示は、近似モデルがデータを適切にあらわすかどうかを確認するために基本情報のいくつかを提供します。

たとえば、線形モデルを 5 つの予測子のうち 2 つが存在せず、切片の項がないデータをモデルにあてはめるには、次のようにします。

X = randn(100,5);
y = X*[1;0;3;0;-1]+randn(100,1);
mdl = fitlm(X,y)
mdl = 

Linear regression model:
    y ~ 1 + x1 + x2 + x3 + x4 + x5

Estimated Coefficients:
                   Estimate     SE          tStat       pValue    
    (Intercept)     0.038164    0.099458     0.38372       0.70205
    x1               0.92794    0.087307      10.628    8.5494e-18
    x2             -0.075593     0.10044    -0.75264       0.45355
    x3                2.8965    0.099879          29    1.1117e-48
    x4              0.045311     0.10832     0.41831       0.67667
    x5              -0.99708     0.11799     -8.4504     3.593e-13

Number of observations: 100, Error degrees of freedom: 94
Root Mean Squared Error: 0.972
R-squared: 0.93,  Adjusted R-Squared 0.926
F-statistic vs. constant model: 248, p-value = 1.5e-52

次の点に注意してください。

  • 表示には Estimate 列に各係数の推定値が含まれます。これらは真の値 [0;1;0;3;0;-1] にかなり近い値です。

  • 係数推定には標準誤差列があります。

  • 予測子 1、3 および 5 に対してレポートされた pValue (標準誤差と仮定して t 統計から導出された) は極端に小さい値です。これらは、応答データ y を作成するために使われた 3 つの予測子です。

  • (Intercept)x2 および x4 に対する pValue は 0.01 よりはるかに大きい値です。これらの 3 つの予測子は、応答データ y を作成するためには使われませんでした。

  • 表示には R2、自由度調整済み決定係数 (R2) および F 統計が含まれます。

分散分析

近似されたモデルの品質を確認するには、ANOVA 表を調べます。たとえば、5 つの予測子をもつ線形モデルで anova を次のように使用します。

X = randn(100,5);
y = X*[1;0;3;0;-1]+randn(100,1);
mdl = fitlm(X,y);
tbl = anova(mdl)
tbl = 
             SumSq      DF    MeanSq     F          pValue    
    x1        106.62     1     106.62     112.96    8.5494e-18
    x2       0.53464     1    0.53464    0.56646       0.45355
    x3        793.74     1     793.74     840.98    1.1117e-48
    x4       0.16515     1    0.16515    0.17498       0.67667
    x5        67.398     1     67.398      71.41     3.593e-13
    Error     88.719    94    0.94382

この表には既定の表示以外の結果も含まれます (モデル表示を参照)。x2x4 の効果が重要でないことが、この表で明確にわかります。目的によってはこのモデルから x2x4 を削除することを検討してください。

診断プロット

診断プロットによって外れ値を特定でき、モデルや近似で他の問題を確認できます。たとえば、carsmall データを読み込み、Cylinders (ノミナル) と Weight の関数として MPG のモデルを作成します。

load carsmall
tbl = table(Weight,MPG,Cylinders);
tbl.Cylinders = ordinal(tbl.Cylinders);
mdl = fitlm(tbl,'MPG ~ Cylinders*Weight + Weight^2');ˋ

データとモデルのてこ比のプロットを作成します。

plotDiagnostics(mdl)

てこ比の高い点がいくつかあります。ただし、てこ比の高い点が外れ値かどうかはこのプロットではわかりません。

クックの距離が大きい点を探します。

plotDiagnostics(mdl,'cookd')

クックの距離が大きい点が 1 つあります。この点を特定してモデルから削除します。データ カーソルを使用して、外れ値をクリックして特定するか、プログラムを使用して次のように外れ値を特定します。

[~,larg] = max(mdl.Diagnostics.CooksDistance);
mdl2 = fitlm(tbl,'MPG ~ Cylinders*Weight + Weight^2',...
    'Exclude',larg);

残差 — 学習データのモデル品質

モデルまたはデータ内の誤差、外れ値または相関を検出できるいくつかの残差プロットが存在します。最もシンプルな残差プロットは既定のヒストグラム プロットであり、これは残差の範囲と頻度を示します。また、確率プロットは残差の分布が正規分布と一致する分散を比較する方法を示します。

carsmall データを読み込み、Cylinders (ノミナル) と Weight の関数として MPG のモデルを作成します。

load carsmall
tbl = table(Weight,MPG,Cylinders);
tbl.Cylinders = ordinal(tbl.Cylinders);
mdl = fitlm(tbl,'MPG ~ Cylinders*Weight + Weight^2');

残差の検査

plotResiduals(mdl)

12 を超える観測は外れ値の可能性があります。

plotResiduals(mdl,'probability')

外れ値の可能性がある 2 つの値もこのプロットに表示されます。そうでない場合は、確率プロットはかなり直線に見え、これは正規分布している残差に適切に近似していることを意味します。

次のように入力すると、2 つの外れ値を特定してデータから削除することができます。

outl = find(mdl.Residuals.Raw > 12)
outl =

    90
    97

外れ値を削除するには、次の Exclude 名前と値のペアを使用します。

mdl2 = fitlm(tbl,'MPG ~ Cylinders*Weight + Weight^2',...
    'Exclude',outl);

mdl2 の残差プロットを調べます。

plotResiduals(mdl2)

新しい残差プロットはかなり左右対称で、目立った問題はないように見えます。しかし、残差の間にいくらか系列相関が見られます。新しいプロットを作成して、そのような効果が存在するかどうかを確認します。

plotResiduals(mdl2,'lagged')

散布図は、第 1 象限と第 3 象限に他の 2 つの象限より多くのプロットが表示され、残差の間に正の系列相関があることを示しています。

可能性のある別の問題は、広範囲の観測で残差が大きい場合です。現在のモデルにこの問題がないか確認してください。

plotResiduals(mdl2,'fitted')

近似値が比較的大きい場合に残差が大きくなる傾向があります。モデルの誤差はおそらく測定値に比例します。

予測子の効果を理解するためのプロット

次の例では、さまざまなプロットを使用して各予測子が回帰モデルに与える影響を理解する方法について説明します。

  1. carsmall データ内のいくつかの予測子から燃費のモデルを作成します。

    load carsmall
    tbl = table(Weight,MPG,Cylinders);
    tbl.Cylinders = ordinal(tbl.Cylinders);
    mdl = fitlm(tbl,'MPG ~ Cylinders*Weight + Weight^2');
  2. 応答のスライス プロットを調査します。これは、各予測子の効果を個別に表示します。

    plotSlice(mdl)

    青い縦の破線によって表される個々の予測子値をドラッグできます。また、赤い破線の曲線で表されている、同時信頼限界と非同時信頼限界を選択することもできます。

  3. 効果のプロットを使用して、応答に対する予測子の効果を別に画面に表示します。

    plotEffects(mdl)

    このプロットで、Weight が約 2,500 から 4,732 に変わると MPG が約 30 低下する (上の青い円の位置) ことがわかります。また、気筒の数が 8 から 4 に変化すると、MPG は約 10 増加します (下の青い円の位置)。青い水平ラインはこれらの予測子の信頼区間を表しています。この予測は、別の予測子が変化したときの 1 つの予測子の平均をもとにしています。このように 2 つの予測子が相関関係にある場合、結果の解釈を慎重に行ってください。

  4. 別の予測子が変化したときの 1 つの予測子の平均の効果を表示するのではなく、相互作用プロットの結合部の相互作用を調べます。

    plotInteraction(mdl,'Weight','Cylinders')

    相互作用プロットは、別の予測子を固定して、1 つの予測子を変化させたときの効果を表します。このプロットのほうがはるかに多くの情報を提供します。たとえば、比較的軽量の自動車 (Weight = 1,795) では、気筒の数が減少すると燃費は向上しますが、比較的重い自動車 (Weight = 4,732) の場合は、気筒の数が減少すると燃費も低下します。

  5. 相互作用についてさらに詳しく理解するには、予測子を含む相互作用プロットを調べます。1 つの予測子を固定し、別の予測子を変化させたこのプロットでは、その効果が曲線で示されます。気筒をさまざまな数に固定して、相互作用を調べます。

    plotInteraction(mdl,'Cylinders','Weight','predictions')

    次に、重量のさまざまな固定レベルで相互作用を調べます。

    plotInteraction(mdl,'Weight','Cylinders','predictions')

項の効果を理解するためのプロット

次の例では、回帰モデルにおける各項の影響を理解する方法について、さまざまなプロットを使用して説明します。

  1. carsmall データ内のいくつかの予測子から燃費のモデルを作成します。

    load carsmall
    tbl = table(Weight,MPG,Cylinders);
    tbl.Cylinders = ordinal(tbl.Cylinders);
    mdl = fitlm(tbl,'MPG ~ Cylinders*Weight + Weight^2');
  2. 追加変数としての Weight^2 で、追加変数プロットを作成します。

    plotAdded(mdl,'Weight^2')

    このプロットは、Weight^2MPG の両方を Weight^2 以外の項に近似した結果を示します。plotAdded を使用する理由は、Weight^2 を追加することによって、どのような追加の改善がモデルに期待できるかを理解するためです。これらの点を近似する線の係数は、フル モデルの Weight^2 の係数です。係数表の表示から理解できるように、Weight^2 予測子は有意性 (pValue < 0.05) の境目を超えたあたりにあります。これはプロットでも確認できます。信頼限界に水平ライン (定数 y) を含むことができないと思われるため、ゼロ勾配モデルはこのデータと整合しません。

  3. モデルの追加変数プロットを全体として作成します。

    plotAdded(mdl)

    モデル全体は非常に有意性が高いため、信頼限界は水平ラインを包含するには至りません。そのラインの勾配は、最良近似の方向に投影された予測子への近似の勾配、すなわち、係数ベクトルのノルムです。

モデルの変更

モデルを変更するには、次の 2 つの方法があります。

  • step — 一度に 1 つずつ項を加算または減算します。ここで step は追加または削除する最も重要な項を選択します。

  • addTerms および removeTerms — 指定された項を追加または削除します。単一モデルまたはモデル範囲の選択 で説明したいずれかの型を使用して、項を指定します。

stepwiselm を使用してモデルを作成した場合、 異なる上限モデルまたは下限モデルを指定したときにのみ step の効果があります。RobustOpts を使用してモデルを近似した場合は、step は機能しません。

たとえば、carbig データをもとに燃費の線形モデルから開始します。

load carbig
tbl = table(Acceleration,Displacement,Horsepower,Weight,MPG);
mdl = fitlm(tbl,'linear','ResponseVar','MPG')
mdl = 


Linear regression model:
    MPG ~ 1 + Acceleration + Displacement + Horsepower + Weight

Estimated Coefficients:
                    Estimate      SE            tStat       pValue    
    (Intercept)         45.251         2.456      18.424    7.0721e-55
    Acceleration     -0.023148        0.1256     -0.1843       0.85388
    Displacement    -0.0060009     0.0067093    -0.89441       0.37166
    Horsepower       -0.043608      0.016573     -2.6312      0.008849
    Weight          -0.0052805    0.00081085     -6.5123    2.3025e-10


Number of observations: 392, Error degrees of freedom: 387
Root Mean Squared Error: 4.25
R-squared: 0.707,  Adjusted R-Squared 0.704
F-statistic vs. constant model: 233, p-value = 9.63e-102

step を使用してモデルの改良を最大 10 ステップ試みます。

mdl1 = step(mdl,'NSteps',10)
1. Adding Displacement:Horsepower, FStat = 87.4802, pValue = 7.05273e-19

mdl1 = 


Linear regression model:
    MPG ~ 1 + Acceleration + Weight + Displacement*Horsepower

Estimated Coefficients:
                               Estimate      SE            tStat      pValue    
    (Intercept)                    61.285        2.8052     21.847    1.8593e-69
    Acceleration                 -0.34401       0.11862       -2.9     0.0039445
    Displacement                -0.081198      0.010071    -8.0623    9.5014e-15
    Horsepower                   -0.24313      0.026068    -9.3265    8.6556e-19
    Weight                     -0.0014367    0.00084041    -1.7095      0.088166
    Displacement:Horsepower    0.00054236    5.7987e-05     9.3531    7.0527e-19


Number of observations: 392, Error degrees of freedom: 386
Root Mean Squared Error: 3.84
R-squared: 0.761,  Adjusted R-Squared 0.758
F-statistic vs. constant model: 246, p-value = 1.32e-117

1 回のみの変更後に step が停止しました。

モデルの簡略化を試みるには、Acceleration および Weight の項を mdl1 から削除します。

mdl2 = removeTerms(mdl1,'Acceleration + Weight')
mdl2 = 


Linear regression model:
    MPG ~ 1 + Displacement*Horsepower

Estimated Coefficients:
                               Estimate      SE           tStat      pValue     
    (Intercept)                    53.051        1.526     34.765    3.0201e-121
    Displacement                -0.098046    0.0066817    -14.674     4.3203e-39
    Horsepower                   -0.23434     0.019593     -11.96     2.8024e-28
    Displacement:Horsepower    0.00058278    5.193e-05     11.222     1.6816e-25


Number of observations: 392, Error degrees of freedom: 388
Root Mean Squared Error: 3.94
R-squared: 0.747,  Adjusted R-Squared 0.745
F-statistic vs. constant model: 381, p-value = 3e-115

mdl2DisplacementHorsepower のみを使用し、Adjusted R-Squared メトリックの mdl1 とほぼ同じくらい良好なデータの近似を得ます。

新しいデータに対する応答を予測またはシミュレート

新しいデータに対する応答を予測またはシミュレートするために、3 つの手法を使用できます。

predict

この例では、predict メソッドを使用して予測で信頼区間を予測して取得する方法を示します。

  1. carbig データを読み込んで、応答 MPG の既定の線形モデルを AccelerationDisplacementHorsepower および Weight の予測子にします。

    load carbig
    X = [Acceleration,Displacement,Horsepower,Weight];
    mdl = fitlm(X,MPG);
  2. 最小値、平均値、最大値から予測子の 3 行の配列を作成します。いくつかの NaN 値が存在するため、NaN 値を無視する関数を使用します。

    Xnew = [nanmin(X);nanmean(X);nanmax(X)]; % new data
  3. 予測モデル応答と予測子の信頼区間を求めます。

    [NewMPG NewMPGCI] = predict(mdl,Xnew)
    NewMPG =
       34.1345
       23.4078
        4.7751
    
    NewMPGCI =
       31.6115   36.6575
       22.9859   23.8298
        0.6134    8.9367

    平均応答の信頼限界は最大応答または最小応答の信頼限界よりも範囲が狭くなりますが、これは妥当な結果です。

feval

テーブルまたはデータセット配列からモデルを構築するときに、多くの場合、feval のほうが predict よりも平均応答を予測するのに便利です。ただし、feval は信頼限界を提供しません。

この例は feval メソッドを使用して平均応答を予測する方法を示しています。

  1. carbig データを読み込んで、応答 MPG の既定の線形モデルを AccelerationDisplacementHorsepower および Weight の予測子にします。

    load carbig
    tbl = table(Acceleration,Displacement,Horsepower,Weight,MPG);
    mdl = fitlm(tbl,'linear','ResponseVar','MPG');
  2. 最小値、平均値、最大値から予測子の 3 行の配列を作成します。いくつかの NaN 値が存在するため、NaN 値を無視する関数を使用します。

    X = [Acceleration,Displacement,Horsepower,Weight];
    Xnew = [nanmin(X);nanmean(X);nanmax(X)]; % new data

    配列 Xnew には予測に適した列数があるため、feval はこれを使用して予測を実行できます。

  3. 予測モデル応答を求めます。

    NewMPG = feval(mdl,Xnew)
    NewMPG =
       34.1345
       23.4078
        4.7751

random

random メソッドは新しいランダム応答値をシミュレートします。これは、平均予測に学習データと同じ分散値のランダム外乱を加えたものと同じです。

この例は random メソッドを使用して応答をシミュレートする方法を示しています。

  1. carbig データを読み込んで、応答 MPG の既定の線形モデルを AccelerationDisplacementHorsepower および Weight の予測子にします。

    load carbig
    X = [Acceleration,Displacement,Horsepower,Weight];
    mdl = fitlm(X,MPG);
  2. 最小値、平均値、最大値から予測子の 3 行の配列を作成します。いくつかの NaN 値が存在するため、NaN 値を無視する関数を使用します。

    Xnew = [nanmin(X);nanmean(X);nanmax(X)]; % new data
  3. ランダム性を含む新しい予測モデル応答を生成します。

    rng('default') % for reproducibility
    NewMPG = random(mdl,Xnew)
    NewMPG =
    
       36.4178
       31.1958
       -4.8176
  4. 負の値の MPG は妥当ではないと思われるため、予測を 2 回試行してください。

    NewMPG = random(mdl,Xnew)
    NewMPG =
    
       37.7959
       24.7615
       -0.7783
    NewMPG = random(mdl,Xnew)
    NewMPG =
    
       32.2931
       24.8628
       19.9715

    Xnew の 3 行目の予測 (maximal) は明らかに信頼できません。

近似モデルの共有

以下のコマンドを使用して mdl のような線形回帰モデルが生成されたとします。

load carbig
tbl = table(Acceleration,Displacement,Horsepower,Weight,MPG);
mdl = fitlm(tbl,'linear','ResponseVar','MPG');

このモデルを他のユーザーと共有するには、以下の手順が可能です。

  • モデルの表示を提供する。

    mdl
    mdl = 
    
    
    Linear regression model:
        MPG ~ 1 + Acceleration + Displacement + Horsepower + Weight
    
    Estimated Coefficients:
                        Estimate      SE            tStat       pValue    
        (Intercept)         45.251         2.456      18.424    7.0721e-55
        Acceleration     -0.023148        0.1256     -0.1843       0.85388
        Displacement    -0.0060009     0.0067093    -0.89441       0.37166
        Horsepower       -0.043608      0.016573     -2.6312      0.008849
        Weight          -0.0052805    0.00081085     -6.5123    2.3025e-10
    
    
    Number of observations: 392, Error degrees of freedom: 387
    Root Mean Squared Error: 4.25
    R-squared: 0.707,  Adjusted R-Squared 0.704
    F-statistic vs. constant model: 233, p-value = 9.63e-102
  • モデル定義と係数のみを提供する。

    mdl.CoefficientNames
    ans = 
    
        '(Intercept)'    'Acceleration'    'Displacement'    'Horsepower'    'Weight'
    mdl.Coefficients.Estimate
    ans =
    
       45.2511
       -0.0231
       -0.0060
       -0.0436
       -0.0053
    mdl.Formula
    ans = 
    
    MPG ~ 1 + Acceleration + Displacement + Horsepower + Weight
この情報は役に立ちましたか?