重回帰分析、変数選択と Closing＋α

統計解析 05 クラス : 第13回 (01/14/09)

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

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

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

プログラム : les1301.sas

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

出力結果 : les1301.lst
SAS システム 2 21:39 Tuesday, January 13, 2009 Model: MODEL1 Dependent Variable: TAIJYUU Analysis of Variance Sum of Mean Source DF Squares Square F Value Prob>F Model 2 9021.64522 4510.82261 93.735 0.0001 Error 116 5582.25764 48.12291 C Total 118 14603.90286 Root MSE 6.93707 R-square 0.6178 Dep Mean 59.02857 Adj R-sq 0.6112 C.V. 11.75205 SAS システム 3 21:39 Tuesday, January 13, 2009 Parameter Estimates Parameter Standard T for H0: Variable DF Estimate Error Parameter=0 Prob > |T| INTERCEP 1 -109.979309 12.42720051 -8.850 0.0001 SHINTYOU 1 0.831468 0.07722550 10.767 0.0001 KYOUI 1 0.345779 0.08166801 4.234 0.0001 SAS システム 4 21:39 Tuesday, January 13, 2009 S H T K C I A J O A T R N I K I D R S P E T J Y T U R U R S O S Y Y O A K Y U E I B E O U U K A E W D D S X U U I U I R A 1 1 1 F 145.0 38.0 . J 10000 . . . 2 F 146.7 41.0 85 J 10000 Vodafone 6000 41.3883 -0.38828 3 F 148.0 42.0 . J 50000 . . . 4 F 148.0 43.0 80 J 50000 DoCoMo 4000 40.7403 2.25971 5 F 148.9 . . J 60000 . . . 6 F 149.0 45.0 . G 60000 . . . 7 F 150.0 46.0 86 40000 . 44.4779 1.52210 8 F 150.0 . . J 10000 softbank 80 . . 9 F 151.0 45.0 . J 20000 docomo 5000 . . 10 F 151.0 50.0 . G 60000 J-PHONE . . . 11 F 151.7 41.5 80 J 35000 . 43.8167 -2.31672 12 F 152.0 35.0 77 J 60000 DoCoMo 2000 43.0288 -8.02883 13 F 152.0 43.0 . J 20000 au 3500 . . 14 F 152.0 44.0 . 45000 DoCoMo 4000 . . 15 F 153.0 41.0 . J 125000 No . . . SAS システム 6 21:39 Tuesday, January 13, 2009 プロット : TAIJYUU*SHINTYOU. 凡例: A = 1 OBS, B = 2 OBS, ... 100 + A | A A A TAIJYUU | A A A A A | B CABAB BACAA B B A AA | A A B A B BBA BAGBC ACAA AABBA 50 + A A AEB CDEAD BBACB A | A A B A A | | | 0 + --+-----------+-----------+-----------+-----------+-----------+- 140 150 160 170 180 190 SHINTYOU SAS システム 7 21:39 Tuesday, January 13, 2009 プロット : TAIJYUU*KYOUI. 凡例: A = 1 OBS, B = 2 OBS, ... 100 + A | A B TAIJYUU | A A AA A | A ACACGACCA A A | A A CCAAFBFKAAAA A 50 + A A AA FEHJCB | AA B B | | | 0 + ---+-----------+-----------+-----------+-----------+-- 40 60 80 100 120 KYOUI SAS システム 8 21:39 Tuesday, January 13, 2009 プロット : TAIJYUU*PRED1. 凡例: A = 1 OBS, B = 2 OBS, ... 100 + A | A A A TAIJYUU | A A A A A | A CBABBACAAB AAAC A | A A BA ACBBAB BCCEACA C BA 50 + B CBEBCCDBEAAC B | AA AB A | | | 0 + ---+-----------+-----------+-----------+-----------+-- 40 50 60 70 80 Predicted Value of TAIJYUU SAS システム 9 21:39 Tuesday, January 13, 2009 プロット : RESID1*PRED1. 凡例: A = 1 OBS, B = 2 OBS, ... | R 50 + e | s | A i 25 + A A d | A A u | A A A B A A A CAABB A A A a 0 +---AA--BB---CBEACCC-BABBBABAAABDBCBCCB-B---AB----A------- l | A A ABB C B A AAA B ACAA A | A -25 + ---+------------+------------+------------+------------+-- 40 50 60 70 80 Predicted Value of TAIJYUU SAS システム 10 21:39 Tuesday, January 13, 2009 プロット : RESID1*SHINTYOU. 凡例: A = 1 OBS, B = 2 OBS, ... | R 50 + e | s | A i 25 + A A d | A A u | A A B A A DABAB B A A a 0 +----------A-A-A-A-AAAEB-CDEAB-BACBB-BAGAC-CCCAA-B-A-A------------ l | A A C AA CB A A A BACAA AA | A -25 + ---+-----------+-----------+-----------+-----------+-----------+-- 140 150 160 170 180 190 SHINTYOU SAS システム 11 21:39 Tuesday, January 13, 2009 プロット : RESID1*KYOUI. 凡例: A = 1 OBS, B = 2 OBS, ... | R 50 + e | s | A i 25 + A A d | A A u | A BB C ABEA B A a 0 +----------------------B--C--GCENBGBJBBE--B--A---A-------- l | AAACBAH BACA A | A -25 + ---+------------+------------+------------+------------+-- 40 60 80 100 120 KYOUI SAS システム 12 21:39 Tuesday, January 13, 2009 プロット : RESID1*TAIJYUU. 凡例: A = 1 OBS, B = 2 OBS, ... | R 50 + e | s | A i 25 + A A d | A A u | AAAAAB D AABAA AA A a 0 +----------------BABDDFGDAE-CFEECABAE----A---------------- l | A A BE AB AAAD CAB | A -25 + ---+------------+------------+------------+------------+-- 20 40 60 80 100 TAIJYUU SAS システム 13 21:39 Tuesday, January 13, 2009 Univariate Procedure Variable=RESID1 Residual Moments N 119 Sum Wgts 119 Mean 0 Sum 0 Std Dev 6.878028 Variance 47.30727 Skewness 1.864791 Kurtosis 6.275161 USS 5582.258 CSS 5582.258 CV . Std Mean 0.630508 T:Mean=0 0 Pr>|T| 1.0000 Num ^= 0 119 Num > 0 45 M(Sign) -14.5 Pr>=|M| 0.0100 Sgn Rank -527 Pr>=|S| 0.1632 W:Normal 0.8811 Pr<W 0.0001 SAS システム 17 21:39 Tuesday, January 13, 2009 Univariate Procedure Variable=RESID1 Residual Histogram # Boxplot 35+* 1 * .* 2 * .*** 5 0 .******************* 37 +--+--+ .************************************ 72 *-----* -15+* 2 | ----+----+----+----+----+----+----+- * may represent up to 2 counts SAS システム 18 21:39 Tuesday, January 13, 2009 Univariate Procedure Variable=RESID1 Residual Normal Probability Plot 35+ * | * * | *****++++++ | +++************* | *** ******************** -15+*++*+++++++ +----+----+----+----+----+----+----+----+----+----+ -2 -1 0 +1 +2

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

