統計ER

統計ソフトRの使い方を中心に、統計解析方法の解説をするブログ。ありそうでなかなか見つからないサンプルサイズ計算などニッチな方法について紹介しています。

統計ソフトRのISLRパッケージWeeklyデータの分析例

ブログランキングに参加しています。
まずはぽちぽちっとお願いします。
↓↓↓↓↓↓↓↓↓↓↓↓↓↓↓↓↓↓
にほんブログ村 科学ブログ 数学へ

機械学習で、よりよく推測できるモデルを選ぶ。

統計ソフトRのISLRパッケージのWeeklyデータで基礎的な機械学習を行ってみた。

 

 

データの準備

ISLRパッケージのWeeklyデータは、S&P500指数の週当たりのリターンのデータ。

www.google.co.jp

9つの変数、1089行のデータから成り立っている。

始めに一回だけインストールする。

install.packages("ISLR")

 

使うときには呼び出す。

library(ISLR)

 

Weeklyデータを確認してみる。

> str(Weekly)
'data.frame':   1089 obs. of  9 variables:
 $ Year     : num  1990 1990 1990 1990 1990 1990 1990 1990 1990 1990 ...
 $ Lag1     : num  0.816 -0.27 -2.576 3.514 0.712 ...
 $ Lag2     : num  1.572 0.816 -0.27 -2.576 3.514 ...
 $ Lag3     : num  -3.936 1.572 0.816 -0.27 -2.576 ...
 $ Lag4     : num  -0.229 -3.936 1.572 0.816 -0.27 ...
 $ Lag5     : num  -3.484 -0.229 -3.936 1.572 0.816 ...
 $ Volume   : num  0.155 0.149 0.16 0.162 0.154 ...
 $ Today    : num  -0.27 -2.576 3.514 0.712 1.178 ...
 $ Direction: Factor w/ 2 levels "Down","Up": 1 1 2 2 2 1 2 2 2 1 ...

 

ロジスティック回帰モデル

Lag1からLag5、Volumeの6つの変数でDirectionを予測するロジスティック回帰モデルを作成する。

Directionはポジティブなリターンかネガティブなリターンか。

Lag2だけ統計学的有意に関連するという結果。Lag2は2週間前からみたパーセンテージリターン。ポジティブリターンに関連している。

> res1 <- glm(Direction ~ Lag1 + Lag2 + Lag3 + Lag4 + Lag5 + Volume, family=binomial, data=Weekly)
> 
> summary(res1)

Call:
glm(formula = Direction ~ Lag1 + Lag2 + Lag3 + Lag4 + Lag5 + 
    Volume, family = binomial, data = Weekly)

Deviance Residuals: 
    Min       1Q   Median       3Q      Max  
-1.6949  -1.2565   0.9913   1.0849   1.4579  

Coefficients:
            Estimate Std. Error z value Pr(>|z|)   
(Intercept)  0.26686    0.08593   3.106   0.0019 **
Lag1        -0.04127    0.02641  -1.563   0.1181   
Lag2         0.05844    0.02686   2.175   0.0296 * 
Lag3        -0.01606    0.02666  -0.602   0.5469   
Lag4        -0.02779    0.02646  -1.050   0.2937   
Lag5        -0.01447    0.02638  -0.549   0.5833   
Volume      -0.02274    0.03690  -0.616   0.5377   
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for binomial family taken to be 1)

    Null deviance: 1496.2  on 1088  degrees of freedom
Residual deviance: 1486.4  on 1082  degrees of freedom
AIC: 1500.4

Number of Fisher Scoring iterations: 4

 

学習データとテストデータに分割して正答率を見る

機械学習手法の基本は、学習データで作ったモデルが、テストデータでどのくらい予測に使えるかを確認するもの。

ばっちり予測できるなら、作ったモデルは使える可能性が高い。

学習データでは予測できるが、テストデータになったら途端に予測できなくなるようなら、汎用性は低く使えるモデルではない。

1990 年から2008 年までを学習データ、残り(2009 と2010 年)をテストデータとする。

Direction をLag2だけで推測するロジスティック回帰モデルで分析してみる。正答率は混合行列で計算する。

まず、training用とtest用のData frameに分ける。

Weekly.tr <- subset (Weekly, Weekly$Year < 2009)

Weekly.te <- subset (Weekly, Weekly$Year >= 2009)

 

training用のデータでlogistic regression modelを作成する。

> res.tr <- glm(Direction ~ Lag2, family=binomial, data=Weekly.tr)
> 
> summary(res.tr)

Call:
glm(formula = Direction ~ Lag2, family = binomial, data = Weekly.tr)

Deviance Residuals: 
   Min      1Q  Median      3Q     Max  
-1.536  -1.264   1.021   1.091   1.368  

Coefficients:
            Estimate Std. Error z value Pr(>|z|)   
(Intercept)  0.20326    0.06428   3.162  0.00157 **
Lag2         0.05810    0.02870   2.024  0.04298 * 
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for binomial family taken to be 1)

    Null deviance: 1354.7  on 984  degrees of freedom
Residual deviance: 1350.5  on 983  degrees of freedom
AIC: 1354.5

Number of Fisher Scoring iterations: 4

 

test用データを用いていま作ったモデルから予測値を算出する。

