랜덤포리스트(Random Forest) 총정리(배깅 포함)

앙상블(Ensemble) 입문

드디어 앙상블(Ensemble) 모델이다. 앙상블이란 말 그대로 각기 다른 모델을 조합해 예측력을 향상시키는 것을 말한다. 책, Predictive Ananlytics에는 다양한 예측 사업들이 소개되는데 가장 높은 예측력을 가지는 모델은 단연 앙상블이었다. 어떤 모델이든 장점과 단점을 가지게 되는데 앙상블은 각각의 장점은 취하면서 단점은 보완할 수 있다. 때문에 앙상블은 최적의 단일 모델을 뛰어 넘는 새로운 예측 모델을 제시한다.

높은 예측력은 분명 훌륭한 모델의 가장 중요한 조건중 하나이다. 그러나 이 기준 역시 분명히 단점이 있다. 예를 들어, 1. 앙상블과 같이 예측력이 높은 모델은 해석에 어려움을 가진다. 또한 2. 다른 단일 모델에 비해 예측 시간이 많이 걸린다. 존경하는 Data Smart의 저자 John Foreman은 예측력을 최대로 하는 여러 모델의 조합을 찾는 시간에 차라리 한가지 모델에서 그 현상을 최적으로 설명할 수 있는 설명변수를 찾을 것을 조언한다. 그의 의견에 전적으로 동의한다. 그러나 앙상블 모델의 대표인 랜덤포리스트(Random Forests)와 부스팅(Boosting)은 예측 분석에 있어 반드시 익혀야할 최소한이라 생각한다. 각각은 회귀트리 또는 분류트리와의 앙상블을 통해 최적 모델을 구현한다. 회귀트리와 분류트리는 이전글, 회귀트리(regression tree)와 분류트리(classification tree) 예측을 참조바란다. 랜덤포리스트를 이해하기 위해서는 먼저 배깅을 이해해야한다. 랜덤포리스트는 일종의 배깅이며 배깅 방법의 대부분은 랜덤포리스트와 같기 때문이다.

배깅(Bagging) 알고리즘

배깅은 Bootstrap Agregating의 준말이다. 직역하면 “주머니 통합”이다. 이것은 배깅 알고리즘에 대한 비유로 배깅이 작동하는 알고리즘을 이해하면 비유가 이해가 된다. 배깅을 구현하는 순서는 아래와 같다.(상기 그림 참조)
1. D개의 전체데이터가 있다.
2. D와 같은 갯수의 데이터를 전체데이터에서 복원추출하고 이 집단(주머니)을 D1이라 하자.
3. D1 데이터를 통해 트리를 만든다.
4. 위 2번, 3번과 같은 방법으로 m개의 또 다른 Dm과 각각의 트리를 만든다. 각각의 트리들은 크기도 크고 prune되지 않은 것이다.
5. m개의 트리의 평균을 통해 예측을 한다. 회귀트리(양적반응변수)가 아닌 분류트리(질적반응변수)라면 평균이 아닌 각 트리가 제시하는 값의 다수결(voting)로 예측한다.

Key Questions: 어떻게 트리가 가진 문제점을 보완할까?

회귀트리와 분류트리는 확보한 데이터의 1/2을 훈련데이터로 사용하기 때문에 최종적으로 구한 pruning 트리 자체도 훈련데이터 샘플에 따라 변동이 있을 뿐아니라 검정MSE의 경우 큰 변동을 가진다. 나는 이 문제점에 대해 이전 글회귀트리(regression tree)와 분류트리(classification tree) 예측의 “의문점 및 개선점”에서 개선점과 함께 언급했었다. 1. 어떻게 하면 훈련셋 변동의 문제를 해결하면서 트리의 예측력을 유지할 수 있을까? 회귀분석에서는 이 답을 교차검증(10 folds CV)에서 찾았다. 그러나 트리의 경우 pruning할 트리 선정을 위해 훈련셋 내에서 교차검증을 사용하게 된다. 따라서 교차검증을 통해 트리의 갯수를 정했는데 또 다시 교차검증을 통해 새로운 트리의 갯수를 찾는 것은 역설적인 상황이 된다. 2. 그렇다면 교차검증을 완료한 상황에서 훈련셋 검증은 어떻게 해야할까?

