Tree and Forest

今回の分析モデルは 決定木 と ランダム森 です。これらも回帰モデルの一つですが，決定木の特徴は非線形であること(階段関数)や従属変数が間隔尺度(連続)であっても名義尺度(離散)であっても構わないこと，ランダム森の特徴はそれに加えて 機械学習 と呼ばれる手法の一つであるところです。

決定木desicion treeとは

決定木とは，回帰分析のように直線をあてはめるのではなく，分岐を繰り返して対象にたどり着くモデルです。質問にYes/Noで回答して行って，最終的にどこかに分類される占いをイメージしてもらえれば分かりやすいでしょう。ただし，占いと違ってこちらはデータに基づき，もっとも有意義な分岐を探すというところが違いますが。

従属変数が連続的な数値であるときは，回帰木regression tree とも言われます。一方，連続変数が離散的な数値であるときは，分類木classification tree とも言われます。

回帰木の例

Rに含まれているcarsというデータを例にしてみましょう。スピードと車が停車するまでの距離のデータセットです。

data(cars)
summary(cars)

##      speed           dist       
##  Min.   : 4.0   Min.   :  2.00  
##  1st Qu.:12.0   1st Qu.: 26.00  
##  Median :15.0   Median : 36.00  
##  Mean   :15.4   Mean   : 42.98  
##  3rd Qu.:19.0   3rd Qu.: 56.00  
##  Max.   :25.0   Max.   :120.00

plot(cars$speed,cars$dist)

停車するまでの距離とスピードの関係について，回帰分析をすると次のようになります。

result.lm <- lm(dist~speed,data=cars)
summary(result.lm)

## 
## Call:
## lm(formula = dist ~ speed, data = cars)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -29.069  -9.525  -2.272   9.215  43.201 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -17.5791     6.7584  -2.601   0.0123 *  
## speed         3.9324     0.4155   9.464 1.49e-12 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 15.38 on 48 degrees of freedom
## Multiple R-squared:  0.6511, Adjusted R-squared:  0.6438 
## F-statistic: 89.57 on 1 and 48 DF,  p-value: 1.49e-12

result.glm <- glm(dist~speed,data=cars,family="Gamma")
summary(result.glm)

## 
## Call:
## glm(formula = dist ~ speed, family = "Gamma", data = cars)
## 
## Deviance Residuals: 
##      Min        1Q    Median        3Q       Max  
## -1.64429  -0.31120  -0.03906   0.14006   1.09571  
## 
## Coefficients:
##               Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  0.0614349  0.0058043  10.584 3.82e-14 ***
## speed       -0.0021315  0.0002772  -7.689 6.52e-10 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for Gamma family taken to be 0.1952369)
## 
##     Null deviance: 22.483  on 49  degrees of freedom
## Residual deviance: 10.954  on 48  degrees of freedom
## AIC: 427.39
## 
## Number of Fisher Scoring iterations: 5

AIC(result.lm)

## [1] 419.1569

AIC(result.glm) #Gaussianモデルで十分

## [1] 427.3911

plot(cars$speed,cars$dist)
abline(result.lm,col="red")

さて，これを回帰木でやるとどうなるでしょうか。treeパッケージをダウンロード・装備(library)して，実行してみましょう。

library(tree)
result.tree <- tree(dist~speed,data=cars)
print(result.tree)

## node), split, n, deviance, yval
##       * denotes terminal node
## 
## 1) root 50 32540.0 42.98  
##   2) speed < 17.5 31  8307.0 29.32  
##     4) speed < 12.5 15  1176.0 18.20  
##       8) speed < 9.5 6   277.3 10.67 *
##       9) speed > 9.5 9   331.6 23.22 *
##     5) speed > 12.5 16  3535.0 39.75 *
##   3) speed > 17.5 19  9016.0 65.26  
##     6) speed < 23.5 14  2847.0 55.71 *
##     7) speed > 23.5 5  1318.0 92.00 *

これは木の形で見た方が分かりやすいところです。

plot(result.tree,type="u")
text(result.tree)

回帰線の形で表現すると，階段関数になっていることが分かります。

plot(cars$speed,cars$dist)
partition.tree(result.tree,add=T,col="red")

mvpartパッケージでも同様のことができます。

library(mvpart)
result.rpart <- rpart(dist~speed,data=cars,method="anova")
result.rpart

## n= 50 
## 
## node), split, n, deviance, yval
##       * denotes terminal node
## 
##  1) root 50 32538.9800 42.98000  
##    2) speed< 17.5 31  8306.7740 29.32258  
##      4) speed< 12.5 15  1176.4000 18.20000  
##        8) speed< 9.5 6   277.3333 10.66667 *
##        9) speed>=9.5 9   331.5556 23.22222 *
##      5) speed>=12.5 16  3535.0000 39.75000 *
##    3) speed>=17.5 19  9015.6840 65.26316  
##      6) speed< 23.5 14  2846.8570 55.71429  
##       12) speed>=18.5 10  1323.6000 52.20000 *
##       13) speed< 18.5 4  1091.0000 64.50000 *
##      7) speed>=23.5 5  1318.0000 92.00000 *

