<4D F736F F F696E74202D BD95CF97CA89F090CD F6489F18B4195AA90CD816A>

Size: px

Start display at page:

Download "<4D F736F F F696E74202D BD95CF97CA89F090CD F6489F18B4195AA90CD816A>"

まいかおいもり
5 years ago
Views:

1 主な多変量解析 9. 多変量解析 1 ( 重回帰分析 ) 目的変数量的説明変数質的あり量的重回帰分析数量化 Ⅰ 類質的判別分析数量化 Ⅱ 類なし主成分分析因子分析多次元尺度構成法数量化 Ⅲ 類数量化 Ⅳ 類その他クラスタ分析共分散構造分析説明変数 : 独立変数予測変数目的変数 : 従属変数基準変数 3 1. 単回帰分析各データの構造 y b ax a α: 1,,, n (nはサンプル数) a: 回帰係数 b: 切片 ε: 誤差 4 回帰式 yˆ b ax 身長体重身長体重身長から体重を予測する height <- c(166,173,177,160,174,175,169,169) weight <- c(66,70,76,65,74,75,74,7) taikei <- data.frame(height,weight) colnames(taikei) <- c(" 身長 ", " 体重 ") rownames(taikei) <- 1:8 taikei.lm <- lm( 体重 ~ 身長, taikei) summary(taikei.lm) 1

2 7 モデル式 8 lm(formula = 体重 ~ 身長, data = taikei) Residuals: Min 1Q Median 3Q Max (Intercept) 身長 ** Signif. codes: 0 *** ** 0.01 * Residual standard error:.78 on 6 degrees of freedom Multiple R-squared: , Adjusted R-squared: F-statistic: 17.1 on 1 and 6 DF, p-value: lm(formula = 体重 ~ 身長, data = taikei) Residuals: Min 残差の基本統計量 1Q Median 3Q Max (Intercept) 身長 ** Signif. codes: 0 *** ** 0.01 * Residual standard error:.78 on 6 degrees of freedom Multiple R-squared: , Adjusted R-squared: par(cex=1.8, pch=16, lwd=3) plot(taikei, ylim=c(60,85), cex=1.5, col="red", xlab=" 身長 (cm)", ylab=" 体重 (Kg) ") abline(taikei.lm) 9 個々の残差の表示 taikei.lm$residuals でも同じ 10 残差 y a yˆ > ( 残差 <- residuals(taikei.lm)) C 標準誤差 t 値 p 値 (Intercept) 身長 ** Signif. codes: 0 *** ** 0.01 * 切片回帰係数 Residual standard error:.78 on 6 degrees of freedom Multiple R-squared: , Adjusted R-squared: F-statistic: 17.1 on 1 and 6 DF, p-value: 回帰係数残差の乗和を最小にする a と b ( x Sxy x)( y y) a S ( x x) xx b y ax 回帰式 Weight Height 1

3 回帰係数の標準誤差 t 値残差平方和残差不偏分散 S e s e y yˆ ( ) Se n k 1 n : 標本数 k : 説明変数の数 13 切片の標準誤差回帰係数の t 値 SE( b) se ( x x) t a a SE(a) 14 回帰係数の標準誤差 SE( a) 1 se n x ( x x) 切片の t 値 t b b SE(b) 推定値の標準誤差 = 残差の標準偏差 Residual standard error:.78 on 6 degrees of freedom Multiple R-squared: , Adjusted R-squared: F-statistic: 17.1 on 1 and 6 DF, p-value: 決定係数 R ( 寄与率 ) 調整済み決定係数 R モデル全体のF 値回帰の変動 0.8<R かなりよい目的変数の変動 0.5<R <0.8 まあよい 15 決定係数 R 調整済み決定係数 R ~ S S R T R ~ S 1 1 e n k R S ( n 1) T n : 標本数 k : 説明変数の数 S T ( y y) ( yˆ y) S R 16 S e ( y yˆ) S T S R S e F 値 R F 1 R n k 1 k 17 予測値残差予測 <- predict(taikei.lm) data.frame(taikei, 予測, 残差 ) 身長体重予測残差

4 19 0 回帰診断図残差と予測値 ( フィット値 ) par(mfrow=c(,)) 残差と予測値残差の正規 Q-Q plot(taikei.lm, col="red", which=1) plot(taikei.lm) 残差の平方根てこ比標準化残差残差の正規 Q-Q 1 残差の平方根 plot(taikei.lm, col="red", which=) plot(taikei.lm, col="red", which=3) Cook の距離 3 てこ比標準化残差 4 plot(taikei.lm, col="red", which=4) plot(taikei.lm, col="red", which=5) 4