결론부터 말하면 Bootstrap 검증은 1번 문제인 훈련셋 변동성을 해결하며 예측력 또한 유지할 수 있다.(Bootstrap 검증은 CV검증과 함께 대표적인 검증법이다.) 또한 Bootstrap 검증시 자연스럽게 사용되지 않는 데이터(Out-of-Bag)를 검정셋으로 사용하면 2번 문제인 검증오차 검증도 할 수 있다.

Bootstrp 검증이란?

표본데이터를 통해 대부분의 경우 알 수 없는 모집단을 추정하는 것은 통계학의 주요 관심사이다. 표본데이터의 평균이 어느 정도의 변동을 가지는가?를 나타낸 것이 평균값의 표준오차이다. 표본데이터가 특정 분포(t분포, z분포 등)를 따른다고 가정하면 앞서 구한 표준오차를 통해 모집단을 추정할 수 있게 된다.

여기서 좀더 확장해보자. 만약 표본데이터를 가지고 표본데이터 갯수만큼 복원추출하면 모집단 뿐만아니라 표본집단과 별 상관이 없어 보일 것이다. 그런데 만약 이 복원추출한 집단을 수백개 만들어 그것의 평균과 평균값의 분산을 모집단의 그것들과 비교하면 어떨까? 결론은 이 비상한 아이디어가 통한다는 것이다. 즉, 모집단의 평균과 평균의 분산은 표본데이터의 수백개의 복원추출로 추정이 가능하게 된다. 우리는 대체로 모집단을 알 수 없다. 그러나 모집단을 추정할 또 하나의 좋은 방법을 가지게 된 것이다.

실제 정말 그런지 구현해 보자. ISLR 패키지의 Portfolio는 두 펀드의 투자수익 X와 Y의 자료 100개로 구성된다. 이 두 펀드에 투자해 최상의 수익을 올리기 위해서는 내 투자 금액의 비율(알파)은 얼마가 되야할까? 알파의 분산이 작다면 사용자는 리스크를 최소화하면서 최적값을 예측할 수 있을 것이다. 최적의 알파값은 (식 5.6)와 같이 각 데이터의 분산과 관련 있다. 모집단의 알파 평균은 0.6이다. 알파 추정값 1000번개의 표준오차는 0.083으로 알려져 있다.

boot 패키지의 boot 함수를 통해 이것을 확인해 보자. boot함수는 (데이터, 함수, bootstrap 갯수)로 이루어진다. 아래와 같이 bootstrap을 1000개를하면 평균은 0.576, 표준편차는 0.0886으로 모집단에서 추출한 표본의 평균과 평균의 오차는 거의 유사함을 알 수 있다.

library(ISLR)
library(boot)

alpha.fn=function(data,index){
  X<-data$X[index]
  Y<-data$Y[index]
  return((var(Y)-cov(X,Y))/(var(X)+var(Y)-2*cov(X,Y)))
}
set.seed(1) #동일한 sample을 1000번 하기 위해 필요
alpha.fn(Portfolio,sample(100,100,replace=T)) #복원추출 명시 필요 

## [1] 0.5963833

boot(Portfolio,alpha.fn,R=1000)

## 
## ORDINARY NONPARAMETRIC BOOTSTRAP
## 
## 
## Call:
## boot(data = Portfolio, statistic = alpha.fn, R = 1000)
## 
## 
## Bootstrap Statistics :
##      original        bias    std. error
## t1* 0.5758321 -7.315422e-05  0.08861826

’’

상기 그림은 실제 위에서 증명한 모집단의 표본(왼쪽)과 하나의 자료를 1000번 복원추출한 Bootstrap(중앙)의 히스토그램이다. boxplot을 보면 결론적으로 두 값의 평균 뿐아니라 분산도 유사함을 알 수 있다.

배깅: Boostrap + Tree 앙상블

