重回帰分析

統計解析 01 クラス : 第10回(06/19/03)

前回は、説明変量が一つである単回帰分析を紹介した。単回帰分析における説明変量が複数になった手法が重回帰分析であるのだが、残差(予測誤差)の二乗和を最小にするという考え方は同じなので、その原理は容易に理解できると期待している。

重回帰分析 : 2変量以上の説明する変量(説明変量)で 1変量を説明(目的変量)

説明変量が複数になる : 単 ===> 重
体重を身長と胸囲で説明したい。予測したい。
[体重]=a[身長]+b[胸囲]+c : 回帰係数を求めたい。
単回帰とアイディアは同じ : 残差(予測誤差)の二乗和を最小にする(最小二乗法)
説明される変量(目的変量)と平行に残差を取る。

プログラム : les1001.sas

/* Lesson 10-1 */ /* File Name = les1001.sas 06/19/03 */ data gakusei; infile 'all03a.prn' firstobs=2; input sex $ height weight chest jitaku $ kodukai carrier $ tsuuwa; proc print data=gakusei(obs=10); run; proc reg data=gakusei; : 回帰分析 model weight=height chest; : 複数変量を指定 output out=outreg1 predicted=pred1 residual=resid1; : 結果項目の保存 run; : proc print data=outreg1(obs=15); run; : proc plot data=outreg1; : 散布図を描く where weight^=. and height^=. and chest^=.; : 解析に使ったデータのみ plot weight*height; : plot weight*chest; : plot weight*pred1; : 観測値と予測値 plot resid1*pred1 /vref=0; : 残差と予測値(残差解析) plot resid1*height/vref=0; : 残差と説明変量(残差解析) plot resid1*chest /vref=0; : 残差と説明変量(残差解析) plot resid1*weight/vref=0; : 残差と目的変量(残差解析) run; : : proc univariate data=outreg1 plot normal; : 残差を正規プロットして確かめる var resid1; : run; :

