多変量解析のプロシジャー、まとめ

統計モデル解析特論I/II : 第14回 (01/25/22)

　最終回の今回は、Q3で紹介した多変量解析の手法について、SASでの実行方法を示すと共に、その他 SAS を利用する上で参考になるサイトを紹介し、最後にまとめを行う。

前回のショート課題: 8名
- a. 演習4(や余力があれば他も)を自力で実行できたであろうか?
  - 実行できました。x6名
  - 無事、実行できました。
    事前に資料をいただけたので、予習が捗りました。
  - 実行することができた．
    しかし，結果を見てどう判断すればよいか少々迷った．自宅生と下宿生で身長に差があるかを調べたとき，正規分布であることは確認できたが，等分散であるかどうかは閾値を1%と5%の間になった．等分散であるとして（閾値5%とした場合）t検定を実行すると2群間に差があると出たが，等分散でないとした場合だと差がないという結果となった．こういう場合はどちらの結果を採用するべきか教えてほしい．（また，wilcoxonの順位和検定の結果は2群間に差がないとなった．）

回帰分析: proc reg

単回帰分析: Q3 第04回 .
重回帰分析: Q3 第04回 .
過去のデータからその構造を把握し、新規に測定されたデータに対する予測を行ないたい。
説明変量の線形結合で目的変量を推定する。
残差は測定誤差であるという考え方に基づいており、目的変量に対して平行に取る。
残差は正規分布をしていると仮定して理論が成り立っている。
体重を、身長と胸囲から予測する等。
- 残差の正規性の確認が必要。

プログラム : Lesson 14-01 : les1401.sas <=== Original File Name = Reg0403.sas
/* Lesson 14-01 */ /* File Name = les1401.sas 01/25/22 */ /* Original File File Name = Reg0403.sas 10/24/21 */ options nocenter linesize=78 pagesize=30; options locale='en_US'; proc printto print = 'StatM21/les1401-Results.txt' new; ods listing gpath='StatM21/SAS_ODS14'; data gakusei; infile 'StatM21/StudAll21d.csv' firstobs=9 dlm=',' dsd missover encoding=sjis termstr=crlf; input sex $ shintyou taijyuu kyoui jitaku : $10. kodukai carryer $ tsuuwa; if shintyou='.' or taijyuu='.' or kyoui='.' then delete; if kyoui<60 then delete; if taijyuu>85 then delete; proc print data=gakusei(obs=10); run; proc reg data=gakusei; : 回帰分析 model taijyuu=shintyou kyoui; : 複数変量を指定(重回帰) output out=out_reg2 predicted=pred2 residual=resid2; : 結果項目の保存 run; : 予測値(pred2)と残差(resid2)をout_reg2に保存 proc print data=out_reg2(obs=15); : 残差を表示 run; : proc plot data=out_reg2; plot taijyuu*shintyou/vaxis=20 to 100 by 20; plot taijyuu*kyoui/vaxis=20 to 100 by 20; plot pred2*taijyuu; plot resid2*shintyou/vref=0; plot resid2*kyoui /vref=0; plot resid2*taijyuu /vref=0; run; proc univariate data=out_reg2 plot normal; var resid2; run;

出力
- [結果タグ] html 形式: les1401-Results.html
- [結果タグ] pdf 形式: les1401-Output.pdf
- テキストファイル: les1401-Results.txt

回帰分析における変数選択 : 説明変数の取捨選択(変数選択)について説明する。

今までのモデル : 指定した説明変数全てを取り込んで計算 : 「フルモデル」と呼ぶ

影響力の大きい説明変数だけをモデルに取り込みたい。

たくさんの説明変数の中から影響力の大きい、いくつかの変数を選別する。

幾つかの選択方法がある :

逐次増減法(stepwise)
前進選択法(forward)
後退選択法(backward)

総当たり法(rsquare)

全部の組合わせを概観してみたい : 総当たり
モデルの探索用 : どの組合わせを使うと説明力が大きいか

例題 : アメリカの大気汚染データ : usair2.txt
- 41都市、7変数。
- 7変量 : 都市名, SO2(二酸化硫黄)濃度, 気温, 製造業数, 人口, 風速, 降雨量, 降雨日数
- SO2 濃度を他の 6変数で説明してみよう。
- 6変数の中には影響力が大きいものがあるはずだ。<==> 影響力が小さいもの
- 参考 : 大気汚染物質広域監視システム「そらまめ君」とその測定項目

プログラム : les1402.sas

/* Lesson 14-02 */ /* File Name = les1402.sas 01/25/22 */ options nocenter linesize=78 pagesize=30; options locale='en_US'; /* options locale='ja_JP'; */ proc printto print = 'StatM21/les1402-Results.txt' new; ods listing gpath='StatM21/SAS_ODS14'; data air; infile 'StatM21/usair2.txt'; 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; title "*** フルモデル ***"; proc reg data=air; : model y=x1 x2 x3 x4 x5 x6; : フルモデル output out=out_reg1 predicted=pred1 residual=resid1; : run; : proc plot data=out_reg1; plot resid1*pred1 /vref=0; : plot resid1*x1 /vref=0; : 指定方法1(ズラズラと列記) 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; title "*** 逐次増減法 ***"; proc reg data=air; : model y=x1-x6 / selection=stepwise; : 逐次増減法 output out=out_reg2 predicted=pred2 residual=resid2; : 連続変数の指定方法(列挙が面倒な場合) run; : proc print data=out_reg2(obs=15); run; proc plot data=out_reg2; plot resid1*pred1 /vref=0; : plot resid1*(x1 x2 x3 x4 x5 x6) /vref=0; : 指定方法2(上と比較) plot resid1*(x1-x6) /vref=0; : 指定方法3(何れも同じ目的の指定) plot resid1*y /vref=0; : run; title "*** 総当たり法 ***"; proc reg data=air; : model y=x1-x6 / selection=rsquare; : 総当たり法 run; :