상기 배깅 알고리즘4와 같이 각 Boostrap은 prune되지 않은 트리이다. 각각의 트리는 분산은 크지만 편향은 작다. 알고리즘5와 같이 m개의 트리들을 평균하는 것은 작은 편향은 유지한 채분산을 줄이는 효과를 가진다. 이것이 배깅이 작동하는 이유이며 대중의 지혜(wisdom of the crowd)라 불리는 이유이다. 각각의 개인은 똑똑치 못하나 개인이 합은 똑똑할 수 있게 된다.

Out-of-Bag(OOB) 오차 추정

복원추출방법은 n개 데이터중 1/3만큼은 사용하지 않게 된다. 어떻게 그것을 아냐고? 간단히 1에서 100까지의 100개 데이터를 복원추출하는 코딩을 짜보자. 여기서는 j=23이라하고 23이 한번도 들어있지 않을 확률이 얼마인지 확인한다.

a<-rep(NA,50000)
for(i in 1:50000){
  a[i]<-sum(sample(1:100,rep=TRUE)==23)>0}
1-mean(a)

## [1] 0.36346

결과는 전체 데이터중 36.5%는 배깅에 사용되지 않는 데이터(OOB)를 의미한다. 실제 수학을 통해 구한 j번째 관측치가 복원추출표본에 있지 않을 확률은 (1-1/n)^n이다. OOB오차는 그것을 사용하지 않은 데이터로 도출한 트리로 예측하기 때문에 검정오차의 유효한 추정치가 된다.

랜덤포리스트란?

랜덤포리스트는 배깅과 같다. 유일한 차이점이라면 데이터 샘플을 복원추출할때 모든 설명변수를 다 사용하는 배깅과 달리 a개의 설명변수만 고려하여 분할을 고려한다. 보통 설명변수 갯수 a는 전체 변수 p의 제곱근을 사용한다.(예: 13개 변수라면 4개 정도 사용) 심증이 강한 설명변수라고 특혜를 주지 않고 랜덤하게 변수를 사용된다. 랜덤포리스트는 일반적으로 배깅보다 나은 성능을 가진다. 설명변수가 많을 경우, 대체로 변수간 상관성이 높은 변수가 섞일 확률이 높은데 랜덤포리스트는 그 가능성을 없애기 때문이다.

결론

’’

상기 그림은 앞서 설명한 배깅과 랜덤포리스트 그리고 각각의 OOB의 검정오차의 비교를 보여준다. 랜덤포리스트의 경우 설명변수 갯수a는 p의 제곱근을 사용했다. 랜덤포리스트의 성능이 배깅보다 우수하고 OOB의 경우 검정오차를 크게 줄일 수 있음을 알 수 있다.

결론에 대한 의문점과 최종 결론

1. bootstrap의 수 m은 어느 정도 커야할까? 답: 결론적으로 m이 100이상이면 충분하지만 검정오차가 안정화될 만큼 큰 값을 사용하자.(예:400개 정도) m이 크다고 과적합을 우려하지 않아도 된다. 2. 랜덤포리스트의 설명변수 갯수 a는 얼마가 적당한가? 답: 회귀트리는 1/3p, 분류트리는 p의 제곱근을 주로 쓰나 실제 이 값의 주변 값을 확인해 볼 필요가 있다. 변수간 상관성에 따라 최적의 a값이 다를 수 있기 때문이다. 3. 배깅이 랜덤포리스트보다 성능이 더 좋지 않다면 왜 쓰는걸까? 답: 현재까지 예를 찾기 어렵다. 굳이 찾자면 설명변수가 5개 이하의 많지 않은 경우 정도이겠다. 따라서 단독 분류트리, 단독 회귀트리뿐 아니라 배깅도 쓸 필요가 없다.
> 결론은 랜덤포리스트의 OOB검정을 사용하되 a값을 조정하여 최선의 a를 찾자.

실전 예제

이전 회귀트리에서 사용한 Boston 집 값 예제를 가져오자(MASS) 또한 randomForest 패키지를 설치한다. 이전 글에서 구현한 회귀트리의 검정MSE 25와 비교해 보자. randomForest함수의 ntree의 default는 500이다.

배깅의 검정MSE