出力結果 : les1001.lst
SAS システム 2 19:52 Thursday, June 12, 2003 Model: MODEL1 Dependent Variable: WEIGHT Analysis of Variance Sum of Mean Source DF Squares Square F Value Prob>F Model 2 6997.07013 3498.53506 86.274 0.0001 Error 83 3365.78475 40.55162 C Total 85 10362.85488 Root MSE 6.36802 R-square 0.6752 Dep Mean 59.54884 Adj R-sq 0.6674 C.V. 10.69377 SAS システム 3 19:52 Thursday, June 12, 2003 Parameter Estimates Parameter Standard T for H0: Variable DF Estimate Error Parameter=0 Prob > |T| INTERCEP 1 -110.841725 13.87776674 -7.987 0.0001 HEIGHT 1 0.678958 0.08661843 7.838 0.0001 CHEST 1 0.647517 0.09237349 7.010 0.0001 SAS システム 4 19:52 Thursday, June 12, 2003 K C H W J O A T R E E C I D R S P E I I H T U R U R S O S G G E A K I U E I B E H H S K A E W D D S X T T T U I R A 1 1 1 F 145.0 38.0 . J 10000 . . . 2 F 148.0 42.0 . J 50000 . . . 3 F 148.0 43.0 80 J 50000 DoCoMo 4000 41.4454 1.5546 4 F 148.9 . . J 60000 . . . 5 F 149.0 45.0 . G 60000 . . . 6 F 150.0 46.0 86 40000 . 46.6884 -0.6884 7 F 151.7 41.5 80 J 35000 . 43.9576 -2.4576 8 F 152.0 35.0 77 J 60000 DoCoMo 2000 42.2187 -7.2187 SAS システム 6 19:52 Thursday, June 12, 2003 プロット : WEIGHT*HEIGHT. 凡例: A = 1 OBS, B = 2 OBS, ... 100 + A | A A WEIGHT | A A A | B BABAB AAAA A A A AA | A A A A B BA BAFBC ABA ABBA 50 + A A ABA AAC C BB CB A | A B A | | | 0 + --+-----------+-----------+-----------+-----------+-----------+-- 140 150 160 170 180 190 HEIGHT SAS システム 7 19:52 Thursday, June 12, 2003 プロット : WEIGHT*CHEST. 凡例: A = 1 OBS, B = 2 OBS, ... 100 + A | A A WEIGHT | B A | A C ABF AA B A A | A A B C A ACBDBI AB A 50 + A A B A CDAHAAA | A B A | | | 0 + -+--------+--------+--------+--------+--------+--------+--------+ 50 60 70 80 90 100 110 120 CHEST SAS システム 8 19:52 Thursday, June 12, 2003 プロット : WEIGHT*PRED1. 凡例: A = 1 OBS, B = 2 OBS, ... 100 + A | A A WEIGHT | A A A | A A CBBBA AA CA A | A B AA B AABB AFEABB AB 50 + AA BAACABBCDA B | AA A A | | | 0 + --+---------+---------+---------+---------+---------+---------+-- 30 40 50 60 70 80 90 Predicted Value of WEIGHT SAS システム 9 19:52 Thursday, June 12, 2003 プロット : RESID1*PRED1. 凡例: A = 1 OBS, B = 2 OBS, ... | R 40 + e | s | A i 20 + A d | A A A A u | A B AA B A A BABA A A a 0 +-------------A--A--BAABAAAAA-AABB-ADDAA--A-BA--------A----------- l | A B AABCA B ABAABBAABA A | -20 + ---+---------+---------+---------+---------+---------+---------+-- 30 40 50 60 70 80 90 Predicted Value of WEIGHT SAS システム 10 19:52 Thursday, June 12, 2003 プロット : RESID1*HEIGHT. 凡例: A = 1 OBS, B = 2 OBS, ... | R 40 + e | s | A i 20 + A d | A A A A u | B A A BBBAB AAAA a 0 +------------A-A-A-A-AAA-A-B-B-BA-BB-A-CAB-BAB---B---A--A--------- l | A A AA B AA CC A BAA A ACBA A | -20 + ---+-----------+-----------+-----------+-----------+-----------+-- 140 150 160 170 180 190 HEIGHT SAS システム 11 19:52 Thursday, June 12, 2003 プロット : RESID1*CHEST. 凡例: A = 1 OBS, B = 2 OBS, ... | R 40 + e | s | A i 20 + A d | A A A A u | B A A A A B ABD B a 0 +-----------------------A---C-ACBGBBAG-A--B---B--------A---------- l | AA B CA FACBC AAB A | -20 + -+--------+--------+--------+--------+--------+--------+--------+- 50 60 70 80 90 100 110 120 CHEST SAS システム 12 19:52 Thursday, June 12, 2003 プロット : RESID1*WEIGHT. 凡例: A = 1 OBS, B = 2 OBS, ... | R 40 + e | s | A i 20 + A d | AA A A u | A B AAB B ABBA AA a 0 +----------------AA-BBBCB-B-AEDBC-A-C----A---------------- l | A A BDABB B AEAAD A | -20 + ---+------------+------------+------------+------------+-- 20 40 60 80 100 WEIGHT SAS システム 13 16:15 Tuesday, June 17, 2003 Univariate Procedure Variable=RESID1 Residual Moments N 86 Sum Wgts 86 Mean 0 Sum 0 Std Dev 6.292652 Variance 39.59747 Skewness 1.176856 Kurtosis 1.64992 USS 3365.785 CSS 3365.785 CV . Std Mean 0.678554 T:Mean=0 0 Pr>|T| 1.0000 Num ^= 0 86 Num > 0 32 M(Sign) -11 Pr>=|M| 0.0230 Sgn Rank -213.5 Pr>=|S| 0.3610 W:Normal 0.914973 Pr< W 0.0001 SAS システム 14 16:15 Tuesday, June 17, 2003 Univariate Procedure Variable=RESID1 Residual Quantiles(Def=5) 100% Max 23.52752 99% 23.52752 75% Q3 4.162278 95% 11.45481 50% Med -1.55236 90% 7.647734 25% Q1 -4.54635 10% -6.88912 0% Min -9.80661 5% -7.21872 1% -9.80661 Range 33.33413 Q3-Q1 8.708626 Mode -6.96171 SAS システム 17 16:15 Tuesday, June 17, 2003 Univariate Procedure Variable=RESID1 Residual Stem Leaf # Boxplot 2 4 1 0 1 57 2 0 1 0113 4 | 0 5556777778888 13 | 0 000111233444 12 +--+--+ -0 44433333333333222222222111111000 32 *-----* -0 998777776666666555555 21 +-----+ -1 0 1 | ----+----+----+----+----+----+-- Multiply Stem.Leaf by 10**+1 SAS システム 18 16:15 Tuesday, June 17, 2003 Univariate Procedure Variable=RESID1 Residual Normal Probability Plot 22.5+ * | * ++ | ***+*++++ | *******+ | +++***** | *********** | * * * ********** -12.5+ ++++++++ +----+----+----+----+----+----+----+----+----+----+ -2 -1 0 +1 +2

