アルゴリズム(多項式回帰)



多項式モデル

多項式モデル

与えられたデータセット (x_i , y_i), i = 1,2,...n, に対して、( X は独立変数、Y は従属変数) 次の形式のモデルを使ってデータに多項式回帰を実行します。

 y_i=\beta _0+\beta _1x_i+\beta_2x_i^2+\beta_3x_i^3+...+\beta_kx_i^k+\varepsilon_i

(1)

ここで、 k は多項式の次数です。Originでは、kは10より小さい正の数です。 パラメータは重み付けした最小二乗法を使って推定されます。 この方法は、独立変数の範囲内で理論曲線と測定データポイント間の差の二乗和を最小化します。フィットの後、仮説検定と残差のプロットを使ってモデルを評価できます。

\beta _0\,\! はy切片で、パラメータ \beta _1\,\!, \beta _2\,\!, ..., \beta _k\,\! は、部分係数(部分勾配)です。 行列形式で記述することができます。

Y=XB+E\,\!

(2)

および、

Y=\begin{bmatrix}
y_1\\
y_2\\
\vdots \\
y_n
\end{bmatrix}
, X=\begin{bmatrix}
1 & x_{11} & x_{1}^2 & \cdots & x_{1}^k\\
1 & x_{21} & x_{2}^2 & \cdots & x_{2}^k\\
\vdots & \vdots  & \vdots & \ddots & \vdots\\
1 & x_{n1} & x_{n}^2 & \cdots & x_{n}^k
\end{bmatrix}  
B=\begin{bmatrix}
\beta_0 \\
\beta_1 \\
\vdots \\
\beta_k
\end{bmatrix} , E=\begin{bmatrix}
\varepsilon _1\\
\varepsilon _2\\
\vdots \\
\varepsilon _n
\end{bmatrix}

\varepsilon_iは、\bar{E}=0かつVar[E]=\sigma^2における正規確率変数として、独立かつ一様に分布していると仮定します。 Bについて、\left \|E\right \|を最小にするには、次の関数を使います。