出力結果:

[結果タグ] html 形式: les1402-Results.html
[結果タグ] pdf 形式: les1402-Output.pdf
テキストファイル: les1402-Results.txt

Thursday, January 20, 2022 05:35:01 PM 85 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 Thursday, January 20, 2022 05:35:01 PM 86 The CORR Procedure 7 Variables: y x1 x2 x3 x4 x5 x6 Simple Statistics Variable N Mean Std Dev Sum Minimum Maximum y 41 30.04878 23.47227 1232 8.00000 110.00000 x1 41 55.76341 7.22772 2286 43.50000 75.50000 x2 41 463.09756 563.47395 18987 35.00000 3344 x3 41 608.60976 579.11302 24953 71.00000 3369 x4 41 9.44390 1.42864 387.20000 6.00000 12.70000 x5 41 36.76902 11.77155 1508 7.05000 59.80000 x6 41 113.90244 26.50642 4670 36.00000 166.00000 Thursday, January 20, 2022 05:35:02 PM 87 The CORR Procedure Pearson Correlation Coefficients, N = 41 Prob > |r| under H0: Rho=0 y x1 x2 x3 x4 x5 x6 y 1.00000 -0.43360 0.64477 0.49378 0.09469 0.05429 0.36956 0.0046 <.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.2340 0.6970 0.0250 0.0126 0.0050 x2 0.64477 -0.19004 1.00000 0.95527 0.23795 -0.03242 0.13183 <.0001 0.2340 <.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 <.0001 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.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.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 *** フルモデル *** Thursday, January 20, 2022 05:35:02 PM 88 The REG Procedure Model: MODEL1 Dependent Variable: y Number of Observations Read 41 Number of Observations Used 41 Analysis of Variance Sum of Mean Source DF Squares Square F Value Pr > F Model 6 14755 2459.10601 11.48 <.0001 Error 34 7283.26641 214.21372 Corrected Total 40 22038 Root MSE 14.63604 R-Square 0.6695 Dependent Mean 30.04878 Adj R-Sq 0.6112 Coeff Var 48.70761 *** フルモデル *** Thursday, January 20, 2022 05:35:02 PM 89 The REG Procedure Model: MODEL1 Dependent Variable: y Parameter Estimates Parameter Standard Variable DF Estimate Error t Value Pr > |t| Intercept 1 111.72848 47.31810 2.36 0.0241 x1 1 -1.26794 0.62118 -2.04 0.0491 x2 1 0.06492 0.01575 4.12 0.0002 x3 1 -0.03928 0.01513 -2.60 0.0138 x4 1 -3.18137 1.81502 -1.75 0.0887 x5 1 0.51236 0.36276 1.41 0.1669 x6 1 -0.05205 0.16201 -0.32 0.7500 *** フルモデル *** Thursday, January 20, 2022 05:35:02 PM 90 Plot of resid1*pred1. Legend: A = 1 obs, B = 2 obs, etc. | 50 + A | | | | R | A e 25 + s | i | A A d | A A u | A A A A a | AA AAA A A l 0 +--------------A-A-------A---A---------A--------------------------A----- | A A BA | A A A | AA B | AB A A | -25 + A -+---------+---------+---------+---------+---------+---------+---------+ -20 0 20 40 60 80 100 120 Predicted Value of y *** フルモデル *** Thursday, January 20, 2022 05:35:02 PM 91 Plot of resid1*x1. Legend: A = 1 obs, B = 2 obs, etc. | 50 + A | | | | R | A e 25 + s | i | A A d | A A u | A A A A a | A A A A A A A l 0 +-----------------AAAA--------A-----------------------------A--------- | B A A A | A A A | A A A A | A A A A A | -25 + A ---+-------+-------+-------+-------+-------+-------+-------+-------+-- 40 45 50 55 60 65 70 75 80 x1 *** フルモデル *** Thursday, January 20, 2022 05:35:02 PM 92 Plot of resid1*x2. Legend: A = 1 obs, B = 2 obs, etc. | 50 + A | | | | R | A e 25 + s | i | A A d | A A u | A A A A a | AAA AA A A l 0 +---AAA--A------------A---------------------------------------------A--- | A B AA | A A A | A AA A | A B A A | -25 + A -+---------+---------+---------+---------+---------+---------+---------+ 0 500 1000 1500 2000 2500 3000 3500 x2 *** フルモデル *** Thursday, January 20, 2022 05:35:02 PM 93 Plot of resid1*x3. Legend: A = 1 obs, B = 2 obs, etc. | 50 + A | | | | R | A e 25 + s | i | A A d | A A u | AA A A a | A AA A B A l 0 +---A--AA---A------------------A------------------------------------A--- | A AA AA | A A A | A AA A | A AA AA | -25 + A -+---------+---------+---------+---------+---------+---------+---------+ 0 500 1000 1500 2000 2500 3000 3500 x3 *** フルモデル *** Thursday, January 20, 2022 05:35:02 PM 94 Plot of resid1*x4. Legend: A = 1 obs, B = 2 obs, etc. | 50 + A | | | | R | A e 25 + s | i | A A d |A A u | A A AA a | A A AA A AA l 0 +--------------------------A---B----------A--A----------------------A--- | A A A A A | A A A | A A A A | A A A A A | -25 + A -+---------+---------+---------+---------+---------+---------+---------+ 6 7 8 9 10 11 12 13 x4 *** フルモデル *** Thursday, January 20, 2022 05:35:02 PM 95 Plot of resid1*x5. Legend: A = 1 obs, B = 2 obs, etc. | 50 + A | | | | R | A e 25 + s | i | A A d | A A u | A AA A a | A A A AA A A l 0 +------------------------------------B---A--A------A-----------------A-- | A A A A A | A A A | A A A A | A A AA A | -25 + A ---+----------+----------+----------+----------+----------+----------+-- 0 10 20 30 40 50 60 x5 *** フルモデル *** Thursday, January 20, 2022 05:35:02 PM 96 Plot of resid1*x6. Legend: A = 1 obs, B = 2 obs, etc. | 50 + A | | | | R | A e 25 + s | i | A A d | A A u | A A B a | A A A A AA A l 0 +---------------------------A---------------A-BA-A-------------------- | A B A A | A A A | A A A A | B A A A | -25 + A ---+-------+-------+-------+-------+-------+-------+-------+-------+-- 20 40 60 80 100 120 140 160 180 x6 *** フルモデル *** Thursday, January 20, 2022 05:35:02 PM 97 Plot of resid1*y. Legend: A = 1 obs, B = 2 obs, etc. | 50 + A | | | | R | A e 25 + s | i | A A d | A A u | A A A A a | AA BA A A l 0 +------AA-------A----A---------A--------------------------A------- | A AA AA | AA A | B B | ABA A | -25 + A ---+---------+---------+---------+---------+---------+---------+-- 0 20 40 60 80 100 120 y *** 逐次増減法 *** Thursday, January 20, 2022 05:35:02 PM 98 The REG Procedure Model: MODEL1 Dependent Variable: y Number of Observations Read 41 Number of Observations Used 41 Stepwise Selection: Step 1 Variable x2 Entered: R-Square = 0.4157 and C(p) = 23.1089 Analysis of Variance Sum of Mean Source DF Squares Square F Value Pr > F Model 1 9161.74469 9161.74469 27.75 <.0001 Error 39 12876 330.15789 Corrected Total 40 22038 *** 逐次増減法 *** Thursday, January 20, 2022 05:35:02 PM 99 The REG Procedure Model: MODEL1 Dependent Variable: y Stepwise Selection: Step 1 Parameter Standard Variable Estimate Error Type II SS F Value Pr > F Intercept 17.61057 3.69159 7513.50474 22.76 <.0001 x2 0.02686 0.00510 9161.74469 27.75 <.0001 Bounds on condition number: 1, 1 ------------------------------------------------------------------------------ Stepwise Selection: Step 2 Variable x3 Entered: R-Square = 0.5863 and C(p) = 7.5586 *** 逐次増減法 *** Thursday, January 20, 2022 05:35:02 PM 100 The REG Procedure Model: MODEL1 Dependent Variable: y Stepwise Selection: Step 2 Analysis of Variance Sum of Mean Source DF Squares Square F Value Pr > F Model 2 12921 6460.63359 26.93 <.0001 Error 38 9116.63526 239.91145 Corrected Total 40 22038 Parameter Standard Variable Estimate Error Type II SS F Value Pr > F Intercept 26.32508 3.84044 11273 46.99 <.0001 x2 0.08243 0.01470 7548.02378 31.46 <.0001 x3 -0.05661 0.01430 3759.52248 15.67 0.0003 Bounds on condition number: 11.434, 45.735 ------------------------------------------------------------------------------ Stepwise Selection: Step 3 *** 逐次増減法 *** Thursday, January 20, 2022 05:35:02 PM 101 The REG Procedure Model: MODEL1 Dependent Variable: y Stepwise Selection: Step 3 Variable x6 Entered: R-Square = 0.6174 and C(p) = 6.3610 Analysis of Variance Sum of Mean Source DF Squares Square F Value Pr > F Model 3 13606 4535.41173 19.90 <.0001 Error 37 8431.66725 227.88290 Corrected Total 40 22038 Parameter Standard Variable Estimate Error Type II SS F Value Pr > F Intercept 6.96585 11.77691 79.72552 0.35 0.5578 x2 0.07433 0.01507 5547.32154 24.34 <.0001 x3 -0.04939 0.01454 2628.36952 11.53 0.0016 x6 0.16436 0.09480 684.96801 3.01 0.0913 *** 逐次増減法 *** Thursday, January 20, 2022 05:35:02 PM 102 The REG Procedure Model: MODEL1 Dependent Variable: y Stepwise Selection: Step 3 Bounds on condition number: 12.65, 78.633 ------------------------------------------------------------------------------ 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. Summary of Stepwise Selection Variable Variable Number Partial Model Step Entered Removed Vars In R-Square R-Square C(p) F Value Pr > F 1 x2 1 0.4157 0.4157 23.1089 27.75 <.0001 2 x3 2 0.1706 0.5863 7.5586 15.67 0.0003 3 x6 3 0.0311 0.6174 6.3610 3.01 0.0913 *** 逐次増減法 *** Thursday, January 20, 2022 05:35:03 PM 103 Obs id y x1 x2 x3 x4 x5 x6 pred2 resid2 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 *** 逐次増減法 *** Thursday, January 20, 2022 05:35:03 PM 104 Plot of resid2*pred2. Legend: A = 1 obs, B = 2 obs, etc. | | 50 + A | | | A R 25 + e | B s | A AAA A i | A A AAA A B d 0 +----------------A--AA-----------------------------------------A-------- u | CAAA A A A A a | B AA l | AA A A A -25 + | A | | -50 + | ---+----------+----------+----------+----------+----------+----------+-- 0 20 40 60 80 100 120 Predicted Value of y *** 逐次増減法 *** Thursday, January 20, 2022 05:35:03 PM 105 Plot of resid2*x1. Legend: A = 1 obs, B = 2 obs, etc. | | 50 + A | | | A R 25 + e | A A s | A A A A A i | A B A A AA A d 0 +------------------AA------A---A-------------------------------------- u | A A B B A AA A a | A A A A l | A A A A A -25 + | A | | -50 + | ---+-------+-------+-------+-------+-------+-------+-------+-------+-- 40 45 50 55 60 65 70 75 80 x1 *** 逐次増減法 *** Thursday, January 20, 2022 05:35:03 PM 106 Plot of resid2*x2. Legend: A = 1 obs, B = 2 obs, etc. | | 50 + A | | | A R 25 + e | A A s | A A AA A i | AAB A AA A d 0 +-A--A----------------A---------------------------------------------A--- u | AAA ABAA A A a | A A AA l | A A A AA -25 + | A | | -50 + | -+---------+---------+---------+---------+---------+---------+---------+ 0 500 1000 1500 2000 2500 3000 3500 x2 *** 逐次増減法 *** Thursday, January 20, 2022 05:35:03 PM 107 Plot of resid2*x3. Legend: A = 1 obs, B = 2 obs, etc. | | 50 + A | | | A R 25 + e | A A s | A B A A i | AAAA BA A d 0 +-A----A-----------------------A------------------------------------A--- u | A AA AAAA AA A a | A A B l | A A AA A -25 + | A | | -50 + | -+---------+---------+---------+---------+---------+---------+---------+ 0 500 1000 1500 2000 2500 3000 3500 x3 *** 逐次増減法 *** Thursday, January 20, 2022 05:35:03 PM 108 Plot of resid2*x4. Legend: A = 1 obs, B = 2 obs, etc. | | 50 + A | | | A R 25 + e | A A s |A A AA A i | AAAAB A A d 0 +-----A----------A------------------------A--A-------------------------- u | A A AA A AB A A a | A A A A l | A A A A A -25 + | A | | -50 + | -+---------+---------+---------+---------+---------+---------+---------+ 6 7 8 9 10 11 12 13 x4 *** 逐次増減法 *** Thursday, January 20, 2022 05:35:03 PM 109 Plot of resid2*x5. Legend: A = 1 obs, B = 2 obs, etc. | | 50 + A | | | A R 25 + e | A A s | A AA A A i | A A A A AA A A d 0 +------------------------------------A---A------A-A--------------------- u | A AB A A A BA a | A A A A l | A A A A A -25 + | A | | -50 + | ---+----------+----------+----------+----------+----------+----------+-- 0 10 20 30 40 50 60 x5 *** 逐次増減法 *** Thursday, January 20, 2022 05:35:03 PM 110 Plot of resid2*x6. Legend: A = 1 obs, B = 2 obs, etc. | | 50 + A | | | A R 25 + e | A A s | A B AA i | A A A A B A A d 0 +----------------------------------------A--A--A------A--------------- u | A A A BA A A A A a | A B A l | A A A A A -25 + | A | | -50 + | ---+-------+-------+-------+-------+-------+-------+-------+-------+-- 20 40 60 80 100 120 140 160 180 x6 ≪中略≫ *** 逐次増減法 *** Thursday, January 20, 2022 05:35:03 PM 117 Plot of resid2*y. Legend: A = 1 obs, B = 2 obs, etc. | | 50 + A | | | A R 25 + e | AA s | A AAA A i | AA AAA A B d 0 +---------------A--A-A------------------------------------A------- u | ABAA B A A A a | A B A l | B A A A -25 + | A | | -50 + | ---+---------+---------+---------+---------+---------+---------+-- 0 20 40 60 80 100 120 y *** 総当たり法 *** Thursday, January 20, 2022 05:35:03 PM 118 The REG Procedure Model: MODEL1 Dependent Variable: y R-Square Selection Method Number of Observations Read 41 Number of Observations Used 41 Number in Model R-Square Variables in Model 1 0.4157 x2 1 0.2438 x3 1 0.1880 x1 1 0.1366 x6 1 0.0090 x4 1 0.0029 x5 ------------------------------------------- 2 0.5863 x2 x3 2 0.5161 x1 x2 2 0.4981 x2 x6 2 0.4214 x2 x5 2 0.4194 x2 x4 2 0.4066 x1 x3 2 0.3657 x3 x6 2 0.2483 x3 x5 2 0.2458 x1 x5 2 0.2439 x3 x4 2 0.2291 x1 x6 2 0.1917 x1 x4 2 0.1587 x5 x6 2 0.1378 x4 x6 2 0.0120 x4 x5 ------------------------------------------- 3 0.6174 x2 x3 x6 3 0.6125 x1 x2 x3 3 0.5930 x2 x3 x5 3 0.5930 x2 x3 x4 3 0.5622 x1 x2 x5 3 0.5452 x1 x2 x6 3 0.5452 x1 x2 x4 3 0.5083 x2 x4 x6 3 0.5047 x2 x5 x6 3 0.4649 x1 x3 x5 3 0.4446 x1 x3 x6 3 0.4320 x1 x3 x4 3 0.4250 x2 x4 x5 3 0.3808 x3 x5 x6 3 0.3702 x3 x4 x6 3 0.2550 x1 x4 x5 3 0.2484 x3 x4 x5 3 0.2462 x1 x5 x6 3 0.2332 x1 x4 x6 3 0.1590 x4 x5 x6 ------------------------------------------- 4 0.6396 x1 x2 x3 x5 4 0.6329 x1 x2 x3 x4 4 0.6291 x1 x2 x3 x6 4 0.6285 x2 x3 x4 x6 4 0.6176 x2 x3 x5 x6 4 0.6028 x1 x2 x4 x5 4 0.5997 x2 x3 x4 x5 4 0.5747 x1 x2 x4 x6 4 0.5622 x1 x2 x5 x6 4 0.5164 x2 x4 x5 x6 4 0.5035 x1 x3 x4 x5 4 0.4708 x1 x3 x4 x6 4 0.4649 x1 x3 x5 x6 4 0.3871 x3 x4 x5 x6 4 0.2550 x1 x4 x5 x6 ------------------------------------------- 5 0.6685 x1 x2 x3 x4 x5 5 0.6501 x1 x2 x3 x4 x6 5 0.6396 x1 x2 x3 x5 x6 5 0.6290 x2 x3 x4 x5 x6 5 0.6040 x1 x2 x4 x5 x6 5 0.5043 x1 x3 x4 x5 x6 ------------------------------------------- 6 0.6695 x1 x2 x3 x4 x5 x6