5 例題 ) データセット cars を使って速度から停止距離を回帰分析せよ > (cars) 5. 重回帰分析各データの構造 ( 説明変数が i 個 ) y b a x a x a x i 1 1 i a 6 speed dist 組み込みのデータセットは > data() で一覧を見ることができるさらにヘルプで詳細が表示される α: 1,,, n (nはサンプル数) a i : 偏回帰係数 b: 切片 ε: 誤差重回帰式 y b a1x1 ax a i x i 誤差の独立性不偏性等分散性正規性を仮定身長体重ウェスト胸囲相関行列 round(cor(taikei),3) 身長体重ウェスト胸囲身長体重ウェスト胸囲 height <- c(166,173,177,160,174,175,169,169) weight <- c(66,70,76,65,74,75,74,7) waist <- c(79,78,80,78,79,81,84,81) chest <- c(83,80,86,81,84,84,91,85) taikei <- data.frame(height,weight,waist,chest) colnames(taikei) <- c(" 身長 ", " 体重 ", " ウェスト ", " 胸囲 ") rownames(taikei) <- 1:8 対散布図 9 30 pairs(taikei, pch=1, bg="red", cex=, panel=panel.smooth) 5

6 31 3 重回帰分析 (taikei.lm <- lm( 体重 ~.,data=taikei)) summary(taikei.lm) lm(formula = 体重 ~., data = taikei) lm(formula = 体重 ~., data = taikei) (Intercept) 身長ウェスト胸囲 Residuals: (Intercept) 身長 ** ウェスト胸囲 Signif. codes: 0 *** ** 0.01 * Residual standard error: on 4 degrees of freedom Multiple R-squared: 0.956, Adjusted R-squared: F-statistic: 16.6 on 3 and 4 DF, p-value: lm(formula = 体重 ~., data = taikei) Residuals: 個々のデータの残差 (Intercept) 身長 ** ウェスト胸囲 Signif. codes: 0 *** ** 0.01 * Residual standard error: on 4 degrees of freedom 33 偏回帰係数標準誤差 t 値 p 値 (Intercept) 身長 ** ウェスト胸囲 Signif. codes: 0 *** ** 0.01 * 推定値の標準誤差 = 残差の標準偏差 Residual standard error: on 4 degrees of freedom Multiple R-squared: 0.956, Adjusted R-squared: F-statistic: 16.6 on 3 and 4 DF, p-value: モデル全体の F 値決定係数 R 調整済み決定係数 R 34 偏回帰プロット install.packages("car") library(car) avplot(taikei.lm," 身長 ") avplot(taikei.lm," ウェスト ") avplot(taikei.lm," 胸囲 ") 35 説明変数の選択予測精度を低めずに必要最小限の説明変数を選択する AIC= )+(q+1) 36 AIC が小さいほど良いモデル 6

7 37 38 (taikei.lm <- step(taikei.lm)) Start: AIC=8.87 体重 ~ 身長 + ウェスト + 胸囲 - ウェスト胸囲 <none> 身長 Step: AIC=7.4 体重 ~ 身長 + 胸囲 <none> 胸囲身長 lm(formula = 体重 ~ 身長 + 胸囲, data = taikei) (Intercept) 身長胸囲 Start: AIC=8.87 体重 ~ 身長 + ウェスト + 胸囲 - ウェスト胸囲 <none> 身長 Step: AIC=7.4 体重 ~ 身長 + 胸囲 <none> 胸囲身長 summary(taikei.lm) 40 lm(formula = 体重 ~ 身長 + 胸囲, data = taikei) (Intercept) 身長胸囲回帰式体重身長胸囲 lm(formula = 体重 ~ 身長 + 胸囲, data = taikei) Residuals: (Intercept) * 身長 ** 胸囲 * Signif. codes: 0 *** ** 0.01 * Residual standard error: on 5 degrees of freedom Multiple R-squared: 0.903, Adjusted R-squared: F-statistic: 8.89 on and 5 DF, p-value: 予測値残差残差 <- residuals(taikei.lm) 予測 <- predict(taikei.lm) data.frame(taikei, 予測, 残差 ) 身長体重ウェスト胸囲予測残差予測平面 install.packages("scatterplot3d") library(scatterplot3d) s3d <- scatterplot3d(cbind(height,chest,weight), type="h", pch=16) s3d$plane3d(taikei.lm, lty.box ="solid") 4 7