library(randomForest)
library(MASS)
set.seed(1)
train<-sample(1:nrow(Boston),nrow(Boston)/2) #회귀트리와 같아 1/2만 훈련셋
Boston.test<-Boston[-train,"medv"]
bag.boston<-randomForest(medv~.,Boston,subset=train,mtry=13)
yhat.bag<-predict(bag.boston,Boston[-train,])
mean((yhat.bag-Boston.test)^2)

## [1] 13.17605

랜덤포리스트 검정 MSE

set.seed(1)
rf.boston<-randomForest(medv~.,Boston,subset=train,mtry=6)
yhat.rf<-predict(rf.boston,Boston[-train,])
mean((yhat.rf-Boston.test)^2)

## [1] 11.50837

랜덤포리스트 OOB MSE 아래 결과와 같이 변수가 6개일 때 OOB MSE가 최소가 된다.

MSE<-rep(NA,6)
set.seed(1)
rf.boston5<-randomForest(medv~.,data=Boston, mtry=5, importance=TRUE, do.trace=500)

##      |      Out-of-bag   |
## Tree |      MSE  %Var(y) |
##  500 |    9.757    11.56 |

rf.boston6<-randomForest(medv~.,data=Boston, mtry=6, importance=TRUE, do.trace=500)

##      |      Out-of-bag   |
## Tree |      MSE  %Var(y) |
##  500 |    9.205    10.90 |

rf.boston7<-randomForest(medv~.,data=Boston, mtry=7, importance=TRUE, do.trace=500)

##      |      Out-of-bag   |
## Tree |      MSE  %Var(y) |
##  500 |    9.352    11.08 |

• 회귀트리 검정 MSE: 25
• 배깅 검정MSE: 13.17
• 랜덤포리스트 검정MSE: 11.51 (변수 6개일 때)
• 랜덤포리스트 OOB MSE: 9.352 (변수 6개일 때)

각각의 값의 랜덤포리스트의 결과값인 Mean of Squared Residuals는 OOB MSE를 의미한다. 상기 253개의 훈련셋을 사용한 OOB MSE는 전체를 사용하지 않았기에 의미가 없다.

부록: 상대적 변수의 중요도

랜덤포리스트는 회귀모델에 비해 각 변수의 해석이 어렵다. 그러나 RSS(회귀트리)나 이탈도(지니지수, 분류트리)를 사용하면 각 설명변수의 중요도를 알 수 있다. 예를 들어 회귀트리의 경우 주어진 변수가 모델에서 제외될 때 검정RSS의 평균 감소량 크다면 설명변수가 중요하다 말할 수 있다.위에서 얻은 최적의 모델 6개의 트리를 가진 OBB 검정 모델에서 importance를 추가하여 변수의 상대적 중요도를 plot해보자 importance에 나오는 첫열은 검정MSE의 증가량이고 두번째 열은 변수에 대한 분할로 인한 노드 impurity의 평균이다. 이것은 회귀트리의 경우 훈련MSE이고 분류트리는 이탈도이다. lstat(재산수준)과 rm(주택크기)가 가장 중요한 두변수 임을 알 수 있다.

set.seed(1)
rf.boston6<-randomForest(medv~.,data=Boston, mtry=6, importance=TRUE)
importance(rf.boston6)

##           %IncMSE IncNodePurity
## crim    16.179977    2175.26986
## zn       3.412097      71.87389
## indus   10.216919    2040.22314
## chas     3.575822     142.43428
## nox     19.563438    2449.25161
## rm      43.704088   14155.07004
## age     11.959347     869.74253
## dis     21.182821    2350.13055
## rad      5.849383     210.91891
## tax     12.444655     834.90365
## ptratio 15.102880    2042.47655
## black    8.231795     595.72178
## lstat   35.025317   14299.96408

varImpPlot(rf.boston6)

lstat를 100%로 했을때의 상대 크기는 아래와 같다.

rel.importance<-rf.boston6$importance[,1]/max(rf.boston6$importance[,1])*100
rel.importance.dec<-sort(rel.importance, decreasing=T)
barplot(rel.importance.dec,horiz=T,names=names(rel.importance.dec))