結果の見方
- 対象になったのは 86名。
- 説明変量群が予測に役立っているか?
  - 回帰に役立っているか : Prob>F : 小さいと有意
  - 「役立っている」と言える : 0.01% だから 1% で有意
- 決定係数 : R-Square ( 相関係数 : R )
  - 目的変量が説明変量でどの程度説明しているかの割合。
  - 1 に近いほど当てはまりが良いと言える。: 67.5%
  - 説明変量数が増えると大きくなるのが一般的。
- 回帰係数 : Parameter Estimate
- ある特定の説明変量が予測に役立っているか?
  - 回帰係数の検定(帰無仮説:係数=0 か?) : Prob>|T| : 小さいと有意
  - 両方とも(身長も胸囲も)有意
  - 「各係数は 0ではない」と言える : 0.01% だから 1% で有意
- 残差の性質 ===> 正規性 : 残差プロット、残差解析
  - 残差(予測誤差)は正規分布をしていると仮定してモデルが構築されている。
  - この仮定が覆ると、回帰分析として成立していないことになる。
  - 残差が正規分布をしているか確認する必要がある。
  - 均等に散らばっているか?
  - 傾向はないか? : もし傾向があると言うことになれば正規性の仮定が崩れている
  - 体重の大きい 3例が外れ値と考えられるか要確認
  - ...
- ...

特定グループでの解析

「男性のみ」と言う特定のグループに対して、同様の解析を行うには?

プログラム : les1002.sas

/* Lesson 10-2 */ /* File Name = les1002.sas 06/19/03 */ data gakusei; infile 'all03a.prn' firstobs=2; input sex $ height weight chest jitaku $ kodukai carrier $ tsuuwa; proc print data=gakusei(obs=10); run; : proc corr data=gakusei; : 相関係数 where sex='M'; : 男性について run; : : proc reg data=gakusei; : 回帰分析 where sex='M'; : 男性について model weight=height chest; : output out=outreg1 predicted=pred1 residual=resid1; : run; : proc print data=outreg1(obs=15); run; proc plot data=outreg1; where sex='M' and weight^=. and height^=. and chest^=.; : 対象データについて plot weight*height; plot weight*chest; plot weight*pred1; plot resid1*(pred1 height chest weight)/vref=0; : まとめて記述 /* plot resid1*pred1 /vref=0; plot resid1*height/vref=0; plot resid1*chest /vref=0; plot resid1*weight/vref=0; */ run; proc univariate data=outreg1 plot normal; var resid1; run;