結果の見方
- 相関係数(Correlation Coefficient)とゼロと仮定した時の検定結果
- フルモデル
- 逐次選択法(stepwise)
  - 変量増減法。
  - 一度取り込まれても、組合わせによっては削除される。
- 総当たり法(rsquare)
  - 説明変数の組合わせ毎の決定係数(R^2)が表示される : 大きい順に
  - モデルの探索用。
  - 決定係数 : R-Square : 1 に近いほど当てはまりが良いと言える
  - 説明変数が増えると大きくなるのが一般的
  - 興味のある組合わせを見つけ出して、このあと計算させる。残差解析も行うこと。
- 他に、前進選択法(forward)、後退選択法(backward)
- 「数値計算上の最適モデル」と「その分野の知識からの最適モデル」には違いがあることを知っておくこと。
- 残差解析はいつの場合でも必要
  - 残差の性質 ===> 正規性 : 残差プロット、残差解析
  - 均等に散らばっているか?
  - 傾向はないか? : 傾向があると言うことは正規性の仮定が崩れていること
  - ...
- ...

主成分分析: proc princomp
- 主成分分析: 第06回 .
- いくつか(p個)の変量の値を情報の損失をできるだけ少なくして、少数変量(m個、m＜p)の総合的指標(主成分)で代表させる方法。
- いくつかのテストの成績を総合した総合的成績、いろいろな症状を総合した総合的な重症度、種々の財務指標に基づく企業の評価等を求めたい。
- 多次元のデータを少ない次元で表現したい。
- 測定された変量の分散を最大となる軸を見つける。
- プログラム : Lesson 14-03 : les1403.sas <=== Original File Name = PCA0503.sas
  proc princomp data=gakusei out=out_prin; : 相関係数を使って var shintyou taijyuu kyoui; : run; : proc print data=out_prin(obs=15); : 主成分スコア(得点)を表示 run; :
