My Notes

統計学とかR(R言語)とかPython3の覚え書きとか走り書きとか。 座右の銘にしたい: All work and no play makes Jack a dull boy.

『マンガでわかる統計学 回帰分析編』第2章 回帰分析 をR(R言語)で

Rコード

# 『マンガでわかる統計学 回帰分析編』第2章 回帰分析

df <- data.frame(最高気温 = c(29, 28, 34, 31, 25,
                            29, 32, 31, 24, 33,
                            25, 31, 26, 30),
                アイスティーの注文数 = c(77, 62, 93, 84, 59,
                                    64, 80, 75, 58, 91,
                                    51, 73, 65, 84))

df


最高気温 <- df$最高気温
最高気温

アイスティーの注文数 <- df$アイスティーの注文数
アイスティーの注文数




summary(df)




# 文字化けするなら
par(family = "Osaka")
# 散布図
plot(df)


# 相関係数
cor(最高気温, アイスティーの注文数)
round(cor(最高気温, アイスティーの注文数), 4)
# 相関係数の検定
cor.test(最高気温, アイスティーの注文数)


# 回帰式
lm1 <- lm(アイスティーの注文数 ~ 最高気温, data = df)
lm1
summary(lm1)
# 回帰直線
abline(lm1, lwd = 2)


# df、予測値(推定値)、残差
予測値 <- predict(lm1)
残差 <- residuals(lm1)
data.frame(df, 予測値, 残差)


# 回帰診断図
par(mfrow = c(2, 2), oma = c(2, 2, 2, 1), mar = c(4, 4, 3, 2))
plot(lm1)




# 予測してみる
# 回帰方程式 y = ax + b
# 回帰係数a
a <- coef(lm1)[[2]]
a
# 切片b
b <- coef(lm1)[[1]]
b


# 確認 29度の場合
y1 <- a*29 + b
y1
# df、予測値(推定値)、残差
# 予測値 <- predict(lm1)
# 残差 <- residuals(lm1)
# data.frame(df, 予測値, 残差)
#
#    最高気温 アイスティーの注文数   予測値      残差
# 1        29                   77 72.03744  4.962555
#
# 計算は合ってる


max(最高気温)
min(最高気温)
# 24度から34度までなので、この範囲内なら内挿
# 27度から予測
y2 <- a*27 + b
y2
round(y2) # 27度の場合は、アイスティーの注文数は65杯ほどと予測できる


# 外挿してみる(注意事項は、p.96)
y3 <- a*35 + b
y3
round(y3)

R Console

> # 『マンガでわかる統計学 回帰分析編』第2章 回帰分析
> 
> df <- data.frame(最高気温 = c(29, 28, 34, 31, 25,
+                             29, 32, 31, 24, 33,
+                             25, 31, 26, 30),
+                 アイスティーの注文数 = c(77, 62, 93, 84, 59,
+                                     64, 80, 75, 58, 91,
+                                     51, 73, 65, 84))
> 
> df
   最高気温 アイスティーの注文数
1        29                   77
2        28                   62
3        34                   93
4        31                   84
5        25                   59
6        29                   64
7        32                   80
8        31                   75
9        24                   58
10       33                   91
11       25                   51
12       31                   73
13       26                   65
14       30                   84
> 
> 
> 最高気温 <- df$最高気温
> 最高気温
 [1] 29 28 34 31 25 29 32 31 24 33 25 31 26 30
> 
> アイスティーの注文数 <- df$アイスティーの注文数
> アイスティーの注文数
 [1] 77 62 93 84 59 64 80 75 58 91 51 73 65 84
> 
> 
> 
> 
> summary(df)
    最高気温     アイスティーの注文数
 Min.   :24.00   Min.   :51.00       
 1st Qu.:26.50   1st Qu.:62.50       
 Median :29.50   Median :74.00       
 Mean   :29.14   Mean   :72.57       
 3rd Qu.:31.00   3rd Qu.:83.00       
 Max.   :34.00   Max.   :93.00       
> 
> 
> 
> 
> # 文字化けするなら
> par(family = "Osaka")
> # 散布図
> plot(df)
> 
> 
> # 相関係数
> cor(最高気温, アイスティーの注文数)
[1] 0.906923
> round(cor(最高気温, アイスティーの注文数), 4)
[1] 0.9069
> # 相関係数の検定
> cor.test(最高気温, アイスティーの注文数)

    Pearson's product-moment correlation

