알고리즘 소개 : XGBoost

알고리즘 소개 : XGBoost

 

XGBoost(eXtreme Gradient Boosting)는 병렬처리와 최적화를 장점으로 내세우는 Gradient boosting 알고리즘 으로 릴리즈된 이래 Kaggle 대회에서 좋은 성적을 보이며 많은 관심을 끈 방법론입니다. LightGBM, CatBoost 등 다른 gradient boosting 알고리즘이 나오면서 어떤 것이 더 좋은 성능을 보이는지는 계속 논란이 되고 있지만 XGBoost는 여전히 준수한 성능 및 속도와 information gain 기반 importance 산출 기능을 제공해 많은 사랑을 받고 있습니다. 이번 포스팅에선 이 XGBoost에 대해 소개해보고자 합니다.

 

1. CART

분명 시작글과 제목에선 XGBoost를 소개한다고 했는데, 갑자기 CART란 단어가 나옵니다. 왜 나왔을까요? 이유는 XGBoost가 CART(Classification And Regression Trees)를 기반으로 만들어진 모델이기 때문입니다. 알고리즘을 알아보려면 먼저 그 근간부터 소개를 하는 것이 맞겠죠.

그림출처 : XGBoost : A Scalable Tree Boosting System

CART는 이름 그대로 분류 및 회귀 나무’들’로 트리 기반 앙상블 모델을 의미합니다. 앙상블 모델은 다수의 학습방법을 이용하여 결론을 내리는 방법으로, CART는 이 ‘다수의 학습방법’을 여러가지 의사결정나무(Decision Tree)로 정의한 방법론입니다. (앙상블 모델에서 Bagging과 Boosting의 차이, Stacking과 앙상블의 차이 등은 추후 포스팅 해보도록 하겠습니다. 궁금하신 분은 Wikipedia나 정리된 다른 글들을 참고해보세요!)

이 CART의 원리는 Additive learning이라고 정의하곤 합니다. 말이 어려워 보이지만 별 특별한 의미 없이 ‘더해서 + 배운다’입니다. 이를 CART의 상황에 대입해보면 밑의 수식처럼 표현할 수 있습니다.

Gradient boosted tree의 아이디어

Y’는 타겟(Y)에 대한 예측값을, a,b,c…는 각각 트리 A,B,C…에서 나온 가중치들을 말합니다. 이 개념을 XGBoost의 Gradient Boosting Tree로 가져가면 밑과 같이 표현됩니다.

XGBoost의 오브젝트 함수와 y hat

수식이랑 변수 설명으로 덕지덕지 붙여 놨지만 기본 개념은 같습니다. f함수는 알고리즘에서 생성해낸 의사결정나무 모델들을 의미하고 l은 loss function을, regularization term은 나무들이 얼마나 복잡할지에 대해 정의하는 파라미터입니다. 요점은 ‘의사결정 나무 모델을 여러 개 학습시켜서 예측값을 더한 것으로 결정하겠다!’ 인 셈이죠. 이렇게 더해진 예측 점수들을 이용해 결론을 내림으로서 과적합이나 기존 모델이 잘 설명하지 못하는 취약부분에 대한 보완을 할 수 있는 힘을 만듭니다.

 

2. 하이퍼 파라미터와 XGBoost