- 出力
  - [結果タグ] html 形式: les1403-Results.html
  - [結果タグ] pdf 形式: les1403-Output.pdf
  - テキストファイル: les1403-Results.txt
因子分析: proc factor
- 因子分析: 第06回 .
- いくつか(p個)の変量の値を情報の損失をできるだけ少なくして、少数変量(m個、m＜p)の総合的指標(主成分)で代表させる方法。
- いくつかのテストの成績を総合した総合的成績、いろいろな症状を総合した総合的な重症度、種々の財務指標に基づく企業の評価等を求めたい。
- 多次元のデータを少ない次元で表現したい。
- 測定された変量が潜在変量の積として表現できるとのモデルを導入して構成されている。
- 因子数を指定する必要がある。因子数に依って解釈が異なることもある。
- 回転の不定性が存在する。
- 主成分分析と因子分析は背景としている構造が異なる。
- プログラム : Lesson 14-04 : les1404.sas <=== Original File Name = FA0613.sas
  proc factor data=food nfactor=3 rotate=varimax out=out_fscore; : 因子数を3、バリマックス回転を指定 var X01-X10; : run; : proc print data=out_fscore; : 因子スコア(得点)を表示 run; :
- 出力
  - [結果タグ] html 形式: les1404-Results.html
  - [結果タグ] pdf 形式: les1404-Output.pdf
  - テキストファイル: les1404-Results.txt