\frac{\partial E'E}{\partial B}=0

(3)

結果\hat Bは、ベクターデータ B最小二乗推定値値は、線形方程式の解で、次の様に表すことができます。

\hat B=\begin{bmatrix}
\hat \beta_0 \\
\hat \beta_1 \\
\vdots \\
\hat \beta_k
\end{bmatrix}=(X'X)^{-1}X^{\prime }Y

(4)

ここで、X'は、Xの転置で、与えられたXに対するYの推定値は次のようになります。

\hat{Y}=X\hat{B}

(5)

(4)に\hat{B}を置換して、行列Pを定義することができます。

\hat{Y}=[X(X'X)^{-1}X']Y=PY

(6)

残差は次のように定義されます。

res_i=Y-\hat{Y}

(7)

残差平方和は、次の通りです。

RSS=\left \|E \right \|^2={Y}'Y-\hat{B}'X'X\hat{B}

(8)

Note:より高い次数の多項式が従属変数に対して最も効果的であることに注目すべきです。そのため、高次の項(4を超える)を持つモデルは、係数値の精度に非常に敏感であり、係数値のわずかな違いが、計算されたy値に大きな違いをもたらす可能性があります。これに言及するのは、デフォルトでは、多項式フィッティングの結果が小数点以下5桁に丸められるためです。これらの出力されたワークシートの値を手動でフィット曲線に戻すと、丸めで発生するわずかな精度の低下が、高次の項に大きい影響を及ぼし、モデルが間違っていると誤った結論を付ける可能性があります。 最適なパラメータ推定値を使用して手動計算を実行する場合は、丸められた値ではなく、完全精度の値を使用してください。 Originが取得した値を小数点以下5桁(またはその他)に丸める場合、これらの値は表示上のみです。特に指定がない場合は、常に、数値計算では完全精度((double(8))を使用します。 詳細は、Originヘルプファイルの 数値の扱いについてをご覧ください。

一般的に、任意の連続関数を高次の多項式モデルにフィットさせることができます。 ただし、より高次の項にはあまり実用的な意味はありません。

フィット制御

誤差を重みとする

上記のセクションでは、誤差に定数分散があると仮定しています。しかし、実験値をフィッティングする場合、(計測器の確度と精度に影響を及ぼす、)機器誤差を考慮する必要があります。したがって、誤差の定数分散推定は、棄却されます。\varepsilon_iを非定数分散の正規分布であると推定する必要があります。また、誤差は、\sigma^2のようになり、フィッティングで重みとして使用することができます。重みは、次のように定義されます。

W=\begin{bmatrix}
 w_1& 0 & \dots &0 \\ 
0 & w_2 & \dots &0 \\ 
 \vdots& \vdots &\ \ddots &\vdots \\ 
 0& 0 &\dots  & w_n
\end{bmatrix}

フィッティングモデルは、次の式になります。

\sum_{i=1}^n w_i (y_i-\hat y_i)^2=\sum_{i=1}^n w_i [y_i-(\hat{\beta _0}+\hat{\beta _1}x_i)]^2

(9)

重み因子w_i は、3つの式によって与えられます。

重み付けなし

エラーバーは、計算では重みとして取り扱われません。

直接重み付け

w_i=\sigma_i

(10)

機械的

機械的重みとして、値は、機械的誤差に反比例します。大きな誤差がある場合よりも正確であるため、小さな誤差の試行には、大きな重みがあります。

w_i=\frac 1{\sigma_i^2}

(11)

重みとしての誤差は、ワークシートの「YError」として構築されています。

切片固定

固定切片は、y切片\beta_0を設定して、値を固定します。また、 固定切片のため、全ての自由度は、n*=n-1となります。

sqrt(補正カイ二乗値)のスケールエラー

sqrt(補正カイ二乗値)のスケールエラーは、重みを付けたフィットで、使用することができます。このオプションは、フィット処理で出力されるパラメータの誤差だけに影響し、フィット処理やデータには影響しません。 デフォルトで、これにはチェックが入っており、Eの分散\sigma^2は、 パラメータ誤差の計算を考慮しているか、 あるいは、分散Eは、誤差計算を考慮していません。 共分散行列を例に挙げます。

sqrt(補正カイ二乗値)のスケールエラー

Cov(\beta _i,\beta _j)=\sigma^2 (X^{\prime }X)^{-1}

(12)

\sigma^2=\frac{RSS}{n^{*}-k}

sqrt(補正カイ二乗値)のスケールエラーでは無い

Cov(\beta _i,\beta _j)=(X'X)^{-1}\,\!

(13)

重み付けフィットには、(X'X)^{-1}\,\!の代わりに、(X'WX)^{-1}\,\!を使います。

フィット結果

フィットパラメータ

Mfitted-paramater.png

フィット値

式(4)

パラメータの標準誤差

各パラメータにおいて、標準誤差は以下のように得られます。

s_{\hat \beta _j}=s_\varepsilon \sqrt{C_{jj}}

(14)

ここで、C_{jj}は、j番目の(X'X)^{-1}の対角要素です。((X'WX)^{-1}は、重み付けフィットに使われます。)s_\varepsilon は、次式で計算される残差標準偏差です。 (「std dev」、「推定の標準誤差」、「root MSE」のようにも呼びます。)

s_\varepsilon =\sqrt{\frac{RSS}{df_{Error}}}=\sqrt{\frac{RSS}{n^{*}-k}}

(15)

s_\varepsilonは、 Eの分散である\sigma ^2の推定です。

Note:自由度(df)dfErrorについての詳細は、ANOVA 表をご覧ください。

t値と信頼水準

回帰の仮定が成り立つ場合、帰無仮説と対立仮説を使用して回帰係数のt検定を実行できます。

H_0:\beta _j=0\,\!

H_\alpha :\beta _j\neq 0

t値は、次の式で計算できます。

t=\frac{\hat \beta _j-0}{s_{\hat \beta _j}}

(16)

計算されたt値を使って、対応する帰無仮説を棄却するかどうかを決めることができます。通常、与えられたパラメータの信頼水準\alpha\,\! について、|t|>t_{\frac \alpha 2} のときは、H_0\,\! を棄却できます。さらに、 p-値は \alpha\,\! より小さくなります。

Prob>|t|

t 検定の H_0 が真である確率

prob=2(1-tcdf(|t|,df_{Error}))\,\!

(17)

ここで、tcdf(|t|,df_{Error})は、 |t|値におけるスチューデントt 分布の累積分布関数を、誤差の自由度(df_{Error})で計算します。

LCL and UCL

t値から各パラメータの (1-\alpha )\times 100\% 信頼区間を次式で計算することができます。

\hat \beta _j-t_{(\frac \alpha 2,n^{*}-k)}\varepsilon _{\hat \beta _j}\leq \hat \beta _j\leq \hat \beta _j+t_{(\frac \alpha 2,n^{*}-k)}\varepsilon _{\hat \beta _j}

(18)

ここでUCLLCLは、それぞれ上側信頼区間下側信頼区間のことです。

CI 半幅

信頼区間の半値幅は以下の通りです。

CI=\frac{UCL-LCL}2

(19)

フィット統計

いくつかのフィット統計式をここに要約します。

Mfitted-stats.png

自由度

(誤差)の変数に対する自由度。詳細は ANOVA表を参照してください。

自由度あたりカイ二乗

\sigma^2=\frac{RSS}{n^{*}-k}

(20)

残差平方和

残差平方和。式(8)を参照。

R二乗(COD)

フィットの良さは、 決定係数(COD) R^2 で評価でき、次の式になります。

R^2=\frac{Explained\, variation}{Total \, variation}=1-\frac{RSS}{TSS}

(21)

補正R二乗

補正 R^2 は、自由度の R^2 値を調整するのに使用されます。これは次式のように計算されます。

\bar R^2=1-\frac{RSS/df_{Error}}{TSS/df_{Total}}

(22)

R値

相関係数 R値は、R^2 の平方根を使って計算できます。

R=\sqrt{R^2}

(23)

Root MSE(SD)

誤差の平均平方の平方根または、残差標準偏差は、次式に等しくなります。

RootMSE=\sqrt{\frac{RSS}{n^*-k}}

(24)

残差のノルム

RSSの平方根に等しい。

Norm \,of \,Residuals=\sqrt{RSS}

(25)

ANOVA表

多項式フィットのANOVA表は

df 平方和 平均平方 F値 Prob > F
モデル k SS_{reg} = TSS-RSS MS_{reg} = SS_{reg} / k MS_{reg} / MSE p-値
誤差 n* - k RSS MSE = RSS / (n^* - k)
合計 n* TSS
Note: 切片がモデルに含まれてる場合、 n*=n-1 です。それ以外は、 n*=n で平方和の合計は未補正となります。

ここで、平方和の合計TSSは、

TSS =\sum_{i=1}^nw_i(y_i -\frac{\sum_{i=1}^n w_i y_i} {\sum_{i=1}^n w_i})^2 (補正) (26)
TSS=\sum_{i=1}^n w_iy_i^2  (未補正) (27)

F値で、フィットモデルがモデル「y=一定」と、有意に異なるかどうかを検定します。

また、p値、または、有意水準は、F検定と一緒に出力されます。p値が、フィットモデルがモデル「y=一定」と有意に異なっていることを意味する\alpha\,\!よりも小さい場合、帰無仮説を棄却できます。

ある値に切片を固定している場合、F検定のp値には意味が無く、切片一定としない線形多重回帰とは異なります。

適合度検定表

不適合度を実行するには、連結フィットモードが選択されている場合に、少なくともX値がデータセット内や複数データセット内で反復できるように、反復観測、つまり、「複製データ」が必要になります。

複製データでフィットに使われている表記:

y_{ij}は、データセット中のi番目のx値における、j番目の観測値です。
\bar{y}_{i}は、i番目のx値における全てのy値の平均です。
\hat{y}_{ij}は、i番目のx値における、j番目の観測値の予測反応です。

残差平方和は、次の通りです。

RSS=\sum_{i}\sum_{j}(y_{ij}-\hat{y}_{ij})^2
LFSS=\sum_{i}\sum_{j}(\bar{y}_{i}-\hat{y}_{ij})^2
PESS=\sum_{i}\sum_{j}(y_{ij}-\bar{y}_{i})^2

非線形フィッティングの適合度検定表:

DF 平方和 平均平方 F値 Prob > F
不適合度 c-k-1 LFSS MSLF = LFSS / (c - k - 1) MSLF / MSPE p-値
純誤差 n - c PESS MSPE = PESS / (n - c)
誤差 n*-k RSS
Note:

切片がモデルに含まれてる場合、 n*=n-1 です。それ以外は、 n*=n で平方和の合計は未補正となります。勾配が固定の場合、 df_{Model}= 0です。

cは、明確なx値の数を示します。切片が固定である場合、適合度検定のDFは、c-kになります。

共分散行列と相関行列

多重線形回帰の共分散行列は以下によって計算されます。

Cov(\beta _i,\beta _j)=\sigma ^2(X^{\prime }X)^{-1}

(28)

特に、単一線形回帰の共分散行列は以下によって計算されます。


\begin{pmatrix}
Cov(\beta _0,\beta _0) & Cov(\beta _0,\beta _1)\\
Cov(\beta _1,\beta _0) & Cov(\beta _1,\beta _1)
\end{pmatrix}=\sigma ^2\frac 1{SXX}\begin{pmatrix} \sum \frac{x_i^2}n & -\bar x \\-\bar x & 1 \end{pmatrix}

(29)

2つのパラメータ間の相関は、

\rho (\beta _i,\beta _j)=\frac{Cov(\beta _i,\beta _j)}{\sqrt{Cov(\beta _i,\beta _i)}\sqrt{Cov(\beta _j,\beta _j)}}

(30)

信頼帯と推定帯

信頼帯

フィット関数の信頼区間は、フィット関数の値の推定値が独立変数の特定の値でどれほど良いかを示します。フィット関数の正確な値が信頼区間に含まれるように、100\alpha%で指定することができます。ここで、\alpha は指定した信頼水準です。フィット関数に対して、この定義済みの信頼区間は、以下の式で計算できます。

\hat{Y}\pm t_{(\alpha/2,dof)}{\left [\chi^2\mathbf{x}\mathbf{C}\mathbf{x}' \right ]}^{\frac{1}{2}}

(31)

ここで

\mathbf{x}=\left [ \frac{\partial y}{\partial \beta_0},\frac{\partial y}{\partial \beta_1},\frac{\partial y}{\partial \beta_2},...,\frac{\partial y}{\partial \beta_k} \right ]

(32)

そして、 C は共分散行列です。

推定帯

指定した有意水準αに対する推定帯は、特定の独立変数の値で、一連の繰り返し実験によるすべての測定データの100\alpha%がその範囲に含まれるような区間です。フィット関数に対して、この定義済みの推定区間は、以下の式で計算できます。

\hat{Y}\pm t_{(\alpha/2,dof)}{\left [\chi^2(1+\mathbf{x}\mathbf{C}\mathbf{x}') \right ]}^{\frac{1}{2}}

(33)

X/YからY/Xを検索

残差プロット

残差タイプ

作図するには、標準正規化スチューデント化スチューデント化残差から1つの残差タイプを選択します。

残差と独立変数

残差散布図res vs.独立変数x_1,x_2,\dots,x_kでは、それぞれのプロットは別のグラフに配置されます。

残差vs.予測値

残差散布図 res vs. フィット結果\hat{y_i}

残差vs.データ順序

res_i vs. 順番i

残差のヒストグラム

残差のヒストグラムres_i

残差のラグプロット

残差res_i vs. ラグ残差res_{(i–1)}

正規残差確率プロット

残差の正規確率プロットは、分散が正規分布しているかどうかを調べるのに使用します。結果のプロットはおおよそ線形で、誤差範囲は正規分布していると仮定することができます。プロットはパーセンタイル対順序化された残差をベースにしており、パーセンタイルは次のように仮定されます。

\frac{(i-\frac{3}{8})}{(n+\frac{1}{4})}

ここで、n はデータセットの合計数で、ii 番目のデータです。なお、正規確率プロットとQ-Qプロットについてをご覧ください。