特定グループでの解析

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

プログラム : les1302.sas

/* Lesson 13-2 */ /* File Name = les1302.sas 01/14/09 */ data gakusei; infile 'all08ce.prn' firstobs=2; input sex $ shintyou taijyuu kyoui jitaku $ kodukai carryer $ tsuuwa; if sex^='M' & sex^='F' then delete; : 性別不明は除外 if shintyou=. | taijyuu=. | kyoui=. then delete; : 欠損のあるデータは除外 proc print data=gakusei(obs=10); run; proc corr data=gakusei; : 相関係数 where sex='M'; : 男性について run; : : proc reg data=gakusei; : 回帰分析 model taijyuu=shintyou kyoui; : where sex='M'; : 男性について output out=outreg1 predicted=pred1 residual=resid1; : run; : proc print data=outreg1(obs=15); run; proc plot data=outreg1; where sex='M'; : 対象データについて plot taijyuu*shintyou; plot taijyuu*kyoui; plot taijyuu*pred1; plot resid1*(pred1 shintyou kyoui taijyuu)/vref=0; : まとめて記述 /* plot resid1*pred1 /vref=0; plot resid1*shintyou/vref=0; plot resid1*kyoui /vref=0; plot resid1*taijyuu /vref=0; */ run; proc univariate data=outreg1 plot normal; var resid1; run;