全体を通してのデータ分析の手順
　Q3では主に統計手法のアイディアや考え方を、Q4ではそれらを計算する道具立てとしてのSASを習得した。改めて、データを分析する際の、問題設定から報告までの全体的な流れについて紹介しておく。
1. 取り組む問題の理解
2. データ採取計画の立案 <=== 設計は重要
3. データ収集
  - 試しに採ってみる。気付かなかった齟齬はないか。
  - 正式に大規模に採取。
4. データクリーニング: 回答ミス、入力ミス等を洗い出す
5. 基礎統計量
6. 各種分析: 多変量解析等
7. 考察
8. 報告
- 考慮する点: 一連の行動を行う際に考慮する点は以下のようなものがある。
  - 素データ・電子化データに間違いがないかの確認。「人間は間違う動物である」。
  - データの生成背景を知る。自分で測定してみると判ることも多い。
  - 他人が採取したデータは要注意。熟知する必要がある。
  - 自分の経験と照らし合わせて矛盾はないか。

お役立ちサイト: サンプルデータ

いろいろと経験を積む際のサンプルとして役立つ
データを収集する際にも有効ではないか
いくつかのサンプルデータがSASシステム内に保存されている
データの説明やその使い方、プログラムも
- Sashelp Data Sets - SAS Support (PDF ファイル)
- SAS/STAT 9.3 User's Guide
  - 左側の項目: Sashelp Data Sets