plot(result.rpart)
text(result.rpart)

tree関数より少し複雑な木ができました。木はどこまでも複雑にしていけますが，複雑さと分かりやすさはトレードオフになりますので，自分で調節することができます。一つの目安は，複雑さ指標cpをによるもので，データから次のように見積もることができます。

plotcp(result.rpart)

プロットされたオレンジの値が，提案される基準です。そこで剪定(prune)します。

result.prune <- prune(result.rpart,cp=0.044)
plot(result.prune)
text(result.prune)

分類木の例

分類木は，分析対象が名義尺度水準の場合に使われます。タイタニック号のデータ(Rに入っている！)でやってみましょう。ちょっとデータを分かりやすく加工します。

data(Titanic)
Titanic

## , , Age = Child, Survived = No
## 
##       Sex
## Class  Male Female
##   1st     0      0
##   2nd     0      0
##   3rd    35     17
##   Crew    0      0
## 
## , , Age = Adult, Survived = No
## 
##       Sex
## Class  Male Female
##   1st   118      4
##   2nd   154     13
##   3rd   387     89
##   Crew  670      3
## 
## , , Age = Child, Survived = Yes
## 
##       Sex
## Class  Male Female
##   1st     5      1
##   2nd    11     13
##   3rd    13     14
##   Crew    0      0
## 
## , , Age = Adult, Survived = Yes
## 
##       Sex
## Class  Male Female
##   1st    57    140
##   2nd    14     80
##   3rd    75     76
##   Crew  192     20

tmp <- data.frame(Titanic)
Titanic.df <- data.frame(
                          Class = rep(tmp$Class, tmp$Freq),
                          Sex = rep(tmp$Sex, tmp$Freq),
                          Age = rep(tmp$Age, tmp$Freq),
                          Survived = rep(tmp$Survived, tmp$Freq)
              )
summary(Titanic.df)

##   Class         Sex          Age       Survived  
##  1st :325   Male  :1731   Child: 109   No :1490  
##  2nd :285   Female: 470   Adult:2092   Yes: 711  
##  3rd :706                                        
##  Crew:885

head(Titanic.df)

##   Class  Sex   Age Survived
## 1   3rd Male Child       No
## 2   3rd Male Child       No
## 3   3rd Male Child       No
## 4   3rd Male Child       No
## 5   3rd Male Child       No
## 6   3rd Male Child       No

で，rpartを。

result.rpart2 <- rpart(Survived ~ .,data=Titanic.df,method="class")
result.rpart2

## n= 2201 
## 
## node), split, n, loss, yval, (yprob)
##       * denotes terminal node
## 
##  1) root 2201 711 No (0.6769650 0.3230350)  
##    2) Sex=Male 1731 367 No (0.7879838 0.2120162)  
##      4) Age=Adult 1667 338 No (0.7972406 0.2027594) *
##      5) Age=Child 64  29 No (0.5468750 0.4531250)  
##       10) Class=3rd 48  13 No (0.7291667 0.2708333) *
##       11) Class=1st,2nd 16   0 Yes (0.0000000 1.0000000) *
##    3) Sex=Female 470 126 Yes (0.2680851 0.7319149)  
##      6) Class=3rd 196  90 No (0.5408163 0.4591837) *
##      7) Class=1st,2nd,Crew 274  20 Yes (0.0729927 0.9270073) *

plot(result.rpart2)
text(result.rpart2)

ランダムフォレスト

木を一本ではなく，何本も使って予測精度を上げよう，というのがランダムフォレスト(ランダム森)という考え方です。

library(randomForest)

## randomForest 4.6-10
## Type rfNews() to see new features/changes/bug fixes.

result.rf <- randomForest(Survived ~.,data=Titanic.df,ntree=1000)
result.rf

## 
## Call:
##  randomForest(formula = Survived ~ ., data = Titanic.df, ntree = 1000) 
##                Type of random forest: classification
##                      Number of trees: 1000
## No. of variables tried at each split: 1
## 
##         OOB estimate of  error rate: 21.26%
## Confusion matrix:
##       No Yes class.error
## No  1463  27  0.01812081
## Yes  441 270  0.62025316

多くの予測変数があって，変数の重要度を知りたいときはimportanceという指標で見ることができます。またそれをプロットすることもできます。

importance(result.rf)

##       MeanDecreaseGini
## Class        50.142343
## Sex         133.528058
## Age           8.908746

varImpPlot(result.rf)

Tree and Forest

Kosugitti

2014年12月12日

決定木desicion treeとは

回帰木の例

分類木の例

ランダムフォレスト

まとめなど

本日の課題