回帰分析(後編)、主成分分析

統計モデル解析特論I/II : 第05回 (11/09/17)

　今回は、前回に説明できなかった部分を終えた後、多変量解析の1手法である回帰分析について解説する。
　本講義では、数理的な導出方法よりも、当該手法がどのようなモデルやアイディアに基づいて構成されたかの考え方や原理に重点を置いて説明する。その最初に紹介する統計解析手法として、回帰分析を取り上げる。工学系の実験等を行う領域では頻繁に使用される手法であるが、日常的な話題の中でも概念は広く利用されているので、取っ付き易い手法ではないだろうか。
　過去のデータからその構造を把握し、新規に測定されたデータに対する予測を行ないたいと言うときに、回帰分析は有用である。構造のシンプルな単回帰分析でこの手法の原理を理解し、複数の説明変量を用いた重回帰分析に拡張する。残差の取り方や、その二乗和を最少にするという考えは同じである。
　ここでは SAS プログラムとその出力を例示しながら説明するが、手慣れたプログラム言語を使って各自実践してほしい。
　また、残りの時間で主成分分析を紹介する。

　● 目次: 回帰分析
　　 1.　単回帰分析 : 予測等に使う、連続変量の関係 [第03回の資料へジャンプ]
　　 2.　「体重の大きい者を除外」して実行するには? [第03回の資料へジャンプ]
　　 3.　有効桁数に注意せよ : どこまでが「意味ある桁」か? [第03回の資料へジャンプ]
　　 4.　重回帰分析 : 2変量以上の説明する変量(説明変量)で 1変量(目的変量)を説明
　　 5.　特定グループでの解析
　　 6.　[要点] 解析する上での注意点
　　 7.　誤用?!
　　 8.　4つの尺度と回帰分析
　　 9.　回帰分析における変数選択、総当たり法

　● 目次: 主成分分析
　　 10.　2変量の場合の主成分分析 : 理解を助けるため
　　 11.　3変量以上の主成分分析
　　 12.　相関係数を使う理由
　　 13.　主成分の数の決定基準

　● 回帰分析 : 連続変量の予測

重回帰分析 : 2変量以上の説明する変量(説明変量)で 1変量(目的変量)を説明
- 説明変量が複数になる : 単 ===> 重
- 体重を身長と胸囲で説明したい。予測したい。
- [体重]=a[身長]+b[胸囲]+c : 回帰係数を求めたい。
- 単回帰とアイディアは同じ : 残差(予測誤差)の二乗和を最小にする(最小二乗法)
- 説明される変量(目的変量)と平行に残差を取る。
1. プログラム : les0501.sas
  /* Lesson 05-1 */ /* File Name = les0501.sas 11/09/17 */ options linesize=72 pagesize=20; options nocenter linesize=78 pagesize=30; proc printto log = '/folders/myfolders/Kougi/SAS_log0501.txt' print = '/folders/myfolders/Kougi/SAS_out0501.txt' new; ods listing gpath='/folders/myfolders/Kougi/SAS_ODS99'; data gakusei; infile '/folders/myfolders/Kougi/all07ae.prn' firstobs=2; input sex $ shintyou taijyuu kyoui jitaku $ kodukai carryer $ tsuuwa; if sex^='M' & sex^='F' then delete; if shintyou=. | taijyuu=. 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; :