> pre.te <- predict(res.tr, newdata=Weekly.te, type="response")
> 
> summary(pre.te)
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
 0.4488  0.5331  0.5573  0.5555  0.5826  0.6954 

 

混同行列confusion matrixを作成し、正答率 accuracyを計算する。62.5%と計算された。

> (tab.te <- table(Weekly.te$Direction, pre.te>0.5))
      
       FALSE TRUE
  Down     9   34
  Up       5   56
> 
> sum(diag(prop.table(tab.te)))*100
[1] 62.5

 

線形判別分析で行う方法

ロジスティック回帰と似た方法で線形判別分析がある。線形判別分析で同様に予測させるとどうなるか?

線形判別分析の関数lda()はMASSパッケージにあるので、MASSパッケージを呼び出してから使用する。

library(MASS)

res.lda.tr <- lda(Direction ~ Lag2, data=Weekly.tr)

pre.lda.te <- predict(res.lda.tr, newdata=Weekly.te)

(tab.lda.te <- table(Weekly.te$Direction, pre.lda.te$class))

sum(diag(prop.table(tab.lda.te)))*100

 

正答率は、ロジスティック回帰と同じ62.5%となった。

> (tab.lda.te <- table(Weekly.te$Direction, pre.lda.te$class))
      
       Down Up
  Down    9 34
  Up      5 56
> 
> sum(diag(prop.table(tab.lda.te)))*100
[1] 62.5

 

二次判別分析で行う方法

二次判別分析は、線形判別分析で必要だった等分散性が要らないという利点がある。

二次判別分析で正答率を見てみる。MASSパッケージのqda()を使う。

res.qda.tr <- qda(Direction ~ Lag2, data=Weekly.tr)

pre.qda.te <- predict(res.qda.tr, newdata=Weekly.te)

(tab.qda.te <- table(Weekly.te$Direction, pre.qda.te$class))

sum(diag(prop.table(tab.qda.te)))*100 

 

Downに一件も正しく判別できず、正答率は58.7%に下がった。

> (tab.qda.te <- table(Weekly.te$Direction, pre.qda.te$class))
      
       Down Up
  Down    0 43
  Up      0 61
> 
> sum(diag(prop.table(tab.qda.te)))*100
[1] 58.65385

 

参考:

www.statmethods.net

 

k近傍法で行う方法

k近傍法は、判別したい新しいデータの近傍にあるk個データから、どのグループに割り振ったら(判別したら)いいかを考える方法。

kを1個、2個、3個と増やして、判別が正解する率を比較する。

classパッケージを用いる。classパッケージはインストール済みだ。

library(class)

 

タイのデータがあった場合ランダムに分けるため、答えが毎回変わる。

答えが毎回変わるのを防止するために、set.seed()を設定する。

k=1の場合

set.seed(20180902)

res.knn.te1 <- knn(train=data.frame(Weekly.tr$Lag2), test=data.frame(Weekly.te$Lag2), cl=Weekly.tr$Direction, k=1)

summary(res.knn.te1)

(tab.knn.te1 <- table(Weekly.te$Direction,res.knn.te1))

sum(diag(prop.table(tab.knn.te1)))*100

 

正答率は51.0%であった。

> (tab.knn.te1 <- table(Weekly.te$Direction,res.knn.te1))
      res.knn.te1
       Down Up
  Down   21 22
  Up     29 32
> 
> sum(diag(prop.table(tab.knn.te1)))*100
[1] 50.96154

 

k=2の場合

set.seed(20180902)

res.knn.te2 <- knn(train=data.frame(Weekly.tr$Lag2), test=data.frame(Weekly.te$Lag2), cl=Weekly.tr$Direction, k=2)

summary(res.knn.te2)

(tab.knn.te2 <- table(Weekly.te$Direction,res.knn.te2))

sum(diag(prop.table(tab.knn.te2)))*100

 

正答率は53.8%に上昇した。

> (tab.knn.te2 <- table(Weekly.te$Direction,res.knn.te2))
      res.knn.te2
       Down Up
  Down   22 21
  Up     27 34
> 
> sum(diag(prop.table(tab.knn.te2)))*100
[1] 53.84615

 

k=3の場合

set.seed(20180902)

res.knn.te3 <- knn(train=data.frame(Weekly.tr$Lag2), test=data.frame(Weekly.te$Lag2), cl=Weekly.tr$Direction, k=3)

summary(res.knn.te3)

(tab.knn.te3 <- table(Weekly.te$Direction,res.knn.te3))

sum(diag(prop.table(tab.knn.te3)))*100

 

正答率55.8%になった。

> (tab.knn.te3 <- table(Weekly.te$Direction,res.knn.te3))
      res.knn.te3
       Down Up
  Down   16 27
  Up     19 42
> 
> sum(diag(prop.table(tab.knn.te3)))*100
[1] 55.76923

 

参考:

stackoverflow.com

まとめ

統計ソフトRのISLRパッケージWeeklyデータを用いて、学習セットで作成したモデルのパフォーマンスをテストデータで検証する機械学習を行った。

ロジスティック回帰、線形判別分析、二次判別分析、k近傍法を比較した。

結果として、ロジスティック回帰と線形判別分析が同率一位(正答率62.5%)であった。