出力結果 : les1302.lst
SAS システム 2 21:39 Tuesday, January 13, 2009 Correlation Analysis 5 'VAR' Variables: SHINTYOU TAIJYUU KYOUI KODUKAI TSUUWA Simple Statistics Variable N Mean Std Dev Sum Minimum Maximum SHINTYOU 74 172.7 6.0782 12776.1 156.0 187.0 TAIJYUU 74 65.1581 9.2875 4821.7 46.0000 100.0 KYOUI 74 88.1757 9.6349 6525.0 46.0000 112.0 KODUKAI 69 56898.6 65702.6 3926000 0 350000 TSUUWA 17 7520.6 4554.3 127850 350.0 15000.0 SAS システム 3 21:39 Tuesday, January 13, 2009 Correlation Analysis Pearson Correlation Coefficients / Prob > |R| under Ho: Rho=0 / Number of Observations SHINTYOU TAIJYUU KYOUI KODUKAI TSUUWA SHINTYOU 1.00000 0.43492 0.19039 0.07982 -0.12917 0.0 0.0001 0.1042 0.5144 0.6212 74 74 74 69 17 TAIJYUU 0.43492 1.00000 0.40420 0.10282 0.18383 0.0001 0.0 0.0004 0.4005 0.4800 74 74 74 69 17 KYOUI 0.19039 0.40420 1.00000 -0.37398 -0.34356 0.1042 0.0004 0.0 0.0015 0.1770 74 74 74 69 17 KODUKAI 0.07982 0.10282 -0.37398 1.00000 0.61874 0.5144 0.4005 0.0015 0.0 0.0106 69 69 69 69 16 TSUUWA -0.12917 0.18383 -0.34356 0.61874 1.00000 0.6212 0.4800 0.1770 0.0106 0.0 17 17 17 16 17 SAS システム 6 21:39 Tuesday, January 13, 2009 Model: MODEL1 Dependent Variable: TAIJYUU Analysis of Variance Sum of Mean Source DF Squares Square F Value Prob>F Model 2 1865.94805 932.97402 14.950 0.0001 Error 71 4430.81209 62.40580 C Total 73 6296.76014 Root MSE 7.89973 R-square 0.2963 Dep Mean 65.15811 Adj R-sq 0.2765 C.V. 12.12395 SAS システム 7 21:39 Tuesday, January 13, 2009 Parameter Estimates Parameter Standard T for H0: Variable DF Estimate Error Parameter=0 Prob > |T| INTERCEP 1 -61.172890 26.51444015 -2.307 0.0240 SHINTYOU 1 0.567545 0.15495093 3.663 0.0005 KYOUI 1 0.321452 0.09775086 3.288 0.0016 SAS システム 10 21:39 Tuesday, January 13, 2009 プロット : TAIJYUU*SHINTYOU. 凡例: A = 1 OBS, B = 2 OBS, ... TAIJYUU | 100 + A | A A | A A 75 + AAA A A A B | B CADAA ABCAA A D A AA | A A A ABA BADBC ABAA AA BA 50 + A B A | | 25 + ---+-----------+-----------+-----------+-----------+-- 150 160 170 180 190 SHINTYOU SAS システム 11 21:39 Tuesday, January 13, 2009 プロット : TAIJYUU*KYOUI. 凡例: A = 1 OBS, B = 2 OBS, ... TAIJYUU | 100 + A | A A | A A 75 + A AA BA A A | A ACABJBCBB A | A A BCAADBEF AA A 50 + A A AA | | 25 + ---+-----------+-----------+-----------+-----------+-- 40 60 80 100 120 KYOUI SAS システム 12 21:39 Tuesday, January 13, 2009 プロット : TAIJYUU*PRED1. 凡例: A = 1 OBS, B = 2 OBS, ... TAIJYUU | 100 + A | A A | A A 75 + A A A AA A A A | A ABACB BAADA A ABBB A | B A A AA B BAACAEAAAAC A 50 + A A AA | | 25 + --+---------+---------+---------+---------+---------+---------+- 50 55 60 65 70 75 80 Predicted Value of TAIJYUU SAS システム 13 21:39 Tuesday, January 13, 2009 プロット : RESID1*PRED1. 凡例: A = 1 OBS, B = 2 OBS, ... | R 50 + e | s | A i 25 + A A d | A A u | A AAAAA A A A a 0 +------A--------A--A-A-A--B-BC-EACAAADA-AA--AC-------A------------ l | A AA A A A AAABAAAAB ACA A | A -25 + ---+---------+---------+---------+---------+---------+---------+-- 50 55 60 65 70 75 80 Predicted Value of TAIJYUU SAS システム 14 21:39 Tuesday, January 13, 2009 プロット : RESID1*SHINTYOU. 凡例: A = 1 OBS, B = 2 OBS, ... | R 50 + e | s | A i 25 + A A d | A A u | A AAA B AA A a 0 +------------------AA---ACE-AFA-ABBADAAAA-A--A-A---------- l | A B AA A BA BAA A B ACA A A | A -25 + ---+------------+------------+------------+------------+-- 150 160 170 180 190 SHINTYOU SAS システム 15 21:39 Tuesday, January 13, 2009 プロット : RESID1*KYOUI. 凡例: A = 1 OBS, B = 2 OBS, ... | R 50 + e | s | A i 25 + A A d | A A u | BD B A a 0 +------------A---------A--A-AC-AE-BAGCBD--B--A---A-------- l | A ABA C EAFA B | A -25 + ---+------------+------------+------------+------------+-- 40 60 80 100 120 KYOUI SAS システム 16 21:39 Tuesday, January 13, 2009 プロット : RESID1*TAIJYUU. 凡例: A = 1 OBS, B = 2 OBS, ... | R 50 + e | s | A i 25 + A A d | A A u | A AB A A AA A a 0 +----------------AAA---E-ED-GAAC-DB------A------------------------ l | A A A BAABAFAA BAA A | A -25 + ---+---------+---------+---------+---------+---------+---------+-- 40 50 60 70 80 90 100 TAIJYUU SAS システム 17 21:39 Tuesday, January 13, 2009 Univariate Procedure Variable=RESID1 Residual Moments N 74 Sum Wgts 74 Mean 0 Sum 0 Std Dev 7.790767 Variance 60.69606 Skewness 1.803304 Kurtosis 4.640393 USS 4430.812 CSS 4430.812 CV . Std Mean 0.905658 T:Mean=0 0 Pr>|T| 1.0000 Num ^= 0 74 Num > 0 30 M(Sign) -7 Pr>=|M| 0.1302 Sgn Rank -259.5 Pr>=|S| 0.1636 W:Normal 0.863389 Pr<W 0.0001 SAS システム 20 21:39 Tuesday, January 13, 2009 Univariate Procedure Variable=RESID1 Residual Stem Leaf # Boxplot 3 2 1 * 2 2 34 2 0 1 6 1 0 1 03 2 | 0 55567889 8 | 0 0111222222233444 16 +--+--+ -0 4444444433333333221111111 25 *-----* -0 9887776666666555 16 +-----+ -1 310 3 | ----+----+----+----+----+ Multiply Stem.Leaf by 10**+1 SAS システム 21 21:39 Tuesday, January 13, 2009 Univariate Procedure Variable=RESID1 Residual Normal Probability Plot 32.5+ * | | * * 17.5+ * ++++++ | ++*++++ | +++***** 2.5+ +++******** | ********** | * ********+ -12.5+ * *+*+++++ +----+----+----+----+----+----+----+----+----+----+ -2 -1 0 +1 +2

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