利用例 : Lesson 14-05 : les1405.sas

Fisher's Iris Data (1936): 多変量解析では非常に有名なデータ
Fisher, R. A. が判別分析の例示データとして提示
Iris flower data set : by Wikipedia
アイリス(あやめ)のデータ。3品種(Setosa, Versicolor, Virginica)。 150サンプル(各品種50サンプル)。萼(がく, Sepal)の幅と長さ、花弁(Petal)の幅と長さの4変量から品種を分類する例。
SAS に保存されている。ファイル名＝sashelp.iris

/* Lesson 14-05 */ /* File Name = les1405.sas 01/25/22 */ options nocenter linesize=78 pagesize=30; options locale='en_US'; /* options locale='ja_JP'; */ proc printto print = 'StatM21/les1405-Results.txt' new; ods listing gpath='StatM21/SAS_ODS14'; title "Sashelp.iris --- Fisher's Iris Data (1936)"; proc contents data=sashelp.iris varnum; : データの変量情報を表示する ods select position; : データの指定方法にも注目 run; title "The First Five Observations Out of 150"; proc print data=sashelp.iris(obs=5) noobs; : 先頭5サンプルを表示 run; title "The Species Variable"; proc freq data=sashelp.iris; : 頻度集計 tables Species; run; proc fastclus data=sashelp.iris out=out_clust maxclusters=3; : クラスター分析 var SepalLength SepalWidth PetalLength PetalWidth; run; proc plot data=out_clust; plot SepalLength*SepalWidth=cluster; plot SepalLength*PetalLength=cluster; plot SepalLength*PetalWidth=cluster; plot SepalWidth*PetalLength=cluster; plot SepalWidth*PetalWidth=cluster; plot PetalLength*PetalWidth=cluster; run; title "Scatterplot Matrix for Iris Data"; proc sgscatter data=sashelp.iris; : [おまけ1] 散布図行列 matrix SepalLength SepalWidth PetalLength PetalWidth / group=Species; run; title "Scatterplot Matrix with histogram for Iris Data"; proc sgscatter data=sashelp.iris; : [おまけ2] ヒストグラム付き散布図行列 matrix SepalLength SepalWidth PetalLength PetalWidth / group=Species diagonal=(kernel histogram); run; title;

出力

[結果タグ] html 形式: les1405-Results.html
[結果タグ] pdf 形式: les1405-Output.pdf
テキストファイル: les1405-Results.txt