8 VIF:Variance Inflation factor( 分散拡大係数 ) 43 VIF 44 VIFが4~10 程度以上は多重共線性が疑われるトレランス ( 許容度 )=1/VIF install.packages("car") library(car) vif(taikei.lm) 多重共線性 (multicolinearlity) 身長ウェスト胸囲説明変数間に関係が強い場合 vif(taikei.lm) 解が求まらないあるいは安定しない身長胸囲相関行列 45 標準回帰係数 46 round(cor(taikei),3) 測定値を標準得点化身長体重ウェスト胸囲身長体重ウェスト胸囲 staikei <- data.frame(scale(taikei)) staikei 身長体重ウェスト胸囲 staikei.lm <- lm( 体重 ~., data=staikei) staikei.lm <- step(staikei.lm) Start: AIC= 体重 ~ 身長 + ウェスト + 胸囲 - ウェスト胸囲 <none> 身長 Step: AIC= 体重 ~ 身長 + 胸囲 <none> 胸囲身長 summary(staikei.lm) lm(formula = 体重 ~ 身長 + 胸囲, data = staikei) Residuals: 標準回帰係数 (Intercept).e e 身長 7.703e e ** 胸囲 4.336e e * Signif. codes: 0 *** ** 0.01 * 0.05 Residual standard error: on 5 degrees of freedo 48 8

9 標準化残差標準化残差 <=3 標準化残差 <- residuals(taikei.lm) 標準化予測 <- predict(taikei.lm) data.frame(taikei, 標準化予測, 標準化残差 ) 身長体重ウェスト胸囲標準化予測標準化残差残差分析誤差分散の推定値 See <- sum(residuals(taikei.lm)^) Ve <- See/(8-1-1) [1] データ数ー説明変数の数ー 1 標準化残差 sr <- residuals(taikei.lm)/sqrt(ve) plot(taikei$height,sr,xlab="height",ylab="standardized residuals") abline(h=0) plot(taikei$waist,sr,xlab="waist",ylab="standardized residuals") abline(h=0) plot(taikei$chest,sr,xlab="chest",ylab="standardized residuals") abline(h=0) 偏回帰プロット avplot(taikei.lm," 身長 ") avplot(taikei.lm," 胸囲 ") 5 回帰診断図 par(mfrow=c(,)) plot(taikei.lm) 残差と予測値残差の正規 Q-Q 53 青木先生の mreg を使う方法目的変数は最右列におく 54 source("../r/all.r", encoding="euc-jp") taikei3 <- data.frame(height, waist, chest, weight) colnames(taikei3) <- c(" 身長 ", " ウェスト ", " 胸囲 ", " 体重 ") rownames(taikei3) <- 1:8 (a <- mreg(taikei3)) 残差の平方根てこ比標準化残差 9

10 偏回帰係数標準誤差 t 値 P 値標準化偏回帰係数トレラ身長ウェスト胸囲定数項回帰の分散分析表平方和自由度平均平方 F 値 P 値回帰残差全体重相関係数 = 重相関係数の二乗 = 自由度調整済重相関係数の二乗 = 対数尤度 = AIC = (sreg(taikei3)) 有効ケース数 : 8 従属変数 : 体重平均値不偏分散標準偏差身長ウェスト胸囲体重 ***** 相関係数行列 ***** 身長ウェスト胸囲体重身長ウェスト胸囲体重変数編入基準 Pin: 0.05 変数除去基準 Pout: 0.05 編入候補変数 : 身長 ***** ステップ 1 ***** 編入変数 : 身長 P : ***** 編入されました 56 偏回帰係数標準誤差 t 値 P 値標準化偏回帰係数トレランス分散拡大要因途中省略 ===================== 結果 ===================== 偏回帰係数標準誤差 t 値 P 値標準化偏回帰係数トレランス分散拡大要因身長胸囲定数項平方和自由度平均平方 F 値 P 値回帰残差全体重相関係数決定係数 ( 重相関係数の二乗 ) 自由度調整済み重相関係数の二乗対数尤度 AIC $ 分析結果偏回帰係数標準誤差 t 値 P 値標準化偏回帰係数トレランス分散拡大要因身長胸囲定数項 $ 分散分析表平方和自由度平均平方 F 値 P 値回帰残差全体 $ 決定係数重相関係数決定係数 ( 重相関係数の二乗 ) 自由度調整済み重相関係数の二乗対数尤度 AIC $ 予測観察値予測値残差標準化残差重回帰分析をする前の注意点説明変数は 10 個未満とするデータ数は説明変数の数より 0 以上多い 60 説明変数の数 < データ数 /3 各説明変数は目的変数と直線相関がある説明変数間の相関はあまり大きすぎない 10