[演習1] : 「男性のみ」で、かつ「体重の大きい 3名を除外」して実行してみよ。
- プログラム : les1303.sas、出力結果 : les1303.lst
  where sex='M' and taijyuu<80;
- 当てはまりは良くなったか? : 異常値と外れ値の意味するもの
- 残差の正規性はどのように変化したか?
[演習2] : 「小遣い額」と「月平均通話料」の間に関係があるか、また、「自宅生と下宿生」に依って違いがあるか等を調べてみよ。なお、先週の例でも判るように、社会人学生と思われる者の小遣い額は他の学生と性質が異なるので、除外して考える必要があるかもしれないことに注意せよ。

回帰分析における変数選択 :

今までのモデル : 指定した説明変数全てを取り込んで計算 : フルモデル
影響力の大きい説明変数だけをモデルに取り込みたい。
たくさんの説明変数の中から影響力の大きい、いくつかの変数を選別する。
幾つかの選択方法がある : 逐次増減法(stepwise)、総当り法(rsquare)、前進選択法(forward)、後退選択法(backward)、...
例題 : アメリカの大気汚染データ : usair2.prn
- 41都市、7変量。
- SO2 濃度(y、二酸化硫黄)を他の 6変量で説明してみよう。
- 「J:\コンピュータによる統計解析05(林篤裕)\usair2.prn」に置いておく。
- 参考 : 大気汚染物質広域監視システム「そらまめ君」とその測定項目

プログラム : les1304.sas

/* Lesson 13-4 */ /* File Name = les1304.sas 01/14/09 */ data air; infile 'usair2.prn'; input id $ y x1 x2 x3 x4 x5 x6; /* label y='SO2 of air in micrograms per cubic metre' x1='Average annual temperature in F' x2='Number of manufacturing enterprises employing 20 or more workers' x3='Population size (1970 census); in thousands' x4='Average annual wind speed in miles per hour' x5='Average annual precipitation in inches' x6='Average number of days with precipitation per year' ; */ proc print data=air(obs=10); run; proc corr data=air; run; proc reg data=air; : model y=x1 x2 x3 x4 x5 x6; : フルモデル output out=outreg1 predicted=pred1 residual=resid1; : run; : proc plot data=outreg1; plot resid1*pred1 /vref=0; : plot resid1*x1 /vref=0; : ズラズラと列記 plot resid1*x2 /vref=0; : plot resid1*x3 /vref=0; : plot resid1*x4 /vref=0; : plot resid1*x5 /vref=0; : plot resid1*x6 /vref=0; : plot resid1*y /vref=0; : run; proc reg data=air; : model y=x1-x6 / selection=stepwise; : 逐次増減法 output out=outreg1 predicted=pred1 residual=resid1; : 連続変数の指定方法 run; : proc print data=outreg1(obs=15); run; proc plot data=outreg1; plot resid1*pred1 /vref=0; : plot resid1*(x1 x2 x3 x4 x5 x6) /vref=0; : 簡略形(上と比較せよ) plot resid1*(x1-x6) /vref=0; : 簡略形(これも同じ意味) plot resid1*y /vref=0; : run; proc reg data=air; : model y=x1-x6 / selection=rsquare; : 総当り法 run; :