Sashelp.iris --- Fisher's Iris Data (1936) 73 Thursday, January 20, 2022 05:27:13 PM The CONTENTS Procedure Variables in Creation Order # Variable Type Len Label 1 Species Char 10 Iris Species 2 SepalLength Num 8 Sepal Length (mm) 3 SepalWidth Num 8 Sepal Width (mm) 4 PetalLength Num 8 Petal Length (mm) 5 PetalWidth Num 8 Petal Width (mm) The First Five Observations Out of 150 74 Thursday, January 20, 2022 05:27:13 PM Sepal Sepal Petal Petal Species Length Width Length Width Setosa 50 33 14 2 Setosa 46 34 14 3 Setosa 46 36 10 2 Setosa 51 33 17 5 Setosa 55 35 13 2 The Species Variable Thursday, January 20, 2022 05:27:13 PM 75 The FREQ Procedure Iris Species Cumulative Cumulative Species Frequency Percent Frequency Percent --------------------------------------------------------------- Setosa 50 33.33 50 33.33 Versicolor 50 33.33 100 66.67 Virginica 50 33.33 150 100.00 The Species Variable Thursday, January 20, 2022 05:27:13 PM 76 The FASTCLUS Procedure Replace=FULL Radius=0 Maxclusters=3 Maxiter=1 Initial Seeds Cluster SepalLength SepalWidth PetalLength PetalWidth --------------------------------------------------------------------------- 1 77.00000000 38.00000000 67.00000000 22.00000000 2 57.00000000 44.00000000 15.00000000 4.00000000 3 49.00000000 25.00000000 45.00000000 17.00000000 Criterion Based on Final Seeds = 3.7097 Cluster Summary Maximum Distance RMS Std from Seed Radius Nearest Cluster Frequency Deviation to Observation Exceeded Cluster ----------------------------------------------------------------------------- 1 33 3.8831 12.9226 3 2 50 2.7803 12.4803 3 3 67 4.1797 18.5320 1 The Species Variable Thursday, January 20, 2022 05:27:13 PM 77 The FASTCLUS Procedure Replace=FULL Radius=0 Maxclusters=3 Maxiter=1 Cluster Summary Distance Between Cluster Cluster Centroids ----------------------------- 1 18.3409 2 34.2516 3 18.3409 Statistics for Variables Variable Total STD Within STD R-Square RSQ/(1-RSQ) --------------------------------------------------------------------- SepalLength 8.28066 4.48242 0.710915 2.459187 SepalWidth 4.35866 3.24819 0.452092 0.825123 PetalLength 17.65298 4.29764 0.941527 16.101961 PetalWidth 7.62238 2.38707 0.903243 9.335201 OVER-ALL 10.69224 3.70171 0.881751 7.456709 Pseudo F Statistic = 548.07 The Species Variable Thursday, January 20, 2022 05:27:13 PM 78 The FASTCLUS Procedure Replace=FULL Radius=0 Maxclusters=3 Maxiter=1 Approximate Expected Over-All R-Squared = 0.62728 Cubic Clustering Criterion = 24.559 WARNING: The two values above are invalid for correlated variables. Cluster Means Cluster SepalLength SepalWidth PetalLength PetalWidth --------------------------------------------------------------------------- 1 69.00000000 30.96969697 58.27272727 21.27272727 2 50.06000000 34.28000000 14.62000000 2.46000000 3 59.47761194 27.61194030 44.52238806 14.53731343 Cluster Standard Deviations Cluster SepalLength SepalWidth PetalLength PetalWidth --------------------------------------------------------------------------- 1 5.012484414 2.909948974 4.577613511 2.401467354 2 3.524896872 3.790643691 1.736639965 1.053855894 3 4.831582365 2.953966126 5.360795421 3.011736428 The Species Variable Thursday, January 20, 2022 05:27:13 PM 79 Plot of SepalLength*SepalWidth. Symbol is value of CLUSTER. | S | e 84 + p | a | 1 1 1 1 l | 1 1 72 + 1 1 1 1 L | 3 1 3 3 e | 1 3 3 3 3 3 1 n | 3 3 3 3 1 3 1 3 3 1 g 60 + 3 3 3 3 3 3 3 3 t | 3 3 3 3 3 3 2 2 2 h | 3 3 3 3 3 2 2 2 2 2 |3 3 3 3 2 2 2 2 2 2 2 2 2 ( 48 + 3 3 2 2 2 2 2 m | 2 2 2 2 2 2 2 m | 2 ) | 36 + | -+-------------+-------------+-------------+-------------+-------------+- 20 25 30 35 40 45 Sepal Width (mm) NOTE: 64 obs hidden. The Species Variable Thursday, January 20, 2022 05:27:13 PM 80 Plot of SepalLength*PetalLength. Symbol is value of CLUSTER. | S | e 84 + p | a | 1 1 1 1 l | 1 1 72 + 1111 1 L | 333 1 11 1 1 e | 3 33 331 1111 n | 333 33333 1111 1 g 60 + 3 3 3333333 3 t | 2 2 2 33 3333 3 333 h | 222 2 33 3 33 | 222222 2 3 3 3 3 ( 48 + 2222 2 3 3 m | 2 222 m | 2 ) | 36 + | ---+---------+---------+---------+---------+---------+---------+-- 10 20 30 40 50 60 70 Petal Length (mm) NOTE: 53 obs hidden. The Species Variable Thursday, January 20, 2022 05:27:13 PM 81 Plot of SepalLength*PetalWidth. Symbol is value of CLUSTER. | S | e 84 + p | a | 1 1 1 l | 1 1 72 + 1 1 1 1 L | 3 3 1 1 e | 3 3 3 3 1 1 1 1 1 1 1 n | 3 3 3 3 3 1 1 1 1 1 g 60 + 3 3 3 3 3 3 3 t | 2 2 2 3 3 3 3 3 3 3 3 h | 2 2 3 3 3 3 3 | 2 2 2 2 2 2 3 3 3 ( 48 + 2 2 2 3 3 m | 2 2 m | 2 ) | 36 + | ---+---------+---------+---------+---------+---------+-- 0 5 10 15 20 25 Petal Width (mm) NOTE: 74 obs hidden. The Species Variable Thursday, January 20, 2022 05:27:13 PM 82 Plot of SepalWidth*PetalLength. Symbol is value of CLUSTER. | 50 + S | e | p | a | 2 l | 22 40 + 2 W | 2 222 2 1 1 i | 2 22 1 d | 2222 t | 2222 2 3 3 1 11 1 h | 22222 33 333 1 11111 11 30 + 2 22 2 33 333 33331 1 11 1 1 ( | 2 3 3333 33333 3 1 1 1 1 m | 3 3333 3 3 3 3 3 1 m | 3 33 3 33 1 ) | 2 3 33 3 3 | 3 3 3 20 + 3 ---+---------+---------+---------+---------+---------+---------+-- 10 20 30 40 50 60 70 Petal Length (mm) NOTE: 39 obs hidden. The Species Variable Thursday, January 20, 2022 05:27:13 PM 83 Plot of SepalWidth*PetalWidth. Symbol is value of CLUSTER. | 50 + S | e | p | a | 2 l | 2 2 40 + 2 W | 2 2 2 1 1 i | 2 2 2 1 d | 2 2 2 t | 2 2 2 2 3 1 1 1 1 h | 2 2 3 3 3 3 1 1 1 30 + 2 2 2 3 3 3 3 1 3 1 1 1 1 1 ( | 2 3 3 3 3 3 1 1 1 1 3 m | 3 3 3 3 3 3 3 1 m | 3 3 3 3 1 3 3 ) | 2 3 3 3 | 3 3 20 + 3 ---+---------+---------+---------+---------+---------+-- 0 5 10 15 20 25 Petal Width (mm) NOTE: 69 obs hidden. The Species Variable Thursday, January 20, 2022 05:27:13 PM 84 Plot of PetalLength*PetalWidth. Symbol is value of CLUSTER. | P 72 + e | 1 1 1 1 t | 1 1 a 60 + 1 1 1 1 1 1 1 l | 3 1 1 1 1 1 1 | 3 3 3 3 3 1 1 3 L 48 + 3 3 3 3 3 3 3 e | 3 3 3 3 3 3 n | 3 3 3 3 3 g 36 + 3 3 t | 3 3 h | 24 + ( | 2 2 m | 2 2 2 2 2 2 m 12 + 2 2 2 2 ) | | 0 + ---+---------+---------+---------+---------+---------+-- 0 5 10 15 20 25 Petal Width (mm) NOTE: 88 obs hidden.

