重回帰で回帰係数の符号が逆転している場合の意味 [データサイエンス、統計モデル]
個々の変数の相関が正（負）なのに、重回帰式において、その係数が負（正）になることがあります。
そんな時に、どう考えるか？
試しにネットに落ちている中古車データを拾ってきました。
相関係数を計算すると
となります。
走行距離や年数が古くなれば、中古車価格は下落するという事が分かります。
ここで、中古車価格を目的変数に、走行距離と年数を説明変数にすると
中古車価格 = -10.0 × 走行距離 - 10.7 × 年数 + 252.0
という重回帰式が得られます。
走行距離が1万km増えると、中古車価格は10.0万円下落し、
同様に年数が1年増えると、中古車価格は10.7万円下落することが分かります。
さて、次に、年数を目的変数に、走行距離と中古車価格を説明変数にすると
年数 = -0.0499 × 価格 - 0.154× 走行距離 + 13.86
という結果になりました。
年数と走行距離の相関係数は、0.562と正の値なのですが、
回帰式のパラメータは、マイナスとなっています。
つまり、「走行距離が1万キロ増えると、年数が0.154年減少する」となっていて感覚と値が違う結果になってしまいます。
これをどう考えるか？です。
重回帰の場合、解釈としては、注意が必要で、「他の変数が同じ値だった場合（他の変数を影響を統制した場合）の効果」をみていることになります。
つまり、「走行距離が同じであったとしたら、走行距離が1万キロ増えた場合、年数はどうなるか？」を表していることになります。
例として、ある車種が、
【車種 A】走行距離：10万キロ、価格：100万円、年数：10年
だったとした場合、
【車種 B】走行距離：20万キロ、価格：100万円、年数：？？？年
車種Bは、20万キロも走っているけど、車種Aと同じです。
つまり、年数が10年ではなく、もっと小さい数字であることが期待されます。
この構造を表しているのが重回帰のパラメータということになります。
一方、一般的には、走行距離が大きくなれば、年数も同時に大きくなる（正の相関がある）という事を表しているのが、相関係数となります。
そんな時に、どう考えるか？
試しにネットに落ちている中古車データを拾ってきました。
相関係数を計算すると
となります。
走行距離や年数が古くなれば、中古車価格は下落するという事が分かります。
ここで、中古車価格を目的変数に、走行距離と年数を説明変数にすると
中古車価格 = -10.0 × 走行距離 - 10.7 × 年数 + 252.0
という重回帰式が得られます。
走行距離が1万km増えると、中古車価格は10.0万円下落し、
同様に年数が1年増えると、中古車価格は10.7万円下落することが分かります。
さて、次に、年数を目的変数に、走行距離と中古車価格を説明変数にすると
年数 = -0.0499 × 価格 - 0.154× 走行距離 + 13.86
という結果になりました。
年数と走行距離の相関係数は、0.562と正の値なのですが、
回帰式のパラメータは、マイナスとなっています。
つまり、「走行距離が1万キロ増えると、年数が0.154年減少する」となっていて感覚と値が違う結果になってしまいます。
これをどう考えるか？です。
重回帰の場合、解釈としては、注意が必要で、「他の変数が同じ値だった場合（他の変数を影響を統制した場合）の効果」をみていることになります。
つまり、「走行距離が同じであったとしたら、走行距離が1万キロ増えた場合、年数はどうなるか？」を表していることになります。
例として、ある車種が、
【車種 A】走行距離：10万キロ、価格：100万円、年数：10年
だったとした場合、
【車種 B】走行距離：20万キロ、価格：100万円、年数：？？？年
車種Bは、20万キロも走っているけど、車種Aと同じです。
つまり、年数が10年ではなく、もっと小さい数字であることが期待されます。
この構造を表しているのが重回帰のパラメータということになります。
一方、一般的には、走行距離が大きくなれば、年数も同時に大きくなる（正の相関がある）という事を表しているのが、相関係数となります。
「君の名は。」7/26発売予定 [テレビ / 映画]
ついに、予約販売が開始されました。
結局、映画観に行けなかったので、今から楽しみです！
結局、映画観に行けなかったので、今から楽しみです！
「君の名は。」Blu-rayスタンダード・エディション(早期購入特典:特製フィルムしおり付き)
- 出版社/メーカー: 東宝
- メディア: Blu-ray
映画『モアナと伝説の海』 [Disney / ディズニー]
7月5日に販売。
![モアナと伝説の海 MovieNEX [ブルーレイ+DVD+デジタルコピー(クラウド対応)+MovieNEXワールド] [Blu-ray]](https://images-fe.ssl-images-amazon.com/images/I/61nwiTiaixL._SL160_.jpg "モアナと伝説の海 MovieNEX [ブルーレイ+DVD+デジタルコピー(クラウド対応)+MovieNEXワールド] [Blu-ray]")
モアナと伝説の海 MovieNEX [ブルーレイ+DVD+デジタルコピー(クラウド対応)+MovieNEXワールド] [Blu-ray]
- 出版社/メーカー: ウォルト・ディズニー・ジャパン株式会社
- メディア: Blu-ray
かけっこ教室 2日目 [ファミリー]
先週に続き、今日が2回目。
毎日する宿題が出て、それなりに頑張っていました。
せっかくなので、前回よりもタイムが伸びてくれればと思っていたところ、
0.7秒ほど速くなっていました。
スタートダッシュがまだ甘いので、そこを鍛えれば、8秒台に行くんじゃないかなって感じました。
2回だけのトレーニングだと、限界があり、実際は、これからもトレーニングを続けていくことで、さらに速く走れるようになると思いますが、他の習い事との兼ね合いですね。。。
これ以上、増やすのは可愛そうだし。。。
毎日する宿題が出て、それなりに頑張っていました。
せっかくなので、前回よりもタイムが伸びてくれればと思っていたところ、
0.7秒ほど速くなっていました。
スタートダッシュがまだ甘いので、そこを鍛えれば、8秒台に行くんじゃないかなって感じました。
2回だけのトレーニングだと、限界があり、実際は、これからもトレーニングを続けていくことで、さらに速く走れるようになると思いますが、他の習い事との兼ね合いですね。。。
これ以上、増やすのは可愛そうだし。。。
TOEICテスト [英語 / ドイツ語]
公式のTOEICテストを受験しました。
SNSの会話のやり取りが問題に追加されいたりと、ずいぶんと傾向が変わった感じです。
テストの結果は、1か月くらいかかるみたいなので、6月中旬くらいですかね。。。
SNSの会話のやり取りが問題に追加されいたりと、ずいぶんと傾向が変わった感じです。
テストの結果は、1か月くらいかかるみたいなので、6月中旬くらいですかね。。。
第2期 電王戦 [将棋]
PONANZAが名人に連勝という結果に終わりました。
毎年、劇的に強くなるコンピュータですが、ここまで強くなっているとは想像できませんでした。
終盤、詰むか詰まないかの段階になると、コンピュータはメチャクチャ強いですが、
今回は、中盤の仕掛けの段階に置いても名人のヨミを超えていました。
電王戦は、今回が最後となりましたが、個人的には、羽生さんとコンピュータとの対局が観たかったです。
毎年、劇的に強くなるコンピュータですが、ここまで強くなっているとは想像できませんでした。
終盤、詰むか詰まないかの段階になると、コンピュータはメチャクチャ強いですが、
今回は、中盤の仕掛けの段階に置いても名人のヨミを超えていました。
電王戦は、今回が最後となりましたが、個人的には、羽生さんとコンピュータとの対局が観たかったです。
コンジョイント分析 [データサイエンス、統計モデル]
8章 コンジョイント分析について
よくやる間違い
属性水準（説明変数）として
・バッテリ時間
・保証期間
・色
・価格
を考えたときに、「価格」が入っているのは間違い
ダメな理由は、コンジョイント分析は、属性間の独立性を担保する必要があり、
価格は他の属性との相関があるため使えない。
現代マーケティング・リサーチ -- 市場を読み解くデータ分析
- 作者: 照井 伸彦
- 出版社/メーカー: 有斐閣
- 発売日: 2013/11/22
- メディア: 単行本（ソフトカバー）
よくやる間違い
属性水準（説明変数）として
・バッテリ時間
・保証期間
・色
・価格
を考えたときに、「価格」が入っているのは間違い
ダメな理由は、コンジョイント分析は、属性間の独立性を担保する必要があり、
価格は他の属性との相関があるため使えない。
Pokémon GO アドベンチャーウィーク [ゲーム]
オムナイト、カブト、ヨーギラスなどが大量に出ていましたね。
ヨーギラスは、普段より増えているとは思うのですが、
それでも出現する率は少なくて、結局1匹だけしか捕まえることができませんでした。。。
ヨーギラスは、普段より増えているとは思うのですが、
それでも出現する率は少なくて、結局1匹だけしか捕まえることができませんでした。。。
怪獣酒場 新橋蒸溜所 [グルメ / クッキング]
怪獣酒場 新橋蒸溜所
http://kaiju-sakaba.com/shimbashi/
川崎にあった怪獣酒場ですが、ついに！東京都内にも出来ました。
新橋ということで、ずいぶんと近いですね。
これは是非行かねば。
ただ、現在のところ予約はできないみたい。
http://kaiju-sakaba.com/shimbashi/
川崎にあった怪獣酒場ですが、ついに！東京都内にも出来ました。
新橋ということで、ずいぶんと近いですね。
これは是非行かねば。
ただ、現在のところ予約はできないみたい。
主人公はベリアルの息子！新TVシリーズ『ウルトラマンジード』 [ウルトラマン]
主人公はベリアルの息子！新TVシリーズ『ウルトラマンジード』http://m-78.jp/news/n-4686/
見るからに悪そうなウルトラマンジード。
それもそのはずで、ベリアルの息子らしいです。
べリアルって完全な悪なんですが、意外と人気があるので、この組み合わせもうなずけますね。
ウルトラマンゼロ：ウルトラセブンの息子
そして、べリアルと戦ったのが、ウルトラマンゼロ。
そもそもべリアルに息子がいたって設定が驚きですが、
どういう展開になるのか楽しみです。
2017年7月8日(土)～
ちなみに、ウルトラマンジードの主人公は、濱田龍臣。
竜馬伝に出ていた子ですが、16歳になったんですね。。。
見るからに悪そうなウルトラマンジード。
それもそのはずで、ベリアルの息子らしいです。
べリアルって完全な悪なんですが、意外と人気があるので、この組み合わせもうなずけますね。
ウルトラマンゼロ：ウルトラセブンの息子
そして、べリアルと戦ったのが、ウルトラマンゼロ。
そもそもべリアルに息子がいたって設定が驚きですが、
どういう展開になるのか楽しみです。
2017年7月8日(土)～
ちなみに、ウルトラマンジードの主人公は、濱田龍臣。
竜馬伝に出ていた子ですが、16歳になったんですね。。。
極楽湯千葉稲毛店 2回目 [ファミリー]
今日は珍しく大学がない土曜日。
極楽湯千葉稲毛店に行ってきました。
極楽湯 千葉稲毛店トップページ
https://www.gokurakuyu.ne.jp/tempo/chibainage/
行くまでが、それなりに遠いのですが、中は広くてとてもきれいな銭湯です。
中に入っているレストランで食事をしましたが、それなりに美味しく満足です。
会員カード作っていないけど、レストランでの割引などがあるから、次に行ったら作ってみます。
極楽湯千葉稲毛店に行ってきました。
極楽湯 千葉稲毛店トップページ
https://www.gokurakuyu.ne.jp/tempo/chibainage/
行くまでが、それなりに遠いのですが、中は広くてとてもきれいな銭湯です。
中に入っているレストランで食事をしましたが、それなりに美味しく満足です。
会員カード作っていないけど、レストランでの割引などがあるから、次に行ったら作ってみます。
小学校の運動会 2017 [ファミリー]
息子の運動会でした。
かけっこ教室で練習した成果がどうか。
スタートに立った時に、ちゃんとスタートの態勢を整えていました。
効果があったようです！
ただ、もう一度、「気をつけ」の号令がかかり、やり直し。
そして、再度態勢を作ろうとしたら、スタートの合図になりました。
ヨーイからのスタートの合図が早すぎ。。。
息子は大きく出遅れたのですが、そこから良い伸びをして、ゴールしました。
最後は、鼻の差くらいで負けてしまい2位でしたが、運がなかったですね。。。
来年に期待です。
でも、格段に速くなったと思います。
かけっこ教室で練習した成果がどうか。
スタートに立った時に、ちゃんとスタートの態勢を整えていました。
効果があったようです！
ただ、もう一度、「気をつけ」の号令がかかり、やり直し。
そして、再度態勢を作ろうとしたら、スタートの合図になりました。
ヨーイからのスタートの合図が早すぎ。。。
息子は大きく出遅れたのですが、そこから良い伸びをして、ゴールしました。
最後は、鼻の差くらいで負けてしまい2位でしたが、運がなかったですね。。。
来年に期待です。
でも、格段に速くなったと思います。
マルニウメバチ ワインバー [グルメ / クッキング]
マルニウメバチ
https://tabelog.com/tokyo/A1301/A130101/13122481/
店内写真撮影不可のお店。
ただ、一品一品がとても美味しくて、ワインにぴったり合います。
値段はそれなりに高いと思うのですが、、、
（ご馳走様でした！）
https://tabelog.com/tokyo/A1301/A130101/13122481/
店内写真撮影不可のお店。
ただ、一品一品がとても美味しくて、ワインにぴったり合います。
値段はそれなりに高いと思うのですが、、、
（ご馳走様でした！）
Rを使って重回帰分析を色々な方法で解く その1 [データサイエンス、統計モデル]
まずは、パッケージ lm() を使って解く方法と行列を使って解く方法から。
# irisのデータを使用
head(iris)
# Species == "setosa" 50レコードを使用
dat <- subset(iris, Species == "setosa", c(Sepal.Length, Sepal.Width, Petal.Length, Petal.Width ))
# 相関行列
round(cor(dat), 3)
Sepal.Length Sepal.Width Petal.Length Petal.Width
Sepal.Length 1.000 0.743 0.267 0.278
Sepal.Width 0.743 1.000 0.178 0.233
Petal.Length 0.267 0.178 1.000 0.332
Petal.Width 0.278 0.233 0.332 1.000
### lm を使った重回帰式 ###
# Length = β1 * Width + β2 * Length + β3 *Width + β0
fit.lm <- lm(Sepal.Length ~ . , data = dat)
summary(fit.lm)
Call:
lm(formula = Sepal.Length ~ ., data = dat)
Residuals:
Min 1Q Median 3Q Max
-0.40662 -0.17721 0.01222 0.13388 0.49693
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 2.35189 0.39287 5.986 3.03e-07 ***
Sepal.Width 0.65483 0.09245 7.083 6.83e-09 ***
Petal.Length 0.23756 0.20802 1.142 0.259
Petal.Width 0.25213 0.34686 0.727 0.471
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 0.2371 on 46 degrees of freedom
Multiple R-squared: 0.5751, Adjusted R-squared: 0.5474
F-statistic: 20.76 on 3 and 46 DF, p-value: 1.192e-08
### 行列を使って解く方法 ###
solve(t(x_ols) %*% x_ols) %*% t(x_ols) %*% y_ols
[,1]
[1,] 2.3518898
[2,] 0.6548350
[3,] 0.2375602
[4,] 0.2521257
# irisのデータを使用
head(iris)
# Species == "setosa" 50レコードを使用
dat <- subset(iris, Species == "setosa", c(Sepal.Length, Sepal.Width, Petal.Length, Petal.Width ))
# 相関行列
round(cor(dat), 3)
Sepal.Length Sepal.Width Petal.Length Petal.Width
Sepal.Length 1.000 0.743 0.267 0.278
Sepal.Width 0.743 1.000 0.178 0.233
Petal.Length 0.267 0.178 1.000 0.332
Petal.Width 0.278 0.233 0.332 1.000
### lm を使った重回帰式 ###
# Length = β1 * Width + β2 * Length + β3 *Width + β0
fit.lm <- lm(Sepal.Length ~ . , data = dat)
summary(fit.lm)
Call:
lm(formula = Sepal.Length ~ ., data = dat)
Residuals:
Min 1Q Median 3Q Max
-0.40662 -0.17721 0.01222 0.13388 0.49693
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 2.35189 0.39287 5.986 3.03e-07 ***
Sepal.Width 0.65483 0.09245 7.083 6.83e-09 ***
Petal.Length 0.23756 0.20802 1.142 0.259
Petal.Width 0.25213 0.34686 0.727 0.471
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 0.2371 on 46 degrees of freedom
Multiple R-squared: 0.5751, Adjusted R-squared: 0.5474
F-statistic: 20.76 on 3 and 46 DF, p-value: 1.192e-08
### 行列を使って解く方法 ###
solve(t(x_ols) %*% x_ols) %*% t(x_ols) %*% y_ols
[,1]
[1,] 2.3518898
[2,] 0.6548350
[3,] 0.2375602
[4,] 0.2521257
Rを使って重回帰分析を色々な方法で解く その2 [データサイエンス、統計モデル]
Rを使って重回帰分析を色々な方法で解く その1
http://skellington.blog.so-net.ne.jp/2017-05-30
↑
こちらの続き。
汎用最適化関数 optim() を使って解く方法。
### optim を使ったパラメータ推定 ###
b0 <- c(0, 0, 0, 0) # 推定するパラメータの初期化
num_row <- nrow(dat) # レコード数
y_ols = dat$Sepal.Length
x_ols = cbind(1, dat$Sepal.Width, dat$Petal.Length, dat$Petal.Width)
ols <- function (y_ols, x_ols) {
min_rss <- function (param) {
sum((y_ols - x_ols %*% param)^2)
}
k <- ncol(x_ols)
return(optim(par = rep(0, k), fn = min_rss, method = "BFGS", hessian = TRUE))
}
result <- ols(y_ols, x_ols)
result
$par
[1] 2.3518922 0.6548347 0.2375590 0.2521268
$value # 対数尤度値
[1] 2.586648
$counts
function gradient
24 8
$convergence # 収束判定
[1] 0
$message
NULL
$hessian
[,1] [,2] [,3] [,4]
[1,] 100.0 342.80 146.20 24.60
[2,] 342.8 1189.20 502.32 85.24
[3,] 146.2 502.32 216.70 36.56
[4,] 24.6 85.24 36.56 7.14
http://skellington.blog.so-net.ne.jp/2017-05-30
↑
こちらの続き。
汎用最適化関数 optim() を使って解く方法。
### optim を使ったパラメータ推定 ###
b0 <- c(0, 0, 0, 0) # 推定するパラメータの初期化
num_row <- nrow(dat) # レコード数
y_ols = dat$Sepal.Length
x_ols = cbind(1, dat$Sepal.Width, dat$Petal.Length, dat$Petal.Width)
ols <- function (y_ols, x_ols) {
min_rss <- function (param) {
sum((y_ols - x_ols %*% param)^2)
}
k <- ncol(x_ols)
return(optim(par = rep(0, k), fn = min_rss, method = "BFGS", hessian = TRUE))
}
result <- ols(y_ols, x_ols)
result
$par
[1] 2.3518922 0.6548347 0.2375590 0.2521268
$value # 対数尤度値
[1] 2.586648
$counts
function gradient
24 8
$convergence # 収束判定
[1] 0
$message
NULL
$hessian
[,1] [,2] [,3] [,4]
[1,] 100.0 342.80 146.20 24.60
[2,] 342.8 1189.20 502.32 85.24
[3,] 146.2 502.32 216.70 36.56
[4,] 24.6 85.24 36.56 7.14