2. 出力結果 : SAS_out0501.txt, les0501_out.pdf
  2017年11月 9日木曜日 08時39分18秒 2 REG プロシジャモデル : MODEL1 従属変数 : taijyuu 読み込んだオブザベーション数 325 使用されたオブザベーション数 114 欠損値を含むオブザベーション数 211 分散分析要因自由度平方和平均平方 F 値 Pr > F Model 2 8070.70705 4035.35353 85.10 <.0001 Error 111 5263.40733 47.41808 Corrected Total 113 13334 Root MSE 6.88608 R2 乗 0.6053 従属変数の平均 58.79298 調整済み R2 乗 0.5982 変動係数 11.71242 2017年11月 9日木曜日 08時39分18秒 3 REG プロシジャモデル : MODEL1 従属変数 : taijyuu パラメータの推定パラメータ変数自由度推定値標準誤差 t 値 Pr > |t| Intercept 1 -106.30023 12.75197 -8.34 <.0001 shintyou 1 0.80655 0.07854 10.27 <.0001 kyoui 1 0.34947 0.08192 4.27 <.0001 2017年11月 9日木曜日 08時39分24秒 4 OBS sex shintyou taijyuu kyoui jitaku kodukai carryer tsuuwa pred1 resid1 1 F 145.0 38.0 . J 10000 . . . 2 F 146.7 41.0 85 J 10000 Vodafone 6000 41.7256 -0.72559 3 F 148.0 42.0 . J 50000 . . . 4 F 148.0 43.0 80 J 50000 DoCoMo 4000 41.0267 1.97328 5 F 149.0 45.0 . G 60000 . . . 6 F 150.0 46.0 86 40000 . 44.7367 1.26333 7 F 151.0 45.0 . J 20000 docomo 5000 . . 8 F 151.0 50.0 . G 60000 J-PHONE . . . 9 F 151.7 41.5 80 J 35000 . 44.0109 -2.51095 10 F 152.0 35.0 77 J 60000 DoCoMo 2000 43.2045 -8.20449 2017年11月 9日木曜日 08時39分24秒 5 プロット : taijyuu*shintyou 凡例 : A = 1 obs, B = 2 obs, ... 100 + A | | A | A | A | A A 75 + A A | A A BAA A A A A | BB A A BAAA A A A taijyuu | A A A CA C AA B | B A AAAA C AA AA A A A A | A A AA A BABA AAA A 50 + AC B CD AAB B | A A A B A A BA A | AA A A | | A | 25 + --+------------+------------+------------+------------+------------+-- 140 150 160 170 180 190 shintyou 2017年11月 9日木曜日 08時39分24秒 6 プロット : taijyuu*kyoui 凡例 : A = 1 obs, B = 2 obs, ... 100 + A | | A | A | A | A A 75 + A A | BD AB A | A A CAAAAB B A taijyuu | A A AGA B A | BB ADBDB A | A A A AAA CC AB 50 + A A DD G B | B DAA BA | B B | | A | 25 + ---+-------------+-------------+-------------+-------------+-- 40 60 80 100 120 kyoui 2017年11月 9日木曜日 08時39分24秒 7 プロット : taijyuu*pred1 凡例 : A = 1 obs, B = 2 obs, ... 100 + A | | A | A | A | A A 75 + A A | A AB A B AAA | A B AA ABAA A A A A taijyuu | A A A BABB AA B | AAA A A B BAAAA B B A | A A A AA A CABA A A 50 + A B B CBAABAAAB | AA C AA ABA | AA A A | | A | 25 + ---+-------------+-------------+-------------+-------------+-- 40 50 60 70 80 予測値 taijyuu 2017年11月 9日木曜日 08時39分24秒 8 プロット : resid1*pred1 凡例 : A = 1 obs, B = 2 obs, ... 40 + | | A | | A | A 残 20 + 差 | | A A | A A A A A | A A A A A BA | A A A A A A A B A BA B A 0 +----A---AA---B-B-AAA---A-C--A-A-ABAAB-AAA----A--------------- | A B AABBAABAAAAA A A AAAB AA A AAB | A A AB B B A A A B B A | A A A | A | -20 + ---+-------------+-------------+-------------+-------------+-- 40 50 60 70 80 予測値 taijyuu 2017年11月 9日木曜日 08時39分24秒 9 プロット : resid1*shintyou 凡例 : A = 1 obs, B = 2 obs, ... 40 + | | A | | A | A 残 20 + 差 | | A A | A A A AA | A A B A A A A | A BA A A BA A BA AA A 0 +----------A---A---A-AAA-A--BA--A--AAAA-AEA-A--AA-A-A----------------- | A B B DC AAB BA A B BAA A A B A A | A A AA CB A A A ABA A | A A A | A | -20 + --+------------+------------+------------+------------+------------+-- 140 150 160 170 180 190 shintyou 2017年11月 9日木曜日 08時39分24秒 10 プロット : resid1*kyoui 凡例 : A = 1 obs, B = 2 obs, ... 40 + | | A | | A | A 残 20 + 差 | | A A | A AA B | A AA A BB | A BAAAB AAA B B A 0 +-----------------------A---C--A-ABFBBE--AB-----A------------- | EBCAGAEDA AB | AAABBA D AC A | BA | A | -20 + ---+-------------+-------------+-------------+-------------+-- 40 60 80 100 120 kyoui 2017年11月 9日木曜日 08時39分24秒 11 プロット : resid1*taijyuu 凡例 : A = 1 obs, B = 2 obs, ... 40 + | | A | | A | A 残 20 + 差 | | A A | A A A AA | A AA A BB | A AA A A B CAABA A 0 +-----------------A-AAACA-B-C-ACBD-C---A---------------------- | A CADCD AB BCBA B AB | A A ABB B A AC B A | A A A | A | -20 + ---+-------------+-------------+-------------+-------------+-- 20 40 60 80 100 taijyuu
3. 結果の見方
  - 対象になったのは 114名。
  - 説明変量群が予測に役立っているか?
    - 回帰に役立っているか : Prob>F : 小さいと有意
    - 「役立っている」と言える : 0.01% だから 1% で有意
  - 決定係数 : R-Square ( 相関係数 : R )
    - 目的変量が説明変量でどの程度説明しているかの割合。
    - 1 に近いほど当てはまりが良いと言える。: 60.5%
    - 説明変量数が増えると大きくなるのが一般的。
  - 回帰係数 : Parameter Estimate
    - 回帰式: a=0.807, b=0.349, c=-106
  - ある特定の説明変量が予測に役立っているか?
    - 回帰係数の検定(帰無仮説:係数=0 か?) : Prob>|T| : 小さいと有意
    - 両方とも(身長も胸囲も)有意
    - 「各係数は 0ではない」と言える : 0.01% だから 1% で有意
  - 残差の性質 ===> 正規性 : 残差プロット、残差解析
    - 残差(予測誤差)は正規分布をしていると仮定してモデルが構築されている。
    - この仮定が覆ると、回帰分析として成立していないことになる。
    - 残差が正規分布をしているか確認する必要がある。
    - 均等に散らばっているか?
    - 傾向はないか? : もし傾向があると言うことになれば正規性の仮定が崩れている
    - 体重の大きい 3例程度が外れ値と考えられるか要確認 ===> [演習1](第5節)
    - ...
  - ...

特定グループでの解析

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

プログラム : les0502.sas

/* Lesson 12-2 */ /* File Name = les1202.sas 07/05/16 */ data gakusei; infile 'all07ae.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;