お役立ちサイト: マニュアル
- Base SAS 9.2 Procedures Guide
  - 左側のProcedures: 各種手法(主に基礎的手法)の使い方やExampleが掲載されている。
    - [例] means, univariate, freq, plot, sort
- SAS/STAT 9.2 User's Guide, Second Edition
  - 左側のProcedures: 各種手法(上記以外)の使い方やExampleが掲載されている。
    - [例] reg, princomp, factor, fastclus, sgscatter
- Base SAS
- SAS/STAT
- Programming: SAS 9.4 and Viya
- Popular Documentation
  
  SAS_Graphに関して
- SAS University Editionでグラフを描く
- もしも，SASのsgplotとRのggplot2 を比較したら…
課題提出(レポート): 詳細は第13回第10節を参照下さい。
提出期限は「02月10日(木) 夕方」です。厳守。
ファイルを提出できるように設定しておきましたが、私の設定ミス等でうまく提出できない場合は電子メールでその旨、ご連絡ください。
まとめ: 講義を終えるにあたって
　この講義を通して、「データとの接し方」や「統計の考え方」が多少なりとも理解できたであろうか? 大量の数値群からその中に内在する構造を見つけることが「解析」であり「統計の面白味」でもあると思う。近年はビッグデータへの注目もあって統計学が脚光を浴びているが、どのような時代であっても、「データに内在する構造」を明らかにする学問として今後もいろいろな分野で応用されていくと考えている。そのためには、統計理論や解析目的を知っている必要があるのは勿論だが、対象とするデータの背景を知っておくことや、統計ソフトを"道具"として使いこなす技術も習得する必要がある。
　そうした中でSASは長年の歴史を有した信頼できる統計ソフトであり、研究での利用に関してフルスペックの機能が無料で使えるという魅力もあり、知っておいて損のない統計環境だと感じているので、今回講義で利用することにした次第である。
　今後、各自の実験や修論をまとめる際だけでなく、日頃目にする新聞や雑誌と言った生活等のいろいろな場面で、種々のデータに出会うことになると思うが、提示された数値にはどの様な意味(と意図)があり、どう理解して、個々人としてどうアクションを起すかの、一つの判断手段として活用してもらえれば幸いである。
　今後、もし統計に関して何か疑問に出会い、私に連絡・相談してみたいと思った時は、遠慮無くご連絡ください。
　皆さんのご期待にどこまで応えられか心許無い部分もありますが、 Q3、Q4を通しての講義、お疲れ様でした。 COVID-19の動向は依然として不透明ですが、ご健康には留意されてお過ごしください。お元気で。
連絡したければ
- 電子メールアドレス : hayashi.atsuhiro@nitech.ac.jp
- 研究室: 11号館 5階 516室
- 電話: 内線 5119, 052-735-5119
次回は、... : 02月01日 (最終回)
- 各自のレポート作成の時間に充ててください。出席(閲覧)は不問です。
- ご質問があれば研究室を訪問ください。待機しています。16:20-17:50。
- 【おまけ】第13回で飛ばしたTips 2, Tips 3を紹介します。【リモート講義】
  - 今週を最後とするつもりでしたが、詰め込み過ぎになったので、一部を来週に回します。
  - 出席等は求めません。興味があれば視聴ください。
  - 今週の内容について質問があれば、回答することにします。