出力結果 : les1304.lst
SAS システム 1 21:40 Tuesday, January 13, 2009 OBS ID Y X1 X2 X3 X4 X5 X6 1 Phoenix 10 70.3 213 582 6.0 7.05 36 2 Little_R 13 61.0 91 132 8.2 48.52 100 3 San_Fran 12 56.7 453 716 8.7 20.66 67 4 Denver 17 51.9 454 515 9.0 12.95 86 5 Hartford 56 49.1 412 158 9.0 43.37 127 6 Wilmingt 36 54.0 80 80 9.0 40.25 114 7 Washingt 29 57.3 434 757 9.3 38.89 111 8 Jacksonv 14 68.4 136 529 8.8 54.47 116 9 Miami 10 75.5 207 335 9.0 59.80 128 10 Atlanta 24 61.5 368 497 9.1 48.34 115 SAS システム 2 21:40 Tuesday, January 13, 2009 Correlation Analysis 7 'VAR' Variables: Y X1 X2 X3 X4 X5 X6 Simple Statistics Variable N Mean Std Dev Sum Minimum Maximum Y 41 30.0488 23.4723 1232.0 8.0000 110.0 X1 41 55.7634 7.2277 2286.3 43.5000 75.5000 X2 41 463.1 563.5 18987.0 35.0000 3344.0 X3 41 608.6 579.1 24953.0 71.0000 3369.0 X4 41 9.4439 1.4286 387.2 6.0000 12.7000 X5 41 36.7690 11.7715 1507.5 7.0500 59.8000 X6 41 113.9 26.5064 4670.0 36.0000 166.0 SAS システム 3 21:40 Tuesday, January 13, 2009 Correlation Analysis Pearson Correlation Coefficients / Prob > |R| under Ho: Rho=0 / N = 41 Y X1 X2 X3 X4 X5 X6 Y 1.00000 -0.43360 0.64477 0.49378 0.09469 0.05429 0.36956 0.0 0.0046 0.0001 0.0010 0.5559 0.7360 0.0174 X1 -0.43360 1.00000 -0.19004 -0.06268 -0.34974 0.38625 -0.43024 0.0046 0.0 0.2340 0.6970 0.0250 0.0126 0.0050 X2 0.64477 -0.19004 1.00000 0.95527 0.23795 -0.03242 0.13183 0.0001 0.2340 0.0 0.0001 0.1341 0.8405 0.4113 X3 0.49378 -0.06268 0.95527 1.00000 0.21264 -0.02612 0.04208 0.0010 0.6970 0.0001 0.0 0.1819 0.8712 0.7939 X4 0.09469 -0.34974 0.23795 0.21264 1.00000 -0.01299 0.16411 0.5559 0.0250 0.1341 0.1819 0.0 0.9357 0.3052 X5 0.05429 0.38625 -0.03242 -0.02612 -0.01299 1.00000 0.49610 0.7360 0.0126 0.8405 0.8712 0.9357 0.0 0.0010 X6 0.36956 -0.43024 0.13183 0.04208 0.16411 0.49610 1.00000 0.0174 0.0050 0.4113 0.7939 0.3052 0.0010 0.0 SAS システム 5 21:40 Tuesday, January 13, 2009 Model: MODEL1 Dependent Variable: Y Analysis of Variance Sum of Mean Source DF Squares Square F Value Prob>F Model 6 14754.63603 2459.10601 11.480 0.0001 Error 34 7283.26641 214.21372 C Total 40 22037.90244 Root MSE 14.63604 R-square 0.6695 Dep Mean 30.04878 Adj R-sq 0.6112 C.V. 48.70761 SAS システム 6 21:40 Tuesday, January 13, 2009 Parameter Estimates Parameter Standard T for H0: Variable DF Estimate Error Parameter=0 Prob > |T| INTERCEP 1 111.728481 47.31810073 2.361 0.0241 X1 1 -1.267941 0.62117952 -2.041 0.0491 X2 1 0.064918 0.01574825 4.122 0.0002 X3 1 -0.039277 0.01513274 -2.595 0.0138 X4 1 -3.181366 1.81501910 -1.753 0.0887 X5 1 0.512359 0.36275507 1.412 0.1669 X6 1 -0.052050 0.16201386 -0.321 0.7500 SAS システム 14 21:40 Tuesday, January 13, 2009 プロット : RESID1*Y. 凡例: A = 1 OBS, B = 2 OBS, ... | R 50 + A e | s | A i 25 + d | A A AA u | AA AA A A A A a 0 +------AB------AAABA-A---------A--------------------------A------- l | CAA C A | ABA A -25 + A ---+---------+---------+---------+---------+---------+---------+-- 0 20 40 60 80 100 120 Y SAS システム 15 21:40 Tuesday, January 13, 2009 Stepwise Procedure for Dependent Variable Y Step 1 Variable X2 Entered R-square = 0.41572671 C(p) = 23.10893175 DF Sum of Squares Mean Square F Prob>F Regression 1 9161.74469120 9161.74469120 27.75 0.0001 Error 39 12876.15774782 330.15789097 Total 40 22037.90243902 Parameter Standard Type II Variable Estimate Error Sum of Squares F Prob>F INTERCEP 17.61057438 3.69158676 7513.50474182 22.76 0.0001 X2 0.02685872 0.00509867 9161.74469120 27.75 0.0001 Bounds on condition number: 1, 1 ------------------------------------------------------------------------ Step 2 Variable X3 Entered R-square = 0.58632019 C(p) = 7.55859687 DF Sum of Squares Mean Square F Prob>F Regression 2 12921.26717485 6460.63358743 26.93 0.0001 Error 38 9116.63526417 239.91145432 Total 40 22037.90243902 Parameter Standard Type II Variable Estimate Error Sum of Squares F Prob>F INTERCEP 26.32508332 3.84043919 11272.71964000 46.99 0.0001 X2 0.08243410 0.01469656 7548.02378137 31.46 0.0001 X3 -0.05660660 0.01429968 3759.52248365 15.67 0.0003 Bounds on condition number: 11.43374, 45.73494 ------------------------------------------------------------------------ Step 3 Variable X6 Entered R-square = 0.61740155 C(p) = 6.36100514 DF Sum of Squares Mean Square F Prob>F Regression 3 13606.23518823 4535.41172941 19.90 0.0001 Error 37 8431.66725079 227.88289867 Total 40 22037.90243902 Parameter Standard Type II Variable Estimate Error Sum of Squares F Prob>F INTERCEP 6.96584888 11.77690656 79.72552238 0.35 0.5578 X2 0.07433399 0.01506613 5547.32153619 24.34 0.0001 X3 -0.04939437 0.01454421 2628.36952166 11.53 0.0016 X6 0.16435940 0.09480151 684.96801338 3.01 0.0913 Bounds on condition number: 12.65025, 78.63322 ------------------------------------------------------------------------ All variables left in the model are significant at the 0.1500 level. No other variable met the 0.1500 significance level for entry into the model. SAS システム 20 21:40 Tuesday, January 13, 2009 Summary of Stepwise Procedure for Dependent Variable Y Variable Number Partial Model Step Entered Removed In R**2 R**2 C(p) F Prob>F 1 X2 1 0.4157 0.4157 23.1089 27.7496 0.0001 2 X3 2 0.1706 0.5863 7.5586 15.6705 0.0003 3 X6 3 0.0311 0.6174 6.3610 3.0058 0.0913 SAS システム 21 21:40 Tuesday, January 13, 2009 OBS ID Y X1 X2 X3 X4 X5 X6 PRED1 RESID1 1 Phoenix 10 70.3 213 582 6.0 7.05 36 -0.032 10.0316 2 Little_R 13 61.0 91 132 8.2 48.52 100 23.646 -10.6461 3 San_Fran 12 56.7 453 716 8.7 20.66 67 16.285 -4.2849 4 Denver 17 51.9 454 515 9.0 12.95 86 29.410 -12.4103 5 Hartford 56 49.1 412 158 9.0 43.37 127 50.661 5.3392 6 Wilmingt 36 54.0 80 80 9.0 40.25 114 27.698 8.3020 7 Washingt 29 57.3 434 757 9.3 38.89 111 20.079 8.9208 8 Jacksonv 14 68.4 136 529 8.8 54.47 116 10.011 3.9887 9 Miami 10 75.5 207 335 9.0 59.80 128 26.844 -16.8439 10 Atlanta 24 61.5 368 497 9.1 48.34 115 28.673 -4.6731 11 Chicago 110 50.6 3344 3369 10.4 34.44 122 109.181 0.8191 12 Indianap 28 52.3 361 746 9.7 38.74 121 16.840 11.1603 13 Des_Moin 17 49.0 104 201 11.2 30.85 103 21.697 -4.6973 14 Wichita 8 56.6 125 277 12.7 30.58 82 16.053 -8.0528 15 Louisvil 30 55.6 291 593 8.3 43.11 123 19.522 10.4776 SAS システム 35 21:40 Tuesday, January 13, 2009 プロット : RESID1*Y. 凡例: A = 1 OBS, B = 2 OBS, ... 50 + A R | e | A s | AA i | A ABA A A A d 0 +--------BA-A--ABA-A-A---------A--------------------------A------- u | AC C B A A a | B A A A l | A | -50 + ---+---------+---------+---------+---------+---------+---------+-- 0 20 40 60 80 100 120 Y SAS システム 36 21:40 Tuesday, January 13, 2009 N = 41 Regression Models for Dependent Variable: Y Number in R-square Variables in Model Model 1 0.41572671 X2 1 0.24381828 X3 1 0.18800913 X1 1 0.13657727 X6 1 0.00896628 X4 1 0.00294788 X5 -------------------------- 2 0.58632019 X2 X3 2 0.51611499 X1 X2 2 0.49813569 X2 X6 2 0.42138706 X2 X5 2 0.41938296 X2 X4 ≪中略≫ 2 0.01204980 X4 X5 ----------------------------- 3 0.61740155 X2 X3 X6 3 0.61254683 X1 X2 X3 3 0.59304760 X2 X3 X5 3 0.59298732 X2 X3 X4 3 0.56222293 X1 X2 X5 ≪中略≫ 3 0.15899893 X4 X5 X6 -------------------------------- 4 0.63964257 X1 X2 X3 X5 4 0.63287070 X1 X2 X3 X4 4 0.62909408 X1 X2 X3 X6 4 0.62847667 X2 X3 X4 X6 ≪中略≫ 4 0.25499437 X1 X4 X5 X6 ----------------------------------- 5 0.66850854 X1 X2 X3 X4 X5 5 0.65012088 X1 X2 X3 X4 X6 5 0.63964824 X1 X2 X3 X5 X6 5 0.62901313 X2 X3 X4 X5 X6 5 0.60403117 X1 X2 X4 X5 X6 5 0.50433666 X1 X3 X4 X5 X6 -------------------------------------- 6 0.66951181 X1 X2 X3 X4 X5 X6 -----------------------------------------