出力結果 : SAS_out0502.txt, les0502_out.pdf
CORR プロシジャ単純統計量変数 N 平均標準偏差合計 taijyuu 242 62.23884 7.92774 15062 kyoui 71 88.09859 9.68527 6255 kodukai 229 48620 52677 11134000 tsuuwa 88 6422 4521 565098 単純統計量変数最小値最大値 taijyuu 46.00000 100.00000 kyoui 46.00000 112.00000 kodukai 0 350000 tsuuwa 0 30000 2017年11月 9日木曜日 08時58分30秒 20 CORR プロシジャ Pearson の相関係数 H0: Rho=0 に対する Prob > |r| オブザベーション数 shintyou taijyuu kyoui kodukai tsuuwa shintyou 1.00000 0.43758 0.15872 0.07647 -0.03430 <.0001 0.1862 0.2491 0.7510 242 242 71 229 88 taijyuu 0.43758 1.00000 0.40227 0.04119 -0.01583 <.0001 0.0005 0.5352 0.8836 242 242 71 229 88 kyoui 0.15872 0.40227 1.00000 -0.37945 -0.38661 0.1862 0.0005 0.0015 0.1721 71 71 71 67 14 kodukai 0.07647 0.04119 -0.37945 1.00000 0.24685 0.2491 0.5352 0.0015 0.0219 229 229 67 229 86 tsuuwa -0.03430 -0.01583 -0.38661 0.24685 1.00000 0.7510 0.8836 0.1721 0.0219 88 88 14 86 88 2017年11月 9日木曜日 08時58分30秒 23 REG プロシジャモデル : MODEL1 従属変数 : taijyuu 読み込んだオブザベーション数 242 使用されたオブザベーション数 71 欠損値を含むオブザベーション数 171 2017年11月 9日木曜日 08時58分30秒 24 REG プロシジャモデル : MODEL1 従属変数 : taijyuu 分散分析要因自由度平方和平均平方 F 値 Pr > F Model 2 1596.38065 798.19033 13.06 <.0001 Error 68 4155.98301 61.11740 Corrected Total 70 5752.36366 Root MSE 7.81776 R2 乗 0.2775 従属変数の平均 64.72817 調整済み R2 乗 0.2563 変動係数 12.07784 2017年11月 9日木曜日 08時58分30秒 25 REG プロシジャモデル : MODEL1 従属変数 : taijyuu パラメータの推定パラメータ変数自由度推定値標準誤差 t 値 Pr > |t| Intercept 1 -54.72134 27.50850 -1.99 0.0507 shintyou 1 0.52620 0.15946 3.30 0.0015 kyoui 1 0.32534 0.09772 3.33 0.0014 2017年11月 9日木曜日 08時58分32秒 27 プロット : taijyuu*shintyou 凡例 : A = 1 obs, B = 2 obs, ... taijyuu | 100 + A | A A A | A A 75 + A B B A BAA C AAA B A A | B B BC B HAC DF D G D DC C K A B B A | A C BB C CAFGBD M JBKK FAE DBCC A BA AA 50 + A B A C A AACCAA A C C B C | | 25 + ---+--------+--------+--------+--------+--------+--------+--------+-- 155 160 165 170 175 180 185 190 shintyou 2017年11月 9日木曜日 08時58分32秒 28 プロット : taijyuu*kyoui 凡例 : A = 1 obs, B = 2 obs, ... (NOTE: 171 obs が欠損値です。) taijyuu | 100 + A | A A | A 75 + AA C A A | A A CABIBBAD A | A A BCA ADBEF AA A 50 + A A AA | | 25 + ---+-------------+-------------+-------------+-------------+-- 40 60 80 100 120 kyoui 2017年11月 9日木曜日 08時58分32秒 29 プロット : taijyuu*pred1 凡例 : A = 1 obs, B = 2 obs, ... (NOTE: 171 obs が欠損値です。) taijyuu | 100 + A | A A | A 75 + AA AA AA A | A A BACAB BAC AAABAAB A | AA A A A A B B CACDAAAABA A 50 + A A AA | | 25 + --+----------+----------+----------+----------+----------+----------+- 50 55 60 65 70 75 80 予測値 taijyuu 2017年11月 9日木曜日 08時58分32秒 30 プロット : resid1*pred1 凡例 : A = 1 obs, B = 2 obs, ... (NOTE: 171 obs が欠損値です。) | 50 + 残 | 差 | A 25 + A A | A | A AABAA A A 0 +-----A---------A--A-A-A---B-B-CACDB-BAC-AAB-ABB--------A------------- | A AA A A A A BA AAABA BA A | -25 + --+----------+----------+----------+----------+----------+----------+- 50 55 60 65 70 75 80 予測値 taijyuu 2017年11月 9日木曜日 08時58分32秒 31 プロット : resid1*shintyou 凡例 : A = 1 obs, B = 2 obs, ... (NOTE: 171 obs が欠損値です。) | 50 + 残 | 差 | A 25 + A A | A | A A A A B B A 0 +---------------AA-----A-CD-A-G-A-BB-BAC-A-AB-A-B---A--A------------- | A B A A A A A AA A A A A BAA A A | -25 + ---+--------+--------+--------+--------+--------+--------+--------+-- 155 160 165 170 175 180 185 190 shintyou 2017年11月 9日木曜日 08時58分32秒 32 プロット : resid1*kyoui 凡例 : A = 1 obs, B = 2 obs, ... (NOTE: 171 obs が欠損値です。) | 50 + 残 | 差 | A 25 + A A | A | A BD B 0 +-------------A---------A---A-AB-AADACHBABE--B--A---A--------- | A ABB CBCD A A | -25 + ---+-------------+-------------+-------------+-------------+-- 40 60 80 100 120 kyoui 2017年11月 9日木曜日 08時58分32秒 33 プロット : resid1*taijyuu 凡例 : A = 1 obs, B = 2 obs, ... (NOTE: 171 obs が欠損値です。) | 50 + 残 | 差 | A 25 + A A | A | A AAB A A AA 0 +----------------A-AA---FADE--GBAB-DB-------A------------------------- | A A A CAABAE A B A | -25 + --+----------+----------+----------+----------+----------+----------+- 40 50 60 70 80 90 100 taijyuu

結果の見方
- 単変量毎の相関が有意なのは、身長と体重、体重と胸囲の間。
- 対象になったのは 71名。
- 回帰に役立っているか : 役立っている : 0.01% だから 1% で有意
- 決定係数(R-square)は 27.8%
- 個々の説明変量が予測に役立っているか?
  - 係数がゼロか? : 定数項も身長も胸囲も有意(1% で有意)
- 残差の性質 ===> 正規性 : 残差プロット、残差解析
  - 均等に散らばっているか?
  - 傾向はないか? : 傾向があると言うことは正規性の仮定が崩れていること
  - 外れ値? 85Kg より重い 3名程度が吟味対象?
