生物科学研究所のロゴマーク
生物科学研究所 井口研究室
Laboratory of Biology, Okaya, Nagano, Japan

回帰と相関,知っているようで知らない,その本質

Regression and correlation with Excel

Excel の回帰分析を例として

井口豊 (生物科学研究所


更新:

相関係数 0 (無相関)の検定と,回帰直線の傾き(回帰係数)0 の検定が,同等のものであることは,このページ末の注 3 に示した。

相関係数の強弱の基準には,ギルフォードの基準(Guilford's Rule of Thumbが広く用いられる。それに関しては,私の Yahoo! 知恵ノート参照。

日本語の統計学の教科書やウェブサイトでは,相関係数の強弱の区分が示されることは多いが,Guilford の名前や文献が挙がることは稀である。相関係数の強弱区分は慣用的である,ということで済ませている場合も多い。


1. はじめに

このページでは,回帰直線とはどのようなものなのか理論的に考え,最小二乗法による直線回帰の「誤解」について取り上げたい。特に, MicrosoftExcel (エクセル)を用いて回帰分析の具体例を示し,その理論的背景を考えてみたい。

これは,読んで字のごとく,データ点からの二乗和が最小になる近似式(適合式),を求めるものである。ところが,この「データ点からの距離」の取り方が,「くせもの」なのである。

最小二乗法は頻繁に行われる直線近似法である。 Excel などの表計算ソフトでも容易に出来る。しかし,それがゆえに,どんな計算法か知らずに適用されてきているのも事実である。

ここでは,最も単純な2変数x, yのデータを使った最小二乗法による直線(いわゆる回帰直線, regression line)への適合と検定について考える。ここでの計算は,一般的に最小二乗法(Least squares methodと呼ばれている,通常最小二乗法Ordinary Least Square: OLS)である。わざわざ「通常」と呼ばれるからには,通常でない最小二乗法もあることを念頭において欲しい。以下,この通常最小二乗法を単に,最小二乗法と呼ぶ。

一方で,二乗でなく絶対値の和が最小となる回帰式も当然ある。それを最小絶対値法(Least absolute value method)による回帰と呼ぶ。これが意外と知られておらず,回帰や近似と言えば,最小二乗が当然だと思い込んでいる人も多い。それに関しては,私の Yahoo! 知恵ノート参照。

2. 回帰直線は,どのようにして決まるか?

次のような8個の対から成るデータ, x, y を考えよう。

    -------------
       x     y
    -------------
       0      45
       5      25
       6      75
      11      55
      12     105
      17      85
      18     135
      23     115
    --------------

Excel のグラフ作成機能で作った散布図( Scatter Diagram)が図1。 この散布図の場合は,本当は2つの右上がりの直線に適合させるべきであるが,回帰分析の学習だと思って,単一直線への適合を考えてみよう。

左は,例えば小中学生などがやる,最も適合しそうだと直感で(勘で),記入される直線。大人でも,こんな直線を引くのではないかと思うし,私も理解できる。

一方,右は,グラフのデータ点を右クリックして,

グラフの近似曲線の追加 → オプション → 線形近似→ オプション→ 「グラフに数式を表示する」と「R2乗値を表示する」

を選んで示された直線と線形(つまり,一次関数の)回帰式である。これは,最小二乗法で作成された直線でもある。

目測による直線適合(左)と最小二乗法による直線適合(右)

図1. 同一データについて,目測による直線適合(左)と最小二乗法による直線適合(右).

回帰式を求めるだけなら,もとの表で

挿入 → 関数

とたどって, SLOPE 関数と INTERCEPT 関数からも求められる。

左右のグラフの直線の傾きに注目して欲しい。右の最小二乗法では, x が大きくなるほど,直線が下側へズレている。

なぜ,このようになるのだろうか?

これは図2に示すように,最小二乗法が縦方向,つまり 軸方向の距離のみを最小とするように計算されているからである。だから,単純に点と点の中間を通るような直線にはならないのである。

最小二乗法によるデータ点から直線までの距離の計算.y 軸方向のみ最小化される

図 2. 最小二乗法によるデータ点から直線までの距離の計算.y 軸方向のみ最小化される.

一般的には,変数 x は独立変数または説明変数と呼ばれ,変数 y は従属変数または目的変数と呼ばれる。図 2 から分かるように,変数 x は誤差を含まないと仮定されるため,確定変数または固定変数(fixed variableと呼ばれることもある。

3. 最小二乗法の理論的計算

では,最小二乗法とは,どのような計算で,回帰直線は,どのように決定されるのだろうか?

まず,回帰式として1次関数 Y = ax + ba, b はパラメータと呼ばれる定数)を想定する。想定する,とは,逆に言えば,2次関数とか,その他の非線形関数とかを想定する場合もあるのだ。

ある点 xi における観測値を yi として,その時の計算値 Yi との差 E を考える。

E = yi - Yi

これを残差(residualと呼ぶ。最小二乗法とは,全ての点について,この残差の2乗を求め,その総和

Q = Σ E2
 = Σ (yi - Yi)2
 = Σ {yi - (axi + b)}2    (1)

これが最小になるように,パラメータ a, b を決定する方法である。
なお,a は直線の傾きを表すが,回帰では回帰係数と呼ぶ。

この式からも分かるように,y 軸方向の距離だけを最小にする,これが最小二乗法の本質だと言える。

ここで,「本質」と述べたのは,ここで取り上げた1次関数に限らず,2次関数だろうが,その他の非線形関数だろうが,同様の方法でパラメータを決定するからである。

(1) の xi, yi には,測定値が入るので,結局,この式は 変数a, b の式と見ることができる。 したがって,Qの最小値を求めるためには,a, b で偏微分して,その連立方程式を解けば良い。

式を見やすくするために,あらためて, x = xi, y = yi と置いて,以下の計算をする。

具体的には,Q を a, b で偏微分して

∂Q/∂a = -2 Σ x {y - (ax + b)}    (1a)
∂Q/∂b = -2 Σ {y - (ax + b)}    (1b)

それぞれ = 0 と置いて

Σ xy - a Σ x2 - b Σ x = 0
Σ y - a Σ x - bn = 0    (2)

これを,まずa について解くと,

回帰係数の計算

(3)


これを整理すると

回帰直線の傾き

(4)


ここで,mx, my は,それぞれ x, y の平均を表す。(3), (4)のどちらを使っても良いが,(4) の方が,傾き a の性質を直感的に理解しやすい。

(1) より,

定数項の計算

(5)


このy 切片 b の式に注目してほしい。
元の式 y = ax + b に,(5)を代入すると

y = ax + b
 = ax + my - a.mx
 = a ( x - mx ) + my  (6)

すなわち,回帰直線は,x, y の平均の点 (mx, my) を必ず通過する。この点をデータの重心と呼ぶ。

回帰直線がデータの重心を通るという性質は,この後,データの変動や統計量の検定を考えるとき重要となる要素である。

4. 回帰係数と相関係数の違い

では,回帰係数(= 傾き)と相関係数は,どのように違うのだろうか?

(3)と(4)の傾き a の分母の計算に注意して欲しい。x のみで決まり,y は全く関係しないのである。これは x の平均からの変化に伴って,y がその平均からどう変化するかを表した式といっても良い。

「そんなこと y = ax + b の式なんだから,言われなくても分かる」と言うかもしれないが,これが意外と理解されてないのである。

例えば, x, y を入れ替えたデータから出来る回帰直線は,改めて y 軸方向の残差が最小になるようにして得られる直線であり,元の直線の x, y を入れ替えて出来る直線と必ずしも一致しないのである。

あとで具体例を示すが, x, y を入れ替えたデータから出来る回帰直線は,必ずしも元の直線の逆関数とはならないのである。

一方,通常,相関係数と呼ばれる,積率相関(Pearson 相関)係数 r は次のように計算される。

Pearson積率相関係数

(7)


これは, x, y に関して対称式であり, x, y を入れ替えても全く同じ式となる。これが, ExcelCORREL 関数の計算式である。

つまり,回帰係数は, y に対する x の関係,という一方通行的な係数なのに対し,相関係数は,xy 相互間の関係,という双方向な係数なのである。統計データの分析に携わる人の中にも,この点を十分理解せずに,両示数を利用している人がいるので注意してほしい。

ここで,相関係数の分母(偏差平方和)と分子(積和)の大小関係を見てみよう。

n 個の x, y データについて,その各偏差をベクトル xy の成分として考える。

 x = (x1-mx, x2-mx, ... xn-mx)
 y = (y1-my, y2-my, ... yn-my)

すると,積和 Σ(x-mx)(y-my) は,ベクトル x, y の内積
[x, y] = (x1-mx)(y2-my) + (x2-mx)(y2-my) + ……..+ (xn-mx)(yn-my)
と考えられる。

また,同様に, x の偏差平方和 Σ(x-mx)2y の偏差平方和 Σ(y-my)2,は,それぞれベクトル xy の大きさ(絶対値)の2乗と考えられる。

|x|2 = (x1-mx)(x2-mx) + (x2-mx)(x2-mx) + ……..+ (xn-mx)(xn-mx)
|y|2 = (y1-my)(y2-my) + (y2-my)(y2-my) + ……..+ (yn-my)(yn-my)

ここで,コーシー・シュワルツ(Cauchy–Schwarzの不等式より,

| (x, y) |≦|x|・|y|

よって,積和と x, y の偏差平方和の関係は

積和と偏差平方和

となり,相関係数の絶対値が 1 以下であることが証明される。

コーシー・シュワルツの不等式を利用すると,相関係数の分子の絶対値の大きさが,分母の大きさによって抑えられていることが良く分かる。

式 (4) で,回帰の傾き a を,もう一度見て欲しい。分母の式から,y に対する x の関係,という一方通行的な示数であることは既に述べた。

しかし,この分母を

回帰係数と相関係数

と分解し,ひとつのルート内の x を y に入れ替えたのが,相関係数 r なのである。

つまり,相関係数の分母では,x の偏差平方和と y の偏差平方和の相乗平均を取り,x の変化も y の変化も考慮しているのである。

ちなみに,x と y に対し,その積の平方根

相乗平均

を相乗平均と言い,相加平均との関係で高校時代に学んだ人も多いだろう。

ところで,x, y 相互間の結びつきの強さを示すには,相関係数r より,それを2乗したr2 を用いるほうが良い。

決定係数

これを決定係数(coefficient of determinationと呼ぶ。

ただし,注意してほしいのは,これは複数ある決定係数の定義の一つに過ぎないということである。

EXCEL回帰分析のバグ?(注2)(注3)を参考にしてほしい。また,決定係数の定義が複数あることから生じる具体的問題に関しては,決定係数R2の違い: Excel,OpenOffice,および統計解析ソフトRを用いてを参照のこと。

この決定係数は,回帰係数と結びつきが強く,その意味でも統計学上の重要な指標である。

例えば,決定係数を以下のように,変形してみる。

決定係数の変形

(8)


ここで,右辺の左の分数は回帰による y 方向の変動を表すことは以下のように証明できる。

まず,第 3 章冒頭に述べたように,測定値 x を回帰直線に代入して得られる値を Y とする。

すなわち,
Y = a.x + b    (9)

y の計算値(推定値) Yと y の測定値の平均 my との差の平方和 Sr は,
Sr = Σ(Y - my)2    (10)

これが,回帰による(回帰直線上の)y 方向の変動量である。

ここで,回帰直線 y = a.x + b が重心 (mx, my) を通ることに注目すると,
my = a. mx + b    (11)

(9) (11) を (10) に代入し整理すると
Sr = Σ{a.x + b - (a. mx + b)}2
 = Σ{a.(x - mx )}2
 = a2.Σ(x - mx )2
   (12)

ここで,式 (4) の傾き a を再度取り上げる。

回帰直線の傾き

(4)


この a を (12) に代入すると以下のようになる。

平方和の変形

さらにこの Sr を決定係数の式 (8) に代入すると以下のようになる。

決定係数の変形

ここで,最後の式の分母は, データの y 方向の全変動を表すので,これを
Sy = Σ(y - my)2    (13)
と置くと

決定係数は全変動量に占める回帰変動量の割合

(14)


すなわち,決定係数 r2 は,y 軸方向の全変動量 Sy に占める回帰変動量 Sr の割合を表す指標とも言えるのである。

相関係数が,x, y 相互間の双方向な関係を表す示数である,と既に述べた。一方,決定係数は,相関係数を2乗しただけであるが,随分異なる性質を持つことが理解できるだろう。

ところで,相関係数が 0,つまり無相関なら,回帰係数(= 傾き)はどうなるだろうか?傾き a の式(4)と相関係数の式 (7) を再び見よう。

回帰直線の傾き

(4)


Pearson積率相関係数

(7)


相関係数が 0 の時を考えると, r において分母が0にならず,分子が 0 となる時である。その時,回帰係数 a の分母も 0 にならず,分子は 0 となる。すなわち,相関係数が 0 ならば,回帰係数(= 傾き)も 0 なのである。このことを理解していない人がかなり多い。

さらに厳密に,相関係数 0 (無相関)の検定と,回帰直線の傾き(回帰係数)0 の検定が,同等のものであることは,このページ末の注 3 に示した。

なお,相関係数の強弱の基準には,ギルフォードの基準(Guilford's Rule of Thumbが広く用いられる。それに関しては,私の Yahoo! 知恵ノート参照。

日本語の統計学の教科書やウェブサイトでは,相関係数の強弱の区分が示されることは多いが,Guilford の名前や文献が挙がることは稀である。相関係数の強弱区分は慣用的である,ということで済ませている場合も多い。

5. 最小二乗法のデータへの適用

ここでは,第2章冒頭に示した8対の x, y のデータに最小二乗法を適用し,回帰直線を求めてみよう。

まず,x と y の平均,mx, my を求めて
mx = 11.5,my = 80

これら平均と,各 x と y の差(偏差)を求める。 Excel では,表 1 のような計算を,まず行う。

表1. 最小二乗法で回帰直線を求めるときの基本計算

回帰直線を求めるときの基本計算

ここから, x と y の偏差積和 Σ(x-mx)(y-my),および,x の偏差平方和 Σ(x-mx)2 を計算する。

Σ(x-mx)(y-my) = 1600
Σ(x-mx)2 = 410

したがって,式(4)より,直線の傾き a は,

a = 1600/410
 = 3.9024

また,式 (6) より,回帰直線は ,重心を通るので,
my = a.mx + b

よって,
b = my - a.mx
 = 80 - 3.9024*11.5
 = 35.122

となり,結局,回帰直線は,
y = 3.9024x + 35.122    (15)

これが図1右に示したEXCELグラフの回帰直線の式である。

ここで,最小二乗法の定義式として,式 (1) をもう一度取り上げよう。
Q = Σ {yi - (axi + b)}2    (1)

少し面倒だが,実際に,この式に測定値(表1の x, y)を入れて計算してみよう。

Q = {45 - (a*0 + b)}2 + {25 - (a*5 + b)}2 + {75 - (a*6 + b)}2 + {55 - (a*11 + b)}2
  + {105 - (a*12 + b)}2 + {85 - (a*17 + b)}2+ {135 - (a*18 + b)}2 + {115 - (a*23 + b)}2

展開,整理すると
Q = 1468a2-17920a+8b2-1280b+184ab+61000

したがって,Q は,a, b の2次曲面となる。

gnuplotを使って,この関数の 3 次元プロットを見よう(図3)。「すりばち」を横から押しつぶしたような形状であり,最小値が (a, b) = (4, 40) 付近だと視覚的に確認できる。

最小二乗曲面 Q

図3. 最小二乗曲面 Q

図4のように,Q-a 面の方向,つまり, a 軸を真横から見ると,a の最小値が 4 付近にあることが理解できる。

最小二乗曲面 Q の Q-a 面

図4. 最小二乗曲面 Q の Q-a 面

さらに,図5 のように,Q-b 面の方向,つまり,b 軸を真横から見ると,b の最小値が 30 から 40 付近にあると良く分かる。

最小二乗曲面 Q の Q-b 面

図5. 最小二乗曲面 Q の Q-b 面

これらが,a, b の最小値,3.9024 と 35.122 を幾何学的に表現している。

ここで,先ほど述べた,x, yを入れ替えたデータから出来る回帰直線は,改めてy軸方向の残差が最小になるようにして得られる直線であり,元の直線のx, yを入れ替えて出来る直線と必ずしも一致しない,ということを思い出して欲しい。

例えば,今まで見てきた回帰直線
y = 3.9024x + 35.122    (15)
の x と y を入れ替えて,逆関数を考えよう。

x = 3.9024y + 35.122
として,y について整理すると,
y = 0.25625x - 9.000
となる。

一方,データの x と y を入れ替えて,回帰直線を求めると
y = 0.16327x - 1.5612

これは逆関数とは全く異なる式となる。統計を頻繁に扱う人でも,この点を認識不足の人がいる。十分注意してもらいたい点である。

ところで,前章の最後に触れた,相関係数が 0 の時の具体例を見るため,次のような(x, y)のデータを考えよう。

x  1  1/√2   0  -1/√2  -1  -1/√2   0   1/√2
y  0  1/√2   1   1/√2   0   -1/√2  -1  -1/√2

これは,単位円 x2 + y2 = 1 上に,x 軸上から45゜ごとに取った点である。その散布図と回帰直線を示す(図6)。

最小二乗曲面 Q の Q-a 面

図6. 単位円 x2 + y2 = 1 上に 45゜ごとに取った点の分布と回帰直線(赤).

相関係数が 0,すなわち無相関であることは,すぐ分かるだろう。しかし,回帰直線が y = 0 であることが,すぐ分かるだろうか?

直感的には,このような分布だと,回帰直線は,どのようにでも引けるとも言えるし,どのようにも引けないとも言えそうだ。

しかし,回帰直線を求める上で重要な2点を,もう一度考えて欲しい。

  1. データ分布の重心を通り,この場合, (0, 0)。
  2. y 方向の最小二乗のみを考える。

そうすると何となく,y = 0 の回帰直線が見えてくるであろう。この場合の最小二乗法定義式の3 次元プロットも下に示す。

単位円の最小二乗法定義式の3 次元プロット

図7. 単位円 x2 + y2 = 1 上に 45゜ ごとに取ったデータ(図6)の最小二乗法定義式の3 次元プロット.

(a, b) = (0, 0) で最小値を取る2次曲面となる。

6. 回帰の検定

6.1. 係数と定数項がゼロであるかどうかの検定

回帰直線 y = ax + b のパラメータの検定について考える。厳密な証明は専門書に譲り,ここでは,回帰直線の性質をもとに,標本の平均値についての t 検定からの類推で,パラメータの検定に至る。

なお,回帰式における係数や定数項の検定と言えば,通常は,それらがゼロであるという帰無仮説(null hypothesis)に対する検定である。それらがゼロ以外の場合に対する検定は,6.2. 係数と定数項がゼロ以外の場合に対する検定を参照してほしい。

まず,傾き a の検定を考える。回帰直線は y 軸方向の残差のみが問題となることに再び注意しよう。

すると,測定値 y と計算値 Y の差の平方和(残差平方和)
Se =Σ (y-Y)2    (16)
を自由度 n-2 で割った残差分散(residual variance) Ve を考え,それを回帰直線の誤差分散の不偏推定量(不偏誤差分散)と考えるのは不自然ではない。

この残差分散(不偏誤差分散) Ve を,よく知られた標本不偏分散 s2 と対応させると分かりやすいだろう。

残差分散と標本不偏分散

(17)


Veで,自由度 n-2 としたのは,標本サイズ n からパラメータ a, b の2個分引いたものである。

右式では,平均の周りに正規分布に従って観測値が散らばることを利用し,左式では,回帰直線の周りに正規分布に従って観測値が散らばることを利用している。当然ながら,両者ともデータのバラツキの正規性を前提としているのである。

したがって,この残差分散は回帰直線からのデータのバラツキを表す統計量と言える。誤解されがちだが,決定係数 R2 は,そういう意味の指標ではない。それに関しては,私の Yahoo! 知恵ノート参照。

傾き a の式 (4) をもう一度見よう。

回帰係数の計算

(4)


回帰直線は y 軸方向の残差のみが問題となる,に注目すると, a の標準誤差を考える上で,問題となるのは分子だけである。それゆえ,この式の分子を残差分散 Veで置き換え,その平方根を取ったものを,傾きの標準誤差 Sa と見なすのも不自然ではない。

この傾きの標準誤差 Sa を,よく知られた 平均値の標準誤差 SE と対応させると分かりやすいだろう。

傾きと平均の標準誤差

(18)


平均値では,標本サイズ n を大きくするほど SE が小さくなり,平均値がより正確になる,と理解できる。

一方,傾きでは,x の平均 mx から離れたデータが増えるほど,Sa は小さくなる。これは,式 (6) に示したように,回帰直線が重心(x, y の平均)を通ることから,単に全体の標本サイズを大きくすることより,重心から離れたデータ量を増やした方が,より正確な傾きが得られることは明らかであろう。

この Sa を用いて,回帰直線(15)
y = 3.9024x + 35.122
の傾き 3.9024 の検定を進める。

まず,元のデータの x 値を,この回帰直線に代入し,元のデータの y 値と計算値 Y の差(残差)を計算する。

表2. 回帰直線からの残差の計算
Y = 3.9024x + 35.122 として計算

回帰直線からの残差

残差平方和 Σ (y-Y)2 = 3556.098 を,自由度 n-2 = 6 で割った残差分散 Veを求める。

Ve = Σ (y-Y)2 / (n-2) = 592.6829    (19)

さらに, x 偏差平方和を計算。

Σ (x-mx)2 = 410    (20)

(19)と(20)を(18)の Sa に代入して

傾きの標準誤差

結局,
Sa = 1.2023

求めた回帰直線(15)の傾き 3.9024 が,ある傾き b と異なるかどうかは,t 検定で調べられる。

傾きのt値の計算式

(21)


自由度は n-2 = 6

この検定は,標本平均が,ある値と異なるかどうか調べる 1 標本 t 検定と全く同じ形式である。

標本平均が 0 と異なるかどうか検定できるように,傾きが 0 と異なるかどうかも,当然,検定できる。

式 (21) を用いて,傾き 0 に対する検定を行うと,

傾きのt値

結局,
t = 3.246

自由度 6 に対して,
p = 0.018
という結果が出る。

次に,回帰直線 y = ax + b の y 切片 b の検定を考えよう。ここ,式 (6) に示したように,回帰直線が重心(平均),(mx, my) を通ることから,

my = a.mx + b
すなわち
b = my - a.mx と置ける。

したがって,b の分散 V[b] は,
V[b] = V[my - a.mx]
   = V[my] + V[a.mx] + Cov[my, b]

ここで,回帰直線が傾きいかんに関わらず重心を通る,ということは,my と b が独立であることを示し,
Cov[my, b] = 0

また,分散の基本的計算から
V[a.mx] = mx2V[a]

したがって,
V[b] = V[my] + mx2V[a]

ここで,V[my] は, y の平均値の分散である。
式 (17) で定義された残差分散(不偏誤差分散) Ve を用い,式 (17) の標本平均の標準誤差 SE からの類推で,
V[my] = Ve / n
と置ける。

また,V[a] は傾きの分散だから,式 (18) の傾きの標準誤差 Sa の2乗である。

結局,

定数項の分散

したがって, y 切片 b の標準誤差 Sb は,

定数項の標準誤差

ここで,n = 8
x の平均mx = 11.5

式(18)より,残差分散Ve = 592.683
式(19)より,x 偏差平方和 Σ (x-mx)2 = 410

これらを代入して
Sb =16.2869

求めた回帰直線(式7)のy 切片35.122が,0と異なるかどうか,すなわち,回帰直線が原点を通るかどうか,やはり t 検定で調べられる。

定数項のt値

結局,
t = 2.157 自由度 6 に対して,
p = 0.074 という結果が出る。

以上の結果を整理すると,
y 切片 35.122,  Sb =16.2869,  t = 2.157, p = 0.074
傾き 3.9024,  Sa = 1.2023,  t = 3.246,   p = 0.018
となる。

一方,Excelを用いて,

ツール → 分析ツール → 回帰分析

とたどって,与えられたデータ (x, y) を分析した結果は以下の通りである。

表3. EXCEL分析ツールの回帰分析による傾きとy切片の検定

Excel分析ツールによる回帰分析

これは上述の分析結果と同じである。

6.2. 係数と定数項がゼロ以外の場合に対する検定

これまで述べたように,回帰係数や定数項を検定する場合,教科書やウェブページでは,それらが 0 である場合の検定が説明されることが多い。しかしながら実際には,0 以外の任意の係数や定数項に対して検定可能である。係数や定数項 0 に対して検定可能ならば,一般性を失うことなく,任意の係数や定数項に対して検定可能であることは,意外と理解されていないらしい。この問題に関しては,次の ブログ参照。

統計解析ソフト R を用いて,任意の係数や定数項の検定を行なう場合,線形回帰であろうが非線形回帰であろうが, nls 関数を使えば良い。

7.回帰分析とは

前章で述べたことは回帰分析の一部とも言える。しかしながら,通常,回帰分析と言うと,回帰を利用した分散分析のことである。ここでは,その内容を見よう。

既に述べてきたことから察しがつくかもしれないが,データの y 方向への全変動 Sy は,回帰による y 方向への変動 Sr と誤差変動(残差平方和)Se の和であると考えることができる。

すなわち
Sy = Sr + Se    (22)

ただし,ここで言う誤差変動とは,単に測定誤差の意味だけでなく,一次回帰式のx以外の要因も含んでいる。

これらは,それぞれ式 (13), (10), (16)より,

Sy = Σ(y - my)2
Sr = Σ(Y - my)2
Se =Σ (y-Y)2

ここで改めて,y は測定値,Y は回帰式による計算値 ax + b,my は測定値 y の平均であることに注意しよう。いずれも y 方向の値のみを問題としていて,最小二乗回帰の特徴が,ここにも現れている。

では具体的に,第5章で求めた式 (15) の回帰直線の回帰分析を行ってみよう。
y = 3.9024x + 35.122    (15)

これまでと重複した計算があるが,分かりやすくするため,それも含めて以下に Excel の計算方法を記述する。

データ y の平均  my = 80
傾きパラメータ   a = 3.9024
y 切片パラメータ  b = 35.122

表4. データの集計方法

回帰分析データの集計法

式 (22) で述べたとおり,
Sy = Sr + Se    (22)
となっている。

回帰分析の検定を行う前に,まず自由度を考えると次のようになる。

回帰の自由度(DFr) = 説明変数x の数 = 1
残差の自由度(DFe) = 標本サイズ n - パラメータ数(今回は,a, b の 2)= 8 - 2 = 6
全体の自由度(DFy) = 標本サイズ n - 平均を求めたことによる1 = 8 - 1 = 7

当然だが,
DFy = DFr +DFe

次に,全分散 Vy, 回帰分散 Vr, 残差分散 Ve を以下のように考える。

Vy = Sy / DFy
 = 9800 / 7
 = 1400
Vr = Sr / DFr
 = 6244 / 1
 = 6244
Ve = Se / DFe
 = 3556 / 6
 = 592.67

このなかで,回帰分散 Vr と残差分散 Ve の比に対して F 検定を行うと回帰分析である。

F = Vr / Ve
 = 6244 / 592.67
 = 10.5

自由度(1, 6)に対して,p = 0.018

Excelを用いて,

ツール → 分析ツール → 回帰分析

とたどって,与えられたデータ (x, y) を分析した結果は以下の通りである。

表5. EXCEL分析ツールの回帰分析の分散分析

Excel分析ツールによる回帰分析

これは上述の分析結果と同じである。

以下の表6に,Excel分析ツールによる回帰分析の全体像を示した。なお,傾きや切片の上下限95%の部分は省略した。

表6. EXCEL分析ツールの回帰分析表のまとめ。傾きや切片の上下限95%の部分は省略

Excel分析ツールによる回帰分析

ここで,いくつか注意点を述べたい。

まず,決定係数 R2 (重決定 R2)の定義は,回帰変動と合計変動の比である。

重決定 R2

ただし繰り返し注意するが,既に第4章で述べたように,決定係数 R2 には複数の定義が存在する。ここでは,Excel分析ツールによる回帰分析の定義,すなわち,相関係数の2乗となる決定係数について解説していることを念頭に置いてほしい。

また,補正 R2 は自由度調整済み決定係数 R*2 とも呼ばれ,その定義は,残差変動と合計変動を,それぞれ自由度で割ったものの比を 1 から引いたもの,すなわち,残差分散と合計分散の比を1から引いたものである。

自由度調整済み決定係数

ここで注意したいのは,Excel の結果が示す表には,不親切にも,

合計分散 = 合計変動 / 合計自由度
すなわち
Vy = Sy / DFy
が示されていないことである。

また,この自由度調整済み決定係数 R*2 は,値によっては,Excelでバグが出る場合があるので注意が必要である。この点に関しては,私の以下のウェブページを参照。

最後に,観測された分散比は,X 値 1 (傾き)の t 値の 2 乗であることに注意しよう。すなわち,
10.53497942 = 3.24576332

これは自由度 (1, n) の F 検定は,自由度 n の t 検定に等しいからである。当然,確率 p も等しく, 0.017559 となっている。

8. おわりに

粕谷 (1998) は,生物統計の入門書の中で, 最小二乗法がy軸方向の距離のみを最小とするように計算されていることを知らないのは,「恥だ」と述べている。しかし実際には,そのような人も多いのではないかと思う。

また,2 変数のグラフを書くときに,それをy 軸方向に拡大させて,ことさら変化を強調する方法を,Huff (1954) は「針小棒大法」と皮肉ったが,私も同意見である。

パソコンや,そのソフトが普及し,誰でもデータ分析できるようになってきた。しかし,少なくとも,自分が実行しようとしているデータ分析法や検定法の数学的背景も学習してほしいと思っている。そうしないと,分析結果の解釈やグラフ化の過程で,誤った適用や解釈が出てきても,本人は全く気づかないでいる,という困った例も増えてくるのではないかと危惧している。

参考文献

Huff D. (1954) How to lie with statistics.(統計でウソをつく方法,高木秀玄・訳, 1968)ブルーバックス,講談社.

粕谷英一 (1998) 生物学を学ぶ人のための統計のはなし. 文一総合出版.

(注 1)

決定係数 R2 は,パラメータに関して線形結合である場合,つまり線形モデルであることを仮定しており,非線形モデルの適合度は正しく評価できない。しかしながら,それでもなお誤解(誤用?)して使う人が多いようである。それに関しては,私の Yahoo! 知恵ノートも参照。

このことは,例えば,群馬大・青木氏なども繰り返し指摘している。

さらに,学術論文でそれを指摘したものも参照して欲しい。
Spiess and Neumeyer (2010)
An evaluation of R2 as an inadequate measure for nonlinear models in pharmacological and biochemical research: a Monte Carlo approach
BMC Pharmacology 10:6

オープンアクセスで, 全文 見られる。

この論文の第1ページで,非線形回帰では,全変動が回帰変動と残差変動の和にならない,ことを指摘している。

(注 2)

UCLAThe Institute for Digital Research and Education (IDRE)
What are pseudo R-squareds?
で,線形回帰における決定係数 R2 に類似したモデル適合度指標(goodness-of-fit measure)として,非常に多くのR2 (Pseudo-R2) と呼ばれる指標を紹介している。

IDRE では, Nonlinear Regression に関して, SASSPSS における定義や利用上の注意事項も解説しており,非常に参考になる。

特に, SPSS の解説末尾に書かれた

Pseudo-R-squared: Many different measures of pseudo-R-squared exist.・・・however, none of them can be interpreted exactly as R-squared in OLS regression is interpreted.

という指摘に注意すべきである。

なお,関連問題として,
Excel回帰分析のバグ?の(注)

および,内容が重複するが,Yahoo!知恵ノート
決定係数R2の誤解: 必ずしも相関の2乗という意味でなく,負にもなるし,非線形回帰には使えない も参照。

(注 3)

回帰直線の傾き 0 の検定と,相関係数 0 (無相関)の検定が同等であることは,以下のようにして示すことができる。

まず,回帰直線の傾き a を表す式 (4) を取り上げる。
回帰直線の傾きを示す式。

次に,その傾き a を持つ回帰直線が平均値(重心)を通ることを示した式 (6) を考える。
回帰直線が平均を通る式。

それが,式 (17) 残差分散 Ve の中の計算値(予測値) Y となるので,結局,Ve は以下のようになる。
回帰式の残差分散の式。

この分子を展開すると
回帰式の残差分散の展開式。

さらに,式 (18) より,傾きの標準誤差 Sa は以下のようになる。
傾きの標準誤差の式。

その上で,
帰無仮説:傾き a = 0
の検定統計量 t は以下のようになる。
傾きの検定統計量 t の式。

さらに説明を続けるが,以上の式は煩雑なので,x,y の偏差を,それぞれ X,Y と置き換える。すると,以下のようになる。

回帰直線の傾きを示す式を簡略化。

回帰式の残差分散の簡略化。

傾きの標準誤差を示す式を簡略化。

ここでまず。a を Ve に代入して整理すると
回帰式の残差分散を x と y で表示する。

この残差分散 Ve を標準誤差 Sa に代入して整理すると
傾きの標準誤差を x と y で表示する。

この標準誤差 Sa と傾き a を,検定統計量 t の式に代入すると
傾きの検定統計量 t を x と y で表示する。

少し整理すると
傾きの検定統計量 t を x と y で表示し整理する。

さらに変形すると
傾きの検定統計量 t を x と y で表示し変形する。

ここで式 (7) より,相関係数 r は以下のように定義される。
相関係数 r を x と y で表示する

この相関係数 r を,検定統計量 t に代入すると
傾きの検定統計量 t を相関係数 r で表示する。

すなわち,相関係数 r を用いて,検定統計量 t を表現すると
傾きの検定統計量 t を相関係数 r で表示する。

これが,相関係数 0 (無相関)の検定統計量となる t 値の計算式である。つまり,回帰直線の傾き 0 に対する検定と,相関係数 0 に対する検定は,全く同じものなのである。