結果の見方
- フルモデル
- 逐次選択法(stepwise)
  - 変量増減法。
  - 一度取り込まれても、組合わせによっては削除される。
- 総当り法(rsquare)
  - 説明変数の組合わせ毎の決定係数(R^2)を表示する。
  - モデルの探索用。
- 他に、前進選択法(forward)、後退選択法(backward)、...
- 「数値計算上の最適モデル」と「その分野の知識からの最適モデル」には違いがあることを知っておくこと。
- 残差解析はいつの場合でも必要
- ...
SAS の文法 : 簡略な表記
- 連続変量の指定 : x1-x6
- plot をまとめて指定 : plot resid1*(x1-x6);

[要点] 解析する上での注意点
- 残差解析を行なう
- 残差に傾向はないか? : ラッパ状、バナナカーブ、...
- 対象集団を単質のものに絞る : 男性だけ、20歳近傍だけ、...
- 外れ値を持つ個体を除外する : 標準体型の者だけ、...
- ただし、恣意的な除外は謹むべきである
- 「外れ値」と「異常値」は区別して取り扱う・識別する : 「外れて」いても「異常」かどうかは別
  - 外れ値 : 大多数のデータから大きく離れた値
  - 異常値 : 通常ではあり得ない入力ミス等の異常な値