[演習1] : 「男性のみ」で、かつ「体重の大きい 3名を除外」して実行してみよ。
- プログラム : les0503.sas、出力結果 : SAS_out0503.txt, les0503_out.pdf
  where sex='M' and taijyuu<85;
- 当てはまりは良くなったか? : 異常値と外れ値の意味するもの
- 残差の正規性はどのように変化したか?

[要点] 解析する上での注意点
- 残差解析を行なう
- 残差に傾向はないか? : バナナカーブ、ラッパ状、...
- 対象集団を単質のものに絞る : 男性だけ、18歳だけ、...
- 外れ値を持つ個体を除外する : 標準体型の者だけ、...
  - ただし、恣意的な除外は謹むべきである
- 「外れ値」と「異常値」は区別して取り扱う・識別する : 「外れて」いても「異常」かどうかは別
  - 外れ値 : 大多数のデータから大きく離れた値
  - 異常値 : 通常ではあり得ない入力ミス等の異常な値
- 上記以外に「多重共線性(multi-colinearity)」の問題もある。逆行列が取れなくなる。不安定な解。
誤用?!
- 構造がシンプルで理解し易い分、"変な"利用も散見される
- 内挿(観測点の内側)はまだしも、外挿(観測点の外側)を予測するのは難しい or 無謀。
- 予測範囲でデータの構造が一定という"条件"が必要。
- 線形で表現できる関係なのか? 線形の関係式が有効なのか? 非線型?
- 自分で判断できる能力を持つこと。疑うこと。
　[例1] 人間の成長曲線
　[例2] 将来のプログラマ必要数予測 : 21世紀(?)には国民全員がプログラマ ('80s)
　[例3] オリンピック 100m 走の男女記録 : 2156年には女性の方が速い (2004.09.30) :
　　　　　　　 Japan Journal LTD の記事 , 朝日新聞の記事
　　　　　[究極の命題!] 100m に 0.00秒要する(!?)ようになるのは何時?
4つの尺度と回帰分析
- 名義尺度(nominal scale) : 性別、職業、支持政党。分類の項目。
- 順序尺度(ordinal scale) : 成績の優・良・可・不可。順序関係あり。「優と良の差」が「可と不可の差」に等しいわけではない。
- 間隔尺度(interval scale) : セ氏やカ氏で測定された温度。原点は任意。
- 比率尺度(ratio scale) : 長さ、重さ。間隔尺度に加えて絶対原点が存在。
- 質的データ : 名義尺度、順序尺度。離散変量とも。
- 量的データ : 間隔尺度、比率尺度。連続変量とも。
- 離散変量を回帰分析で用いるには注意が必要である。線形関係と言えるのか? 割り当てた値を変更したら影響力が変化する。

回帰分析における変数選択 : 　回帰分析では回帰係数や重相関係数を知ることだけでなく残差解析も重要であることを強調したつもりである。次に説明変数の取捨選択(変数選択)について説明する。

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

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

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

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

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

総当たり法(rsquare)

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

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

プログラム : les0504.sas

/* Lesson 05-4 */ /* File Name = les0504.sas 11/09/17 */ options linesize=72 pagesize=20; options nocenter linesize=78 pagesize=20; proc printto log = '/folders/myfolders/Kougi/SAS_log0504.txt' print = '/folders/myfolders/Kougi/SAS_out0504.txt' new; ods listing gpath='/folders/myfolders/Kougi/SAS_ODS99'; data air; infile '/folders/myfolders/Kougi/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; :

出力結果 : SAS_out0504.txt, les0504_out.pdf
2017年11月 9日木曜日 10時19分59秒 1 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 2017年11月 9日木曜日 10時19分59秒 2 CORR プロシジャ 7 変数 : y x1 x2 x3 x4 x5 x6 単純統計量変数 N 平均標準偏差合計 y 41 30.04878 23.47227 1232 変数 N 平均標準偏差合計 x1 41 55.76341 7.22772 2286 x2 41 463.09756 563.47395 18987 x3 41 608.60976 579.11302 24953 x4 41 9.44390 1.42864 387.20000 x5 41 36.76902 11.77155 1508 x6 41 113.90244 26.50642 4670 2017年11月 9日木曜日 10時20分00秒 5 CORR プロシジャ Pearson の相関係数, N = 41 H0: Rho=0 に対する Prob > |r| 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 2017年11月 9日木曜日 10時20分00秒 8 REG プロシジャモデル : MODEL1 従属変数 : y 読み込んだオブザベーション数 41 使用されたオブザベーション数 41 分散分析要因自由度平方和平均平方 F 値 Pr > F Model 6 14755 2459.10601 11.48 <.0001 Error 34 7283.26641 214.21372 Corrected Total 40 22038 Root MSE 14.63604 R2 乗 0.6695 従属変数の平均 30.04878 調整済み R2 乗 0.6112 変動係数 48.70761 パラメータの推定パラメータ変数自由度推定値標準誤差 t 値 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 2017年11月 9日木曜日 10時20分04秒 11 プロット : resid1*pred1 凡例 : A = 1 obs, B = 2 obs, ... | 50 + A 残 | 差 | A 25 + | A A B | B B A A B 0 +--------------AAA-----ABAAAA---------A-----------------------A------ | A CB BA | C A A -25 + A ---+--------+--------+--------+--------+--------+--------+--------+-- -20 0 20 40 60 80 100 120 予測値 y 2017年11月 9日木曜日 10時20分04秒 19 REG プロシジャモデル : MODEL1 従属変数 : y 読み込んだオブザベーション数 41 使用されたオブザベーション数 41 ステップワイズ法: ステップ 1 変数 x2 の追加 : R2 乗 = 0.4157 C(p) = 23.1089 分散分析要因自由度平方和平均平方 F 値 Pr > F Model 1 9161.74469 9161.74469 27.75 <.0001 Error 39 12876 330.15789 Corrected Total 40 22038 パラメータ Type II 変数推定値標準誤差平方和 F 値 Pr > F Intercept 17.61057 3.69159 7513.50474 22.76 <.0001 x2 0.02686 0.00510 9161.74469 27.75 <.0001 条件数における境界 : 1, 1 ------------------------------------------------------------------------------ ステップワイズ法: ステップ 2 2017年11月 9日木曜日 10時20分04秒 22 REG プロシジャモデル : MODEL1 従属変数 : y 変数 x3 の追加 : R2 乗 = 0.5863 C(p) = 7.5586 分散分析要因自由度平方和平均平方 F 値 Pr > F Model 2 12921 6460.63359 26.93 <.0001 Error 38 9116.63526 239.91145 Corrected Total 40 22038 パラメータ Type II 変数推定値標準誤差平方和 F 値 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 条件数における境界 : 11.434, 45.735 ------------------------------------------------------------------------------ ステップワイズ法: ステップ 3 REG プロシジャモデル : MODEL1 従属変数 : y 変数 x6 の追加 : R2 乗 = 0.6174 C(p) = 6.3610 分散分析要因自由度平方和平均平方 F 値 Pr > F Model 3 13606 4535.41173 19.90 <.0001 Error 37 8431.66725 227.88290 Corrected Total 40 22038 パラメータ Type II 変数推定値標準誤差平方和 F 値 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 条件数における境界 : 12.65, 78.633 ------------------------------------------------------------------------------ モデル内のすべての変数は水準 0.1500 で有意です。モデルへの変数追加で、他の変数は有意水準 0.1500 で満たされていません。ステップワイズ法の要約変数の変数の取り込んだステップ追加削除変数の数偏 R2 乗モデル R2 乗 1 x2 1 0.4157 0.4157 2 x3 2 0.1706 0.5863 3 x6 3 0.0311 0.6174 ステップワイズ法の要約ステップ C(p) F 値 Pr > F 1 23.1089 27.75 <.0001 2 7.5586 15.67 0.0003 3 6.3610 3.01 0.0913 2017年11月 9日木曜日 10時20分06秒 31 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 2017年11月 9日木曜日 10時20分06秒 32 プロット : resid1*pred1 凡例 : A = 1 obs, B = 2 obs, ... 50 + A | 残 | A 差 | B | A AAB A A A 0 +-----A-A--A--BAA-AAA---------A-------------------------------A------- | BAA BAAA A A | AA A A A | A | -50 + --+----------+----------+----------+----------+----------+----------+- 0 20 40 60 80 100 120 予測値 y 2017年11月 9日木曜日 10時20分06秒 46 REG プロシジャモデル : MODEL1 従属変数 : y R2 乗選択法読み込んだオブザベーション数 41 使用されたオブザベーション数 41 取り込んだ変数の数 R2 乗モデルの独立変数 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

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

● 主成分分析
　いくつか(p個)の変量の値を情報の損失をできるだけ少なくして、少数変量(m個、m＜p)の総合的指標(主成分)で代表させる方法として主成分分析(Principal Component Analysis, PCA)と因子分析(Factor Analysis, FA)がある。いくつかのテストの成績を総合した総合的成績、いろいろな症状を総合した総合的な重症度、種々の財務指標に基づく企業の評価等を求めたいといった場合に用いられる。 p変量(p次元)の観測値をm個(m次元)の主成分に縮約させるという意味で、次元を減少させる(reduce)方法と言うこともでき、多変量データを要約する有力な方法である。
　両手法は似た目的に使われるが、元になっている考え方は異なるので利用する場面では注意が必要である。違いに焦点を当てながらまずは主成分分析を説明する。

2変量の場合の主成分分析 : 理解を助けるため

定式化 : 配布資料 54ページ～

重み(係数) : a1、a2
合成変量(線形結合) : z
よく代表するように、a1 と a2 を決める。
より広がって測定できる軸に沿うと情報量が多い。
　　　[参考:立体の測定] ノギス、ノギスの使い方 _
全測定値の分散を最大化する軸を決定する。

身長と体重の総合指標を求めよう : プログラム : les0511.sas

/* Lesson 05-11 */ /* File Name = les0511.sas 11/09/17 */ options linesize=72 pagesize=20; options nocenter linesize=78 pagesize=20; proc printto log = '/folders/myfolders/Kougi/SAS_log0511.txt' print = '/folders/myfolders/Kougi/SAS_out0511.txt' new; ods listing gpath='/folders/myfolders/Kougi/SAS_ODS99'; data gakusei; infile '/folders/myfolders/Kougi/all07ae.prn' firstobs=2; input sex $ shintyou taijyuu kyoui jitaku $ kodukai carryer $ tsuuwa; if sex^='M' & sex^='F' then delete; if shintyou=. | taijyuu=. then delete; proc print data=gakusei(obs=10); run; proc plot data=gakusei; : 散布図 plot shintyou*taijyuu; : 元の変量のプロット run; : proc princomp cov data=gakusei out=outprin; : 主成分分析(分散共分散行列) var shintyou taijyuu; : 2変量 run; : proc print data=outprin(obs=15); : 結果の出力 run; : proc plot data=outprin; : 散布図 plot prin2*prin1/vref=0 href=0; : 主成分得点のプロット run; : : 参考までに、 proc sort data=outprin; : 説明のためにソートしてみる by prin1; : 第一主成分で run; : proc print data=outprin; : 体重がややが効いていることの確認 run; :

出力結果 : SAS_out0511.txt, les0511_out.pdf

身長と体重の散布図
各変量の平均、標準偏差、分散共分散行列
固有値、比率(寄与率)、累積寄与率 : 解釈に使う
固有ベクトル(係数a1とa2) : 解釈に使う
主成分得点 : 各個人の得点(z)、2つある
第1軸と第2軸の主成分得点の散布図

2017年11月 9日木曜日 11時13分08秒 57 プロット : shintyou*taijyuu 凡例 : A = 1 obs, B = 2 obs, ... shintyou | 200 + | | A B A A 180 + A BFCCGCD H B B A A A | CFCJHHUQIKCIBFBB BC A | H DHFMDFBCBGA A AA A A 160 + ABDFEGKADBACB | A EADECAE A A | A BAA 140 + ---+------------+------------+------------+------------+-- 20 40 60 80 100 taijyuu 2017年11月 9日木曜日 11時13分08秒 58 PRINCOMP プロシジャオブザベーション数 325 変数の数 2 単純統計量 shintyou taijyuu Mean 168.6947692 58.78092308 StD 8.0436274 9.34693152 2017年11月 9日木曜日 11時13分08秒 59 PRINCOMP プロシジャ共分散行列 shintyou taijyuu shintyou 64.69994169 52.97829497 taijyuu 52.97829497 87.36512877 総分散 152.06507047 共分散行列の固有値固有値差比率累積 1 130.209351 108.353632 0.8563 0.8563 2 21.855719 0.1437 1.0000 固有ベクトル Prin1 Prin2 shintyou 0.628817 0.777553 taijyuu 0.777553 -.628817 2017年11月 9日木曜日 11時13分08秒 61 s h t k c i a j o a t n i k i d r s P P t j y t u r u r r O s y y o a k y u i i B e o u u k a e w n n S x u u i u i r a 1 2 1 F 145.0 38.0 . J 10000 . -31.0580 -5.35654 2 F 146.7 41.0 85 J 10000 Vodafone 6000 -27.6563 -5.92115 3 F 148.0 42.0 . J 50000 . -26.0613 -5.53915 4 F 148.0 43.0 80 J 50000 DoCoMo 4000 -25.2837 -6.16797 5 F 149.0 45.0 . G 60000 . -23.0998 -6.64805 6 F 150.0 46.0 86 40000 . -21.6934 -6.49931 7 F 151.0 45.0 . J 20000 docomo 5000 -21.8422 -5.09294 8 F 151.0 50.0 . G 60000 J-PHONE . -17.9544 -8.23703 9 F 151.7 41.5 80 J 35000 . -24.1234 -2.34780 10 F 152.0 35.0 77 J 60000 DoCoMo 2000 -28.9889 1.97278 11 F 152.0 43.0 . J 20000 au 3500 -22.7685 -3.05776 12 F 152.0 44.0 . 45000 DoCoMo 4000 -21.9909 -3.68657 13 F 153.0 41.0 . J 125000 No . -23.6948 -1.02257 14 F 153.0 42.0 . G 0 Vodafone 1000 -22.9172 -1.65139 15 F 153.0 46.5 87 G 10000 . -19.4182 -4.48106 2017年11月 9日木曜日 11時13分08秒 63 プロット : Prin2*Prin1 凡例 : A = 1 obs, B = 2 obs, ... 20 + | | | Prin2 | A | A | B B DABBABACABA AAA | AB DDAADDBEEECCHEFCABDD A A 0 +-----------A---AABBABCAFC-FE-BGGBFEEEKFDACCEB--G-AA------A----------- | A AA ABABABBDC BAF B ABAB CBDDBACED A AA | A AA AA A A B B|BA A A A AB A | A | AA A A | | A A -20 + | A A ---+-------+-------+-------+-------+-------+-------+-------+-------+-- -40 -30 -20 -10 0 10 20 30 40 Prin1 2017年11月 9日木曜日 11時13分08秒 64 s h t k c i a j o a t n i k i d r s P P t j y t u r u r r O s y y o a k y u i i B e o u u k a e w n n S x u u i u i r a 1 2 1 F 145.0 38.0 . J 10000 . -31.0580 -5.35654 2 F 152.0 35.0 77 J 60000 DoCoMo 2000 -28.9889 1.97278 3 F 146.7 41.0 85 J 10000 Vodafone 6000 -27.6563 -5.92115 ≪中略≫ 2017年11月 9日木曜日 11時13分09秒 89 OBS sex shintyou taijyuu kyoui jitaku kodukai carryer tsuuwa Prin1 Prin2 292 M 177.0 68 . G 80000 . 12.3908 0.6606 293 M 182.0 64 . G 0 . 12.4247 7.0637 294 M 165.0 78 . G 0 2098 12.6205 -14.9582 295 M 170.0 74 90 J 0 . 12.6544 -8.5551 296 M 175.0 70 95 G 50000 8000 12.6883 -2.1521 297 M 178.0 68 . J 100000 DoCoMo 4000 13.0196 1.4382 298 M 184.0 65 . G 140000 au 10000 14.4599 7.9900 299 M 170.0 78 . 45000 Vodafone 10000 15.7646 -11.0704 300 M 179.9 70 . J 15000 DoCoMo 700 15.7695 1.6579 301 M 175.0 74 . J 0 . 15.7985 -4.6674 302 M 180.0 70 94 G 70000 au 5000 15.8324 1.7357 303 M 180.0 70 . J 40000 au 4000 15.8324 1.7357 304 M 180.0 70 . . . 15.8324 1.7357 305 M 180.0 70 . J 40000 DoCoMo 6500 15.8324 1.7357 306 M 180.0 70 . 5000 3000 15.8324 1.7357 307 M 178.7 71.2 95 0 . 15.9480 -0.0297 308 M 184.0 68.0 85 30000 . 16.7925 6.1035 309 M 173.5 76.5 . G 100000 . 16.7991 -7.4057 310 M 182.0 70.0 90 G 100000 . 17.0900 3.2908 311 M 185.0 68.0 93 J 0 . 17.4213 6.8811 312 M 175.0 77.0 95 G 130000 . 18.1311 -6.5538 313 M 179.1 74.2 . 0 au 4000 18.5321 -1.6052 314 M 175.0 79.0 . J 0 No 0 19.6862 -7.8115 315 M 176.5 78.0 96 J 10000 . 19.8519 -6.0163 316 M 177.0 78.0 . J 40000 . 20.1663 -5.6275 317 M 181.5 74.5 . G 120000 au 3000 20.2746 0.0723 318 M 185.0 72.0 . J 30000 7000 20.5315 4.3658 319 M 178.0 78.0 110 G 50000 . 20.7951 -4.8500 320 M 173.0 84.0 46 G 350000 . 22.3164 -12.5106 321 M 169.3 88.5 94 J 0 . 23.4887 -18.2173 322 M 186 82 . J 0 . 28.9359 -1.1448 323 M 182 90 100 J 40000 . 32.6411 -9.2856 324 M 178 95 . 1000 No . 34.0135 -15.5399 325 M 178 100 112 G 60000 . 37.9013 -18.6840

解釈方法

寄与率 : その軸がどの程度説明力を持っているか : 第1軸だけで十分(85.6%)。第2軸に含まれる説明力は小さい(14.4%)。
固有ベクトル : その軸の特徴を示している : 身長と体重の重みはほぼ同等だが、体重がやや大きめに効いている(第1軸)
主成分得点と散布図 : 各個人がどこに付置されているか
第1軸 : 全体的な体格の指標。身長と体重を足したような指標。
3変量以上の主成分分析
1. 身長、体重、胸囲での総合指標を求めてみよう : プログラム : les1302.sas
  /* Lesson 05-12 */ /* File Name = les0512.sas 11/09/17 */ options linesize=72 pagesize=20; options nocenter linesize=78 pagesize=20; proc printto log = '/folders/myfolders/Kougi/SAS_log0512.txt' print = '/folders/myfolders/Kougi/SAS_out0512.txt' new; ods listing gpath='/folders/myfolders/Kougi/SAS_ODS99'; data gakusei; infile '/folders/myfolders/Kougi/all07ae.prn' firstobs=2; input sex $ shintyou taijyuu kyoui jitaku $ kodukai carryer $ tsuuwa; if sex^='M' & sex^='F' then delete; if shintyou=. | taijyuu=. then delete; proc print data=gakusei(obs=10); run; proc princomp cov data=gakusei out=outprin; : 主成分分析(分散共分散行列) var shintyou taijyuu kyoui; : 3変量 run; : proc print data=outprin(obs=15); : 結果の出力 run; : proc plot data=outprin; : 散布図 plot prin2*prin1/vref=0 href=0; : 主成分得点のプロット plot prin3*prin2/vref=0 href=0; : plot prin3*prin1/vref=0 href=0; : run; :
2. 出力結果 : SAS_out0512.txt, les0512_out.pdf
  - 各変量の平均、標準偏差、共分散行列
  - 固有値、比率(寄与率)、累積寄与率
  - 固有ベクトル
  - 主成分得点
  - 第1軸～第3軸の散布図
```
                               2017年11月 9日 木曜日 11時15分29秒  94
PRINCOMP プロシジャ

オブザベーション数         114
変数の数                          3

                     単純統計量
              shintyou           taijyuu             kyoui
Mean       167.3517544       58.79298246       86.17543860
StD          8.7227627       10.86282708        8.36262822

                       共分散行列
                  shintyou           taijyuu             kyoui
shintyou        76.0865898        69.6653222        23.7439373
taijyuu         69.6653222       118.0010123        43.5906226
kyoui           23.7439373        43.5906226        69.9335507

総分散    264.02115277

               共分散行列の固有値
         固有値           差      比率      累積
   1    189.966471    138.636164      0.7195      0.7195
   2     51.330307     28.605932      0.1944      0.9139
   3     22.724375                    0.0861      1.0000

                固有ベクトル
                 Prin1         Prin2         Prin3
shintyou      0.539085      -.407903      0.736887
taijyuu       0.751825      -.161336      -.639320
kyoui         0.379667      0.898658      0.219698

                               2017年11月 9日 木曜日 11時15分29秒  98
        s
        h     t          k      c
        i     a     j    o      a       t
        n     i   k i    d      r       s      P       P        P
        t     j   y t    u      r       u      r       r        r
 O  s   y     y   o a    k      y       u      i       i        i
 B  e   o     u   u k    a      e       w      n       n        n
 S  x   u     u   i u    i      r       a      1       2        3

  1 F 145.0 38.0  . J  10000             .    .       .       .     
  2 F 146.7 41.0 85 J  10000 Vodafone 6000 -24.9565 10.2382 -4.10085
  3 F 148.0 42.0  . J  50000             .    .       .       .     
  4 F 148.0 43.0 80 J  50000 DoCoMo   4000 -24.6504  4.8920 -5.52002
  5 F 149.0 45.0  . G  60000             .    .       .       .     
  6 F 150.0 46.0 86    40000             . -19.0388  8.9841 -4.64602
  7 F 151.0 45.0  . J  20000 docomo   5000    .       .       .     
  8 F 151.0 50.0  . G  60000 J-PHONE     .    .       .       .     
  9 F 151.7 41.5 80 J  35000             . -23.7835  3.6248 -1.83456
 10 F 152.0 35.0 77 J  60000 DoCoMo   2000 -29.6477  1.8551  1.88299
 11 F 152.0 43.0  . J  20000 au       3500    .       .       .     
 12 F 152.0 44.0  .    45000 DoCoMo   4000    .       .       .     
 13 F 153.0 41.0  . J 125000 No          .    .       .       .     
 14 F 153.0 42.0  . G      0 Vodafone 1000    .       .       .     
 15 F 153.0 46.5 87 G  10000             . -16.6659  8.5784 -2.53532

                               2017年11月 9日 木曜日 11時15分29秒 100

       プロット : Prin2*Prin1   凡例 : A = 1 obs, B = 2 obs, ...
       (NOTE: 211 obs が欠損値です。)
Prin2 |                            |
   20 +                            |
      |           A                |  A    A           A            A
      |           B   BAAABCA   AC |  A AC     A
    0 +-------A----A--ABACCBDD-BCA-AD-CG-ABAFBA-C-AA-A-------A----------------
      |              A  AA      A  A C AAAAAA AA
      |                       AA   | A A A
  -20 +                            |
      |                  A         |
      |                            |
  -40 +                            |    A
      -+-------------+-------------+-------------+-------------+-------------+
      -40           -20            0            20            40            60
                                        Prin1

                               2017年11月 9日 木曜日 11時15分29秒 101
       プロット : Prin3*Prin2   凡例 : A = 1 obs, B = 2 obs, ...
       (NOTE: 211 obs が欠損値です。)
Prin3 |                                                  |
   10 +                                         A   A A  |
      |                                      AA A  AAAB  | B
      |                                         AA   C BCEEDA  B
    0 +-----------------------------------A----B-A---AA-BEBGDCAADA-----AA-----
      |                       A                    A  ACBAAEACA  ABA   A
      |                                          A A     AA   B
  -10 +                                                 A|       A
      |                                                  |           A
      |                                                  |A
  -20 +        A                                         |
      -+---------+---------+---------+---------+---------+---------+---------+
      -50       -40       -30       -20       -10        0        10        20
                                        Prin2

                               2017年11月 9日 木曜日 11時15分29秒 102
       プロット : Prin3*Prin1   凡例 : A = 1 obs, B = 2 obs, ...
       (NOTE: 211 obs が欠損値です。)
Prin3 |                            |
   10 +                            | A   A      A
      |                            |AB AA   CAA
      |       A       A  BB BB  B  |AABD BAA   AB
    0 +---------------ABB-CDBAABBB-AB-AC-BAAB-AA-------A----------------------
      |           BA  AAAC AAA A   A  BAAAAAB     AA
      |           A  A          AA |         A               A
  -10 +                  A       A |
      |                            |                                A
      |                            |                 A
  -20 +                            |    A
      -+-------------+-------------+-------------+-------------+-------------+
      -40           -20            0            20            40            60
                                        Prin1
```
3. 解釈方法
  - 寄与率 : 2軸まで取れば十分のようだ(91.4%)。
  - 第1軸 : 全体的な体格の因子。特に体重が効いている。
  - 第2軸 : 太さの因子(?)。胸囲が正で身長が負。
  - 第3軸 : 華奢さの因子(?)。無視しても良い軸であるが。(8.6%)。

相関行列を使う理由

変量によって測定単位が異なる
使う単位によって解析結果が異なる
分散の影響を受ける
各変量を標準化して用いる : 相関行列 <===> 分散共分散行列

相関行列を用いて体格の総合指標を求めてみよう : プログラム : les1303.sas

 /* Lesson 05-13 */
 /*    File Name = les0513.sas   11/09/17   */
options linesize=72 pagesize=20;
options nocenter linesize=78 pagesize=20;

proc printto log   = '/folders/myfolders/Kougi/SAS_log0513.txt'
             print = '/folders/myfolders/Kougi/SAS_out0513.txt' new;

ods listing gpath='/folders/myfolders/Kougi/SAS_ODS99';

data gakusei;
  infile '/folders/myfolders/Kougi/all07ae.prn'
    firstobs=2;
  input sex $ shintyou taijyuu kyoui 
        jitaku $ kodukai carryer $ tsuuwa;

if sex^='M' & sex^='F' then delete;
if shintyou=. | taijyuu=. then delete;

proc print data=gakusei(obs=10);

run;                                          :
proc princomp data=gakusei out=outprin;       : 相関係数を使って
  var shintyou taijyuu kyoui;                 :
run;                                          :
proc print data=outprin(obs=15);
run;
proc plot data=outprin;
  plot prin2*prin1/vref=0 href=0;
  plot prin3*prin2/vref=0 href=0;
  plot prin3*prin1/vref=0 href=0;
run;

出力結果 : SAS_out0513.txt, les0513_out.pdf

各変量の平均、標準偏差、相関行列
固有値、比率(寄与率)、累積寄与率
固有ベクトル
主成分得点
第1軸～第3軸の散布図


                               2017年11月 9日 木曜日 11時16分02秒 143
PRINCOMP プロシジャ

オブザベーション数         114
変数の数                          3

                     単純統計量
              shintyou           taijyuu             kyoui
Mean       167.3517544       58.79298246       86.17543860
StD          8.7227627       10.86282708        8.36262822

                 相関行列
              shintyou      taijyuu       kyoui
shintyou        1.0000       0.7352      0.3255
taijyuu         0.7352       1.0000      0.4799
kyoui           0.3255       0.4799      1.0000

                相関行列の固有値
         固有値           差      比率      累積
   1    2.04696555    1.33664820      0.6823      0.6823
   2    0.71031735    0.46760025      0.2368      0.9191
   3    0.24271710                    0.0809      1.0000

                固有ベクトル
                 Prin1         Prin2         Prin3
shintyou      0.599200      -.483881      0.637823
taijyuu       0.640769      -.187770      -.744418
kyoui         0.479974      0.854752      0.197544

解釈方法
- 寄与率 : 2軸まで取れば十分のようだ(91.9%)。
- 第1軸 : 全体的な体格の因子。
- 第2軸 : 太さの因子。
- 分散共分散行列を使ったときよりも第1軸の固有ベクトル同士が近い値になった。しかし、軸の解釈に違いはない。その理由はこの例では 3変量のスケールや分散にそれほどの違いがないためと想像される。
- 分散共分散行列と相関行列を使ったときの違いを見てみたければ、 shintyou のみを mm 単位で測定したと考えて、 100倍したものをデータとして両者の出力を比較してみよ。
  プログラム : les0514.sas 、出力結果 : SAS_out0514.txt, les0514_out.pdf

主成分の数の決定基準 : 配布資料 80ページ
明確に決まっているわけではないが、以下のような基準が一般的に用いられている。また、結果の解釈の都合上、多少増減させることもある。
- 累積寄与率がある程度(例えば80%以上)以上大きくなること
- 相関行列を用いている場合、固有値が 1以上
- 固有値のギャップがある程度以上大きいこと
- ...
次回は、... : 11月16日
- 因子分析
- Reasoning Test, 思考テスト
- 学習診断、スコアリングレポート
- ...

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