그렇다면 위의 CART면 모든 것이 해결되느냐? 하면 그렇지 않습니다. 여러가지 모델을 생성할 때 무슨 기준으로 만들어야 하는지, 위에 말한 regularization term(모델 복잡도)은 어떻게 정해야 하는지, gradient boosting의 꽃인 learning rate은 어떻게 설정해야 할지 등 많은 문제가 남아있습니다. XGBoost는 이를 위해 유저가 설정해주는 하이퍼 파라미터가 존재합니다. 복잡한 블랙박스 모델들이 으레 그렇듯, XGBoost도 이 하이퍼 파라미터를 얼마나 잘 조정해주느냐가 성능을 끌어올리는데 중요한 역할을 합니다. XGBoost documentation의 파라미터 페이지 (https://xgboost.readthedocs.io/en/latest/parameter.html)에 모든 인풋값들과 그 의미가 잘 정리되어 있지만, 이 모든 것을 소개하려면 너무 많은 시간이 소요되니 자주 쓰고 중요한 것들 몇가지만 소개해보도록 하겠습니다.

 

-     eta : learning rate, default = 0.3, [0, 1]

eta는 ‘Step size shrinkage’로, 학습 단계별로 가중치를 얼마나 적용할 지 결정하는 숫자입니다. 가중치이므로 0~1 사이의 값을 가지며, 낮을수록 more conservative, 즉 보수적인 모델이 됩니다(다음 단계의 결과물을 적게 반영하기 때문).

-     gamma : min split loss, default = 0, [0, ∞]

Gamma를 이해하기 위해서는 먼저 information gain이라는 것을 이해해야 합니다. Information gain은 의사결정나무가 가지를 칠 때, 즉 새로운 변수를 기준으로 데이터를 분류할 때 타겟 변수에 대해 얼마나 설명할 수 있는지를 측정하는 기준입니다.

XGBoost에서의 Information Gain 수식

여기서 감마는 맨 끝의 숫자로, 계산한 information gain에 페널티를 부여하는 숫자입니다. 즉 이 값이 커질수록 의사결정나무들은 가지를 잘 만들려 하지 않게 되며, 이에 따라 보수적인 모델이 되게 됩니다.

-     max_depth : default = 6, [0, ∞]

max_depth는 말 그대로 의사결정나무의 깊이의 한도입니다. 기본 값은 6으로 커질수록 더 복잡한 모델이 생성되며, 이는 overfitting의 문제를 일으킬 수도 있습니다. 물론 이 값은 어디까지나 한계치이기 때문에 무작정 늘린다고 모든 나무들이 많은 가지를 생성해내지는 않습니다.

-     subsample : default = 1, (0, 1]

training 데이터셋에서 subset을 만들지 전부를 사용할지를 정하는 파라미터입니다. 매번 나무를 만들 때(=iteration) 적용하며 overfitting 문제를 방지하려고 사용합니다.

-     colsample_bytree : default = 1, (0, 1]

나무를 만들 때 칼럼, 즉 변수를 샘플링해서 쓸지에 대한 파라미터입니다. 나무를 만들기 전 한 번 샘플링을 하게 됩니다.

-     scale_pos_weight : default = 1, (0, 1]

분류 모델에서 사용하는 가중치 파라미터로, 극단적으로 적은 타겟값이 존재하는 문제에서 유용합니다.

 

이러한 하이퍼 파라미터는 대부분 값이 정해진 규칙이 없다보니 주로 랜덤서치나 그리드 서치로 최적값을 찾게되며, 근래에는 딥러닝에 대한 관심이 증가하면서 베이지안 최적화를 이용한 방법도 주목을 받고 있습니다. 어떤 방법이든 heuristic한 방법들이기 때문에 많은 연산을 필요로 하며, 따라서 무작정 넓은 값을 찾기보다는 현 데이터 상황에 맞는 중요한 파라미터들을 선별하고 그 의미를 알고 있는 상태에서 찾는 것이 중요합니다.

 

3. XGBoost in R

이제 본격적으로 모델을 어떻게 만드는지 살펴보겠습니다. 데이터는 맨날 쓰던 iris 데이터를 사용할 것이며, 분류 문제와 회귀 문제 둘 다 해보도록 하겠습니다.

3.1 XGBoost를 위한 Data 준비

XGBoost는 매트릭스 기반 인풋을 사용하고, 파이썬과 같이 0부터 시작하는 인덱스를 사용하기 때문에 R 환경에서 사용하려면 데이터 준비 부분에서 주의를 많이 해야합니다. xgboost 패키지에서는 xgb.DMatrix 형태의 객체로 인풋 데이터를 정의하나, 보다 간편히 다루기 위해서 다음과 같이 형식을 취하도록 하겠습니다.

R에서의 XGBoost 데이터 인풋 정의

X는 칼럼 이름을 가지는 데이터 메트릭스 형식이며, Y는 라벨 벡터, 즉 타겟값을 벡터 형태로 넣는 것입니다. 이 형식으로 모델을 생성하기 앞서 객체를 정의하는 것이 편한데, 이유는 R에서 다루는 대부분의 데이터들이 데이터프레임 형식으로 되어있기 때문입니다. 특히, R의 기초 패키지들은 factor변수들을 자동으로 dummy variable(범주형 변수를 0 또는 1로 구성된 값들로 변형)화 시켜주나, xgboost는 매트릭스 형식으로 데이터를 받다보니 그런 기능이 존재하지 않아서 직접 처리를 해주어야 합니다.

3.2 Multi-class classification using XGBoost in R

iris의 Species 변수가 3가지 클래스를 가지고 있기 때문에 binary 분류 모델이 아닌 다중 분류모델을 사용해서 풀어보도록 하겠습니다. 먼저 위에 언급한대로 iris 데이터를 x와 y로 쪼개서 정의하도록 하겠습니다.

iris 데이터셋에서 분류모델을 위한x와 y 정의

이제 이 데이터들을 이용해 5-fold cross validation을 하는 xgboost 모델을 만들어보겠습니다. xgb.cv는 xgboost기반 cross valdiation을 진행하도록 해주는 함수로, 기존의 xgboost, xgb.train 함수가 가지는 인풋 이외에 nfold(몇 개의 폴드(데이터 subset)로 진행할지), prediction(예측값을 만들어낼지 말지) 등의 인풋을 받게 됩니다.

xgboost 5 fold Cross validation with iris data

이제 원하는 대로 모델 인풋을 정합니다. 주의할 점은 label(y데이터)을 넣을때 아직 factor값으로 되있는 y를 숫자로 바꾼후 1을 빼줘야 한다는 것입니다. 이유는 파이썬처럼 xgboost가 0부터 시작하는 값을 받기 때문에 as.numeric을 해서 만든 1,2,3....이 아닌 1을 뺀 0,1,2....를 넣지 않으면 에러가 뜨기 때문입니다. num_class는 복수의 클래스를 가진 분류모델을 짤 때 쓰는 값으로, 클래스가 몇 개가 존재하는지 알려주는 값입니다. levels 함수를 이용해 y값의 층계가 몇 개 존재하는지 확인 후 length 함수를 이용하여 길이, 즉 몇 개 클래스가 있는지를 입력해주면 됩니다.

다음 줄에서 nfold는 몇 번 fold(접다), training/validation으로 나누는 비율을 몇으로 할지 정하는 값으로, 보편적으로 사용하는 5fold를 사용하였습니다. nrounds는 모델에서 iteration을 몇 번 진행할지이며 마지막으로 early_stopping_rounds는 해당 값 이후의 iteration에서 평가지표가 수렴한다고 생각되면 모델을 멈추게 하는 값입니다. 

세번째 줄에서는 objective는 어떤 목적을 가지고 학습을 진행할지 이며 여기서는 multi(멀티 클래스) : softprob(soft probability)를 사용합니다. 원한다면 multi:softmax를 사용하셔도 큰 문제는 없습니다. eval_metric은 모델 평가지표(evaluation metric)로 기본은 merror(multiclass error)로 되있으나 저는 mlogloss(multiclass log loss)로 사용하겠습니다. 마지막의 verbose는 학습마다 평가값들에 대한 메세지를 출력할지이며 prediction은 Cross validation이 모두 끝나면 예측된 값들을 출력할지 입니다.

이렇게 모델을 학습하면 결과물이 list형식으로 저장되며 리스트는 다음과 같은 값들을 가지고 있습니다.

CV 모델의 요소들

여기서 대부분의 값들은 굳이 확인 안하더라도 모델을 세팅하는 과정에서 이미 알 수 있는 것이며, 관심을 가질만한 것들은 evaluation_log와 pred입니다. 둘다 말 그대로 모델 iteration별 평가 로그와 예측값을 나타내며, 다음과 같은 출력값이 나오게 됩니다.

evaluation log 오브젝트

pred오브젝트, multi:softmax 방법으로 진행시 이미 최고 확률로 값을 정하기 때문에 벡터로 나옵니다.

여기서 저는 softprob 방법을 사용했기 때문에 일반적인 분류 모델에서 보는 결과를 보려면 약간의 수정이 필요합니다. 이 수정 방법은 각 행마다 3개 열(클래스) 중 최고 값을 가지는 것을 선택하는 것으로, softmax 방법에서는 모델 내부에서 진행해주는 것입니다. 다음과 같은 코드로 분류 결과를 낼 수 있습니다.

분류 모델 결과 데이터프레임인 pred_df 생성

pred_df를 활용한 결과 리포트

위의 그림에서 보다시피 단순 테이블로 정리 후 결과를 손으로 직접 써서 확인할 수도 있으며, 혹은 caret 패키지의 confusionMatrix 함수를 활용하여 리포트 기능을 강화할 수도 있습니다. iris 데이터가 워낙 깔끔히 정돈된 데이터라서 잘 맞추는 모습을 볼 수 있습니다.

3.3 XGBoost Cross validation 결과 시각화

위와 같이 모델을 생성하고 결과를 간단한 confusion matrix로 만들어 보았으나, 모델의 학습과정이 잘 되었는지, iteration별로 값이 어떻게 변화하였는지는 알기 힘듭니다. 이를 위해서는 evaluation_log 객체를 활용하여 시각화를 하면 쉽게 알아볼 수 있으며, 다음과 같은 코드로 생성 가능합니다.

cross valdiation 결과 시각화를 위한 cvplot 함수

이 함수에 아까 학습한 모델 객체를 넣으면 다음과 같은 그래프가 생성됩니다

Multiclass classification 시각화 결과

train과 test의 mlogloss 값이 일정 라인에서 멈춰있음을 볼 수 있습니다. early_stopping_rounds 옵션에 의해 163 iteration까지 학습 후 모델은 값이 정체되있음을 확인하고 그만 학습했음을 알 수 있습니다.

3.4 Regression using XGBoost in R

이제 xgboost를 이용해 회귀모델을 만들도록 해보겠습니다. 데이터는 iris로 동일하게 사용하며, 타겟은 Sepal.Width를 사용하도록 하겠습니다.

여기서 x데이터에 Species라는 범주형 변수가 포함되어 있기 때문에 이 변수를 one-hot encoding 형태로 변형해줘야 합니다. 이를 위해서 model.matrix라는 함수를 이용하여 해당 데이터를 만든 후 나머지 데이터에 붙이는 식으로 진행하겠습니다.

회귀분석을 위한 데이터 정제

이제 분류모델과 똑같이 모델을 정의하도록 하겠습니다. 크게 달라질 것은 없으나 다중 분류모델이 아니기에 num_class를 정의할 필요가 없으며 objective는 reg:linear(linear regression)을 채택하겠습니다.

5-fold cross validation 회귀모델 정의 및 학습

학습한 결과를 위의 cvplot을 이용해 시각화하면 다음과 같이 결과를 볼 수 있습니다.

회귀 모델 cross validation 결과

3.5 하이퍼 파라미터 튜닝

이제 모델별로 테스트를 해봤으니 하이퍼 파라미터 튜닝에 대해 진행해보도록 하겠습니다. 이 예제에서는 무난하게 사용되는 그리드 서치를 활용하여, eta와  gamma 값을 튜닝하도록 하겠습니다.

그리드 서치를 위해서는 먼저 그리드, 즉 탐색할 변수들에 대한 조합이 필요합니다. R에서는 expand.grid 함수를 이용하면 간편하게 그리드를 생성할 수 있습니다.

그리드 서치를 위한 그리드(서치 범위) 생성

그리드가 생성되면 이제 이 범위만큼 모델을 학습시켜야 합니다. 그냥 for문으로 진행할 수도 있지만 빠른 서치를 위해 foreach와 doParallel 패키지를 활용하여 병렬처리를 하도록 하겠습니다.

그리드 서치 병렬처리문

cl은 병렬처리를 위해 코어들을 할당하는 클러스터이며, registerDoParallel은 해당 클러스터를 작업하기 위해 할당해주는 것, foreach는 할당된 클러스터를 가지고 병렬작업을 진행하는 함수입니다. 각 grid 별로 파라미터 값을 적용하여 모델을 학습시키고, 학습된 모델의 train과 test rmse의 마지막 값(수렴된 값)을 기준으로 가져오게 됩니다.

*R에서 foreach와 병렬처리에 대해서는 기회가 되면 apply, Reduce 등과 묶어서 한번에 포스팅 해보도록 하겠습니다

서치한 결과 중 test 결과값이 최소인 것을 보려면 which.min 함수를 사용하여 값을 가져오면 됩니다. which.min은 해당 벡터 중 최소인 것을 찾아주는 함수로, R과 데이터프레임 포스팅에서 사용했었던 which의 아종이라고 생각하시면 됩니다.

그리드 서치한 결과에서 가장 낮은 test rmse를 가지는 시도와 해당 서치에서 사용한 하이퍼 파라미터

 

3.6 XGBoost와 변수 중요도

흔히 XGBoost와 같은 블랙박스 모델들의 단점은 계산과정을 추적하기가 힘들어 변수 간 중요도를 수치화 하기 어렵다는 것입니다. 하지만 XGBoost에선 이 변수 중요도를 수치화하여 제공하고 있습니다. xgb.importance 함수는 Gain(위의 Information gain), Cover, Frequency를 기준으로 변수들의 중요도를 데이터프레임화 시켜줍니다. 단, xgb.importance는 cross validation 객체(xgb.cv)는 받지 않으므로 xgboost 함수로 다시 모델을 만들어줘야 합니다

모델 재생성과 imprtance dataframe 생성

importance 객체

이 importance 객체는 xgb.plot.importance 함수로 barplot을 그릴 수 있으나, 썩 좋은 plot은 아니기에 직접 그려보도록 하겠습니다.

xgboost importance plot

 

4. 마치며

처음 이 주제를 생각할땐 대학원에서 연구과제를 할 때 공부했던 주제여서 쉬울 것이라 여겼는데 막상 글을 써보니 쉽게 와닿도록 설명하기가 힘든 주제라는 것을 느끼게 된 것 같습니다.

다음 주제는 본문에서 언급했던 R에서의 병렬처리를 해볼까 생각중입니다. 

*주제는 소재 및 상황에 따라 변경될 수도 있습니다

 

[컨텐츠 코드]

	#==========================================================================
	# Topic : XGBoost
	# Date : 2019. 03. 30
	# Author : Junmo Nam
	#==========================================================================
	#==========================================================================
	# Load data and packages
	#==========================================================================
	sapply(c('xgboost','dplyr', #xgboost for modeling, dplyr for data managing
	'foreach','doParallel' #parallel packages for grid search
	),require,character.only = T)
	x = iris %>% select(-Species) %>%
	data.matrix #matrix format
	y = iris$Species
	#==========================================================================
	# Cross validation with whole data : multiclass classification
	#==========================================================================
	#training model
	cv_model1 = xgb.cv(data = x, label = as.numeric(y)-1, num_class = levels(y) %>% length, # claiming data to use
	nfold = 5, nrounds = 200, early_stopping_rounds = 150, # about folds and rounds
	objective = 'multi:softprob', eval_metric = 'mlogloss', # model options
	verbose = F, prediction = T # do not print messages while training, make prediction
	)
	#result : output as list
	attributes(cv_model1)
	#evaluation log
	cv_model1$evaluation_log
	#prediction matrix(for softprob) to dataframe
	cv_model1$pred
	pred_df = cv_model1$pred %>% as.data.frame %>%
	mutate(pred = levels(y)[max.col(.)] %>% as.factor,actual = y)
	#table
	pred_df %>% select(pred,actual) %>% table #table
	caret::confusionMatrix(pred_df$pred,pred_df$actual) #confusion matrix
	#visualizing model
	cvplot = function(model){ #visualizing function
	eval.log = model$evaluation_log
	std = names(eval.log[,2]) %>% gsub('train_','',.) %>% gsub('_mean','',.)
	data.frame(error = c(unlist(eval.log[,2]),unlist(eval.log[,4])),
	class = c(rep('train',nrow(eval.log)),
	rep('test',nrow(eval.log))),
	nround = rep(1:nrow(eval.log),2)
	) %>%
	ggplot(aes(nround,error,col = class))+
	geom_point(alpha = 0.2)+
	geom_smooth(alpha = 0.4,se = F)+
	theme_bw()+
	ggtitle("Apple's Rbox : XGBoost Cross-validation Visualization",
	subtitle = paste0('fold : ',length(model$folds),
	' iteration : ',model$niter
	)
	)+ylab(std)+theme(axis.title=element_text(size=11))
	}
	cvplot(cv_model1)
	#==========================================================================
	# Cross validation with whole data : numeric target regression
	#==========================================================================
	#which one to be a numeric target in iris??
	cor(iris %>% select(-Species)) %>% corrplot::corrplot(method = 'number')
	#reclaiming data
	x_num = iris %>% model.matrix(~0+Species,.) %>% as.data.frame %>% #one-hot encoding matrix
	bind_cols(iris %>% select(-Species,-Sepal.Width)) %>% data.matrix
	y_num = iris$Sepal.Width
	#training model
	cv_model2 = xgb.cv(data = x_num, label = y_num,
	nfold = 5, nrounds = 200, early_stopping_rounds = 150,# about folds and rounds
	objective = 'reg:linear',verbose = F, prediction = T
	)
	#check performance by plot
	cvplot(cv_model2)
	#check performance by prediction (Mean absolute percetage error)
	abs((y_num - cv_model2$pred) / y_num) %>% mean #7%
	#==========================================================================
	# Hyper parameter tuning : grid search
	#==========================================================================
	#make grid
	grid = expand.grid(eta = seq(0.1,0.4,0.05),gamma = seq(0,5,1))
	#make parallel cluster and register
	cl = makeCluster(detectCores()-1) #n-1 core cluster
	registerDoParallel(cl)
	#do search by parallel
	grid_search = foreach(i = 1:nrow(grid),.combine = rbind,.packages = c('dplyr','xgboost')) %dopar% {
	model = xgb.cv(data = x_num, label = y_num,
	nfold = 5, nrounds = 200, early_stopping_rounds = 150,# about folds and rounds
	objective = 'reg:linear',verbose = F, prediction = T,
	params = grid[i,]
	)
	data.frame(train_rmse_last = unlist(model$evaluation_log[,2]) %>% last,
	test_rmse_last = unlist(model$evaluation_log[,4]) %>% last)
	}
	stopCluster(cl) #end parallel cluster
	grid_search[which.min(grid_search$test_rmse_last),]
	grid[which.min(grid_search$test_rmse_last),]
	#==========================================================================
	# Variable importance using xgb.importance
	#==========================================================================
	#xgbooster model using grid search parameter
	model = xgboost(data = x_num, label = y_num,
	nrounds = 200, early_stopping_rounds = 150,# about folds and rounds
	objective = 'reg:linear',verbose = F,
	params = grid[which.min(grid_search$test_rmse_last),])
	#make importance dataframe
	imp = xgb.importance(model = model)
	imp
	xgb.plot.importance(imp)
	#customized visualization
	data.frame(variable = rep(imp$Feature,3),
	value = c(imp$Gain,imp$Cover,imp$Frequency),
	Type = c(rep('Gain',nrow(imp)),rep('Cover',nrow(imp)),rep('Frequency',nrow(imp)))
	) %>% ggplot(aes(variable,value,fill = variable))+
	geom_bar(stat = 'identity')+
	facet_grid(~Type)+
	theme_bw()+
	ggtitle('XGBoost : Customized Importance Plot',
	subtitle = "Author : Apple's R box")

view raw xgboost.r hosted with ❤ by GitHub