11 重回帰分析をした後の注意点疑似相関相関係数が大きいが標準偏回帰係数が小さい多重共線性相関係数と標準偏回帰係数が異符号で大きい説明変数間の相関が高い場合が多い 61 例題 data(attitude) attitude rating 会社に対する総合評価 complaints 従業員の苦情の取り扱い privileges 特別な特権は許さない learning 学習の機会 raises 能力に基づいた昇給 critical 加重 advance 昇進 ratingを他の変数で予測 attitude <- cbind(attitude[,:7],rating=attitude[,1]) 6 相関行列 round(cor(attitude),) rating complaints privileges learning raises critical advance rating complaints privileges learning raises critical advance pairs(attitude,panel=panel.smooth,attitude) , 68, 77, 81, 74, 65, 6, 83, 77, 90, 85, 60 57, 83, 54, 50, 64, 65 4, 45, 7, 7, 69, 75, 55, 59, 79, 60, 79, , 77, 77, 54, 79, 80, 5, 35, 46, 36, 63, attitude <- cbind(attitude[,:7],rating=attitude[,1]) (a <- mreg(attitude)) 65 偏回帰係数標準誤差 t 値 P 値標準化偏回帰係数トレランス complaints privileges learning raises critical advance 定数項回帰の分散分析表平方和自由度平均平方 F 値 P 値回帰 e-05 残差全体重相関係数 = 重相関係数の二乗 = 自由度調整済重相関係数の二乗 = 対数尤度 = AIC =

12 回帰分析の要約情報 attitude.lm1 <- lm(rating ~., data = attitude) summary(attitude.lm1) 67 lm(formula = rating ~., data = attitude) Residuals: Min 1Q Median 3Q Max (Intercept) complaints *** privileges learning raises critical advance Residual standard error: on 3 degrees of freedom Multiple R-squared: 0.736, Adjusted R-squared: F-statistic: 10.5 on 6 and 3 DF, p-value: 1.40e e AIC によるモデルの変数選択 AIC=- ( モデルの最大対数尤度 )+ ( モデルのパラメータ数 ) attitude.lm<-step(attitude.lm1) 70 Start: AIC=13.36 rating ~ complaints + privileges + learning + raises + critical + advance - critical raises privileges advance <none> learning complaints Step: AIC=11.45 rating ~ complaints + privileges + learning + raises + advance - raises privileges advance <none> learning complaints

13 自動車感性評価学 Step: AIC= rating ~ complaints + privileges + learning + advance - privileges advance <none> learning complaints Step: AIC= rating ~ complaints + learning + advance - advance <none> learning complaints Step: AIC=118 rating ~ complaints + learning 74 <none> learning complaints 最終結果 summary(attitude.lm) 75 lm(formula = rating ~ complaints + learning, data = attitude) 76 Residuals: Min 1Q Median 3Q Max (Intercept) complaints e-06 *** learning 回帰診断図 78 par(mfrow=c(,)) Residual standard error: on 7 degrees of freedom Multiple R-squared: 0.708, Adjusted R-squared: F-statistic: 3.74 on and 7 DF, p-value: 6.058e-08 plot(attitude.lm) Residuals Residuals vs Fitted 113 Standardized residuals Normal Q-Q rating = *complaints+0.11*learning Fitted values Theoretical Quantiles Scale-Location Residuals vs Leverage Standardized residuals Standardized residuals Cook's distance Fitted values Leverage 13

14 例題 airquality 年 5 月 ~9 月のニューヨーク Ozone オゾンの量 (ppb) Solar.R 太陽の放射の量 (lang) Wind 風力 (mph) Temp 温度 ( 華氏 F) Month 月 1 1 Day 月の日数 1 31 Ozone を他の変数で予測 14

Microsoft PowerPoint - 資料04 重回帰分析.ppt

Microsoft PowerPoint - 資料04 重回帰分析.ppt 04. 重回帰分析京都大学加納学 Division of Process Control & Process Sstems Engineering Department of Chemical Engineering, Koto Universit manabu@cheme.koto-u.ac.jp http://www-pse.cheme.koto-u.ac.jp/~kano/ Outline