- 上記以外に「多重共線性(multi-colinearity)」の問題もある。逆行列が取れなくなる。不安定な解。
誤用?!
- 構造がシンプルで理解し易い分、"変な"利用も散見される
- 内挿(観測点の内側)はまだしも、外挿(観測点の外側)を予測するのは難しい or 無謀。
- 予測範囲でデータの構造が一定という"条件"が必要。
- 線形で表現できる関係なのか? 線形の関係式が有効なのか? 非線型?
- 自分で判断できる能力を持つこと。疑うこと。
　[例1] 人間の成長曲線
　[例2] 将来のプログラマ必要数予測 : 21世紀(?)には国民全員がプログラマ ('80s)
　[例3] オリンピック 100m 走の男女記録 : 2156年には女性の方が速い (2004.09.30) :
　　　　　　　Japan Journal LTD の記事 , Japan Journal LTD の記事 , 朝日新聞の記事
　　　　　[究極の命題!] 100m に 0.00秒要する(!?)ようになるのは何時?
演習 : 各自のデータに対して回帰分析を行ってみよう
- フルモデル
- 残差解析
- 特定のグループ
- 変数選択
- ...
最終レポート
　半年を通して学んできた SAS の使い方、および統計手法を、自分が興味を持ったデータに適用してみて、興味深い知見を得る体験をしてもらう。
1. 対象データ :
  - 自分で収集したデータ。
  - 一つである必要はない。複数でも良い。
  - 前回用いたものと異なっていても可。
2. 作業内容 :
  1. SAS を使って解析し、興味深い知見を引き出そう。
  2. 以下の点に注意しながらレポートを作成しよう。
  3. 利用するデータ解析手法については、特に制限や指定をしないが、「多変量解析の手法」を使うとより高度なデータ構造が把握できることがある。
  4. 前回のレポートの反省点を踏まえて作成すること。
  5. 興味を持った点や得られた知見に対する考察は人によって個々異なるものであるので、他人と相談することなく自分の力で解析しレポートを作成すること。
3. レポート : 以下に挙げるような項目を含めて作成すること。
  - 所属学部名、学籍番号、氏名
  - データ内容の説明
  - どのような点に興味を持ったか
  - 自分の解析目的
  - 何を知りたいためにどのような手法を使ったのか
  - 得られた知見と考察
  - その他、気付いたこと
  - 講義全体を通しての感想 : 今後の参考にしたいので
4. 提出期限 :
  2009年01月30日(金) 15:00まで
  注意 :
  1. 紙で提出する場合は、事務所の受付終了時刻に注意すること。提出日は事務室の受領印で判断する。
  2. 電子メールで提出する場合に、添付ファイルは使わないこと。また、提出日時はメールヘッダーから判断する。受領確認メールを必ず返すのでこれを受け取って提出完了となる。
  3. レポートを受領した者の学籍番号は、講義の連絡ページに掲載するので、確認すること。ただし、2回とも提出したからと言って単位が認定されるわけではない点には注意されたい。
データやプログラムのバックアップ
　講義を終えると(年度末まで?)、 stat システムにログインできなくなり、それと同時に stat システム内に保存してあるデータやプログラムも呼び出せなくなる(正確には消去されてしまう)。この半年間の勉強成果を残しておきたい人は、Windows 側に転送して、 FD や MO、USB メモリー等に早めにバックアップを取るようにして下さい。
最後に
　この講義を通して、「統計」や「データ解析」と言う言葉に多少なりとも親しみを持っていただけただろうか? 統計手法については、数式よりもその手法の考え方や利用目的に重点をおいて説明したつもりである。また、それらを計算する"道具"として SAS を使った。
　今後、新聞や雑誌と言った生活では勿論のこと、いろいろな場面で、種々の数値列に出会うことになると思うが、提示された数値にはどの様な意味(と意図)があり、どう理解して、個々人としてどうアクションを起すかの、一つの判断手段として活用してもらえれば幸いである。
　なお、今まで紹介していた私のメールアドレスは実は講義用のものであった。今後、もし統計に関して何か疑問に出会い、私に連絡・相談してみたいと思った時は、以下のアドレスを使ってください。
メールアドレス : hayashi@rd.dnc.ac.jp
　皆さんの期待に応えられたか心許無い部分もありますが、半年間ご苦労様でした。
次回、次々回は、... : 01月21日、01月28日
- 自習。レポート作成。
[番外編] 主成分分析と因子分析: 両者の違いを中心に
　いくつか(p個)の変量の値を情報の損失をできるだけ少なくして、少数変量(m個、m＜p)の総合的指標(主成分)で代表させる方法として主成分分析(Principal Component Analysis, PCA)と因子分析(Factor Analysis, FA)がある。いくつかのテストの成績を総合した総合的成績、いろいろな症状を総合した総合的な重症度、種々の財務指標に基づく企業の評価等を求めたいといった場合に用いられる。 p変量(p次元)の観測値をm個(m次元)の主成分に縮約させるという意味で、次元を減少させる(reduce)方法と言うこともでき、多変量データを要約する有力な方法である。
　両者は似た目的に使われるが、元になっている考え方は異なるので利用する場面では注意が必要である。違いに焦点を当てながら説明する。
- 資料 : 東工大能力測定法特論 (第5回, 11/28/08)
主成分分析(PCA)と因子分析(FA): 目的は同じでも異なる手法
- 考え方(基になっているモデル)
  - 主成分分析 : 分散最大
  - 因子分析 : 構造モデルの導入、回転性、一意性に疑問(特に斜交回転)
- 対象データ
  - サンプルサイズのある程度大きいデータが対象となる
  - 変量数もある程度大きいデータが対象となる
  - 当然ながらサンプルサイズの方が変量数よりも大きいこと : 多変量解析全般
- 対象データを熟知している方が解釈しやすい(熟知の必要性)
- 因子の特徴付けはデータのバックグラウンドに深く関係
- 経験を積むとより納得する説明ができる
- 潜在的な構造が仮定できるか? モデルが適用可能な問題か吟味する必要性。
- 単位系を変えて解析してみる。: 主成分分析
- 軸数、因子数を変えて解釈を行ってみる。行きつ戻りつして試行錯誤してみる。
  [裏技] ラインマーカーで絶対値の大きい変量にマークを付ける
  いろいろなデータで経験を積んでみる。
  どっちを使う? : やってみる。解釈してみる。今までの事例と比較してみる。
  ...
  
  データによっては解釈が困難なことも有り得る。また、自分の思い付かない結果を含んでいることもある。

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

重回帰分析、変数選択 と Closing＋α

統計解析 05 クラス : 第13回 (01/14/09)

重回帰分析、変数選択と Closing＋α