data:  最高気温 and アイスティーの注文数
t = 7.4572, df = 12, p-value = 7.661e-06
alternative hypothesis: true correlation is not equal to 0
95 percent confidence interval:
 0.7254041 0.9705020
sample estimates:
     cor 
0.906923 

> 
> 
> # 回帰式
> lm1 <- lm(アイスティーの注文数 ~ 最高気温, data = df)
> lm1

Call:
lm(formula = アイスティーの注文数 ~ 最高気温, data = df)

Coefficients:
(Intercept)     最高気温  
    -36.361        3.738  

> summary(lm1)

Call:
lm(formula = アイスティーの注文数 ~ 最高気温, data = df)

Residuals:
   Min     1Q Median     3Q    Max 
-8.037 -5.693  2.094  4.409  8.225 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept) -36.3612    14.6873  -2.476   0.0292 *  
最高気温      3.7379     0.5012   7.457 7.66e-06 ***
---
Signif. codes:  0***0.001**0.01*0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 5.709 on 12 degrees of freedom
Multiple R-squared:  0.8225, Adjusted R-squared:  0.8077 
F-statistic: 55.61 on 1 and 12 DF,  p-value: 7.661e-06

> # 回帰直線
> abline(lm1, lwd = 2)
> 
> 
> # df、予測値(推定値)、残差
> 予測値 <- predict(lm1)
> 残差 <- residuals(lm1)
> data.frame(df, 予測値, 残差)
   最高気温 アイスティーの注文数   予測値      残差
1        29                   77 72.03744  4.962555
2        28                   62 68.29956 -6.299559
3        34                   93 90.72687  2.273128
4        31                   84 79.51322  4.486784
5        25                   59 57.08590  1.914097
6        29                   64 72.03744 -8.037445
7        32                   80 83.25110 -3.251101
8        31                   75 79.51322 -4.513216
9        24                   58 53.34802  4.651982
10       33                   91 86.98899  4.011013
11       25                   51 57.08590 -6.085903
12       31                   73 79.51322 -6.513216
13       26                   65 60.82379  4.176211
14       30                   84 75.77533  8.224670
> 
> 
> # 回帰診断図
> par(mfrow = c(2, 2), oma = c(2, 2, 2, 1), mar = c(4, 4, 3, 2))
> plot(lm1)
> 
> 
> 
> 
> # 予測してみる
> # 回帰方程式 y = ax + b
> # 回帰係数a
> a <- coef(lm1)[[2]]
> a
[1] 3.737885
> # 切片b
> b <- coef(lm1)[[1]]
> b
[1] -36.36123
> 
> 
> # 確認 29度の場合
> y1 <- a*29 + b
> y1
[1] 72.03744
> # df、予測値(推定値)、残差
> # 予測値 <- predict(lm1)
> # 残差 <- residuals(lm1)
> # data.frame(df, 予測値, 残差)
> #
> #    最高気温 アイスティーの注文数   予測値      残差
> # 1        29                   77 72.03744  4.962555
> #
> # 計算は合ってる
> 
> 
> max(最高気温)
[1] 34
> min(最高気温)
[1] 24
> # 24度から34度までなので、この範囲内なら内挿
> # 27度から予測
> y2 <- a*27 + b
> y2
[1] 64.56167
> round(y2) # 27度の場合は、アイスティーの注文数は65杯ほどと予測できる
[1] 65
> 
> 
> # 外挿してみる(注意事項は、p.96)
> y3 <- a*35 + b
> y3
[1] 94.46476
> round(y3)
[1] 94

散布図のスクリーンショット

f:id:my_notes:20170616062320p:plain

回帰直線のスクリーンショット

f:id:my_notes:20170616062351p:plain

回帰診断図のスクリーンショット

f:id:my_notes:20170616062710p:plain

参考文献

マンガでわかる統計学 回帰分析編

マンガでわかる統計学 回帰分析編

参考にしたRコード

Rによるデータサイエンス データ解析の基礎から最新手法まで

Rによるデータサイエンス データ解析の基礎から最新手法まで

新版ならこち

Rによるデータサイエンス データ解析の基礎から最新手法まで

Rによるデータサイエンス データ解析の基礎から最新手法まで