出力結果 : les1002.lst
SAS システム 3 19:52 Thursday, June 12, 2003 Correlation Analysis Pearson Correlation Coefficients / Prob > |R| under Ho: Rho=0 / Number of Observations HEIGHT WEIGHT CHEST KODUKAI TSUUWA HEIGHT 1.00000 0.44827 0.21465 0.07149 -0.22902 0.0 0.0001 0.1026 0.3752 0.3053 168 168 59 156 22 WEIGHT 0.44827 1.00000 0.66778 -0.02635 0.00972 0.0001 0.0 0.0001 0.7441 0.9657 168 168 59 156 22 CHEST 0.21465 0.66778 1.00000 -0.12322 0.86603 0.1026 0.0001 0.0 0.3701 0.3333 59 59 59 55 3 KODUKAI 0.07149 -0.02635 -0.12322 1.00000 0.03539 0.3752 0.7441 0.3701 0.0 0.8790 156 156 55 159 21 TSUUWA -0.22902 0.00972 0.86603 0.03539 1.00000 0.3053 0.9657 0.3333 0.8790 0.0 22 22 3 21 22 SAS システム 6 19:52 Thursday, June 12, 2003 Model: MODEL1 Dependent Variable: WEIGHT Analysis of Variance Sum of Mean Source DF Squares Square F Value Prob>F Model 2 2670.74675 1335.37337 30.720 0.0001 Error 56 2434.28885 43.46944 C Total 58 5105.03559 Root MSE 6.59314 R-square 0.5232 Dep Mean 64.52034 Adj R-sq 0.5061 C.V. 10.21869 SAS システム 7 19:52 Thursday, June 12, 2003 Parameter Estimates Parameter Standard T for H0: Variable DF Estimate Error Parameter=0 Prob > |T| INTERCEP 1 -67.355995 24.32863116 -2.769 0.0076 HEIGHT 1 0.430400 0.14291427 3.012 0.0039 CHEST 1 0.651915 0.10151936 6.422 0.0001 SAS システム 10 19:52 Thursday, June 12, 2003 プロット : WEIGHT*CHEST. 凡例: A = 1 OBS, B = 2 OBS, ... WEIGHT | 100 + A | A A | 75 + AA C A A | A C AAH AA B A | A A A C A ACBCBE AA A 50 + A A A A | | 25 + -+--------+--------+--------+--------+--------+--------+--------+ 50 60 70 80 90 100 110 120 CHEST SAS システム 11 19:52 Thursday, June 12, 2003 プロット : WEIGHT*PRED1. 凡例: A = 1 OBS, B = 2 OBS, ... WEIGHT | 100 + A | A A | 75 + AA AAAA A | AA A DAAACABA B | A AA AAAAAABBDDB BA 50 + A A A A | | 25 + --+-----------+-----------+-----------+-----------+-----------+-- 40 50 60 70 80 90 Predicted Value of WEIGHT SAS システム 12 19:52 Thursday, June 12, 2003 プロット : RESID1*PRED1. 凡例: A = 1 OBS, B = 2 OBS, ... | R 40 + e | s | i 20 + A A d | A A u | A A A AA A CBA AA a 0 +---------------A---A----A-AAAAABB-AC--A-AA---------A------------- l | A A A AADBB BBB B | -20 + ---+-----------+-----------+-----------+-----------+-----------+-- 40 50 60 70 80 90 Predicted Value of WEIGHT SAS システム 13 19:52 Thursday, June 12, 2003 プロット : RESID1*HEIGHT. 凡例: A = 1 OBS, B = 2 OBS, ... | R 40 + e | s | i 20 + A A d | A A u | A B A C A B A A AA a 0 +----A-------A-----------A-B---A-B-A-A-A--AB-----AA--A-A-A---A---- l | A B A A B C A B A A A BA A A | -20 + ---+---------+---------+---------+---------+---------+---------+-- 155 160 165 170 175 180 185 HEIGHT SAS システム 14 19:52 Thursday, June 12, 2003 プロット : RESID1*CHEST. 凡例: A = 1 OBS, B = 2 OBS, ... | R 40 + e | s | i 20 + A A d | A A u | A A A A B BD B a 0 +-----------------------B---B---AB-BAE-A--B---A--------A---------- l | A B BBBAE AAB A A | -20 + -+--------+--------+--------+--------+--------+--------+--------+- 50 60 70 80 90 100 110 120 CHEST SAS システム 15 19:52 Thursday, June 12, 2003 プロット : RESID1*WEIGHT. 凡例: A = 1 OBS, B = 2 OBS, ... | R 40 + e | s | i 20 + A A d | A A u | A A A B A AC A A AA a 0 +----------A-------A---DACA-D--A-AB------A------------------------ l | A A BA B FABABA A | -20 + ---+---------+---------+---------+---------+---------+---------+-- 40 50 60 70 80 90 100 WEIGHT SAS システム 16 16:15 Tuesday, June 17, 2003 Univariate Procedure Variable=RESID1 Residual Moments N 59 Sum Wgts 59 Mean 0 Sum 0 Std Dev 6.478464 Variance 41.9705 Skewness 1.197453 Kurtosis 1.616956 USS 2434.289 CSS 2434.289 CV . Std Mean 0.843424 T:Mean=0 0 Pr>|T| 1.0000 Num ^= 0 59 Num > 0 23 M(Sign) -6.5 Pr>=|M| 0.1175 Sgn Rank -103 Pr>=|S| 0.4416 W:Normal 0.907993 Pr< W 0.0001 SAS システム 17 16:15 Tuesday, June 17, 2003 Univariate Procedure Variable=RESID1 Residual Quantiles(Def=5) 100% Max 21.70922 99% 21.70922 75% Q3 4.496601 95% 13.83164 50% Med -1.4844 90% 7.306712 25% Q1 -4.78841 10% -6.60478 0% Min -8.95245 5% -8.26289 1% -8.95245 Range 30.66167 Q3-Q1 9.285006 Mode -4.78841 SAS システム 20 16:15 Tuesday, June 17, 2003 Univariate Procedure Variable=RESID1 Residual Stem Leaf # Boxplot 2 2 1 0 1 8 1 | 1 034 3 | 0 555566777 9 | 0 011233444 9 +--+--+ -0 4444333332221111100 19 *-----* -0 99888766655555555 17 +-----+ ----+----+----+----+ Multiply Stem.Leaf by 10**+1 SAS システム 21 16:15 Tuesday, June 17, 2003 Univariate Procedure Variable=RESID1 Residual Normal Probability Plot 22.5+ * | * ++ | *+*++++++ 7.5+ ++*****+ | +++******* | *********** -7.5+ * * * ******** +----+----+----+----+----+----+----+----+----+----+ -2 -1 0 +1 +2

結果の見方
- 単変量毎の相関が有意なのは、身長と体重、体重と胸囲の間。
- 対象になったのは 59名。
- 回帰に役立っているか : 役立っている : 0.01% だから 1% で有意
- 決定係数(R-square)は 0.523
- 個々の説明変量が予測に役立っているか?
  - 係数がゼロか? : 定数項も身長も胸囲も有意(1% で有意)
- 残差の性質 ===> 正規性 : 残差プロット、残差解析
  - 均等に散らばっているか?
  - 傾向はないか? : 傾向があると言うことは正規性の仮定が崩れていること
  - 外れ値? 85Kg より重い 3名程度が吟味対象?

[演習] : 「男性のみ」で、かつ「体重の大きい 3名を除外」して実行してみよ。
- プログラム : les1003.sas、出力結果 : les1003.lst
  where sex='M' and weight<85;
- 当てはまりは良くなったか? : 異常値と外れ値の意味するもの
- 残差の正規性はどのように変化したか?
次回は、... : 06月26日 14:45
- 回帰における変量選択
- ...

講義のホームページへ戻ります