CH 05. 분류
# setting
library(MASS)
library(dplyr)##
## Attaching package: 'dplyr'## The following object is masked from 'package:MASS':
##
## select## The following objects are masked from 'package:stats':
##
## filter, lag## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, unionlibrary(ggplot2)
library(FNN)
library(mgcv)## Loading required package: nlme##
## Attaching package: 'nlme'## The following object is masked from 'package:dplyr':
##
## collapse## This is mgcv 1.8-31. For overview type 'help("mgcv-package")'.library(rpart)
library(klaR)
PSDS_PATH <- file.path('/Users/jinameliachoi/Documents/statistics-for-data-scientists')
loan3000 <- read.csv(file.path(PSDS_PATH, 'data', 'loan3000.csv'))
loan_data <- read.csv(file.path(PSDS_PATH, 'data', 'loan_data.csv'))
loan_data$outcome <- ordered(loan_data$outcome, levels=c('paid off', 'default'))
full_train_set <- read.csv(file.path(PSDS_PATH, 'data', 'full_train_set.csv'))
full_train_set$outcome <- ordered(full_train_set$outcome, levels=c('paid off', 'default'))나이브 베이즈 (Naive Bayes)
주어진 결과에 대해 예측변수 값을 관찰할 확률을 사용하여,
예측변수가 주어졌을 때, 결과 Y = i를 관찰할 확률을 추정한다.
나이브 베이즈 계산 방법:
이진 응답(Y가 0 아니면 1)에 대해, 각 예측변수에 대한 조건부확률을 구한다.
이 확률은 training dataset에서 Y = i인 레코드 중 X의 값의 비율로 구할 수 있다.각 확률값을 곱한다음, Y = i에 혹한 레코드들의 비율을 곱한다.
모든 클래스에 대해 1~2단계를 반복한다.
2단계에서 모든 클래스에 대해 구한 확률 값을 모두 더한 값으로 클래스 i의 확률을 나누면 결과 i의 확률을 구할 수 있다.
이 예측변수에 대해 가장 높은 확률을 갖는 클래스를 해당 레코드에 할당한다.
서로 모든 X들이 독립이라고 가정한 것: 예측변수 벡터의 정확한 조건부확률은 각 조건부확률의 곱으로도 충분히 잘 추정할 수 있다고 가정함
나이브 베이즈 모델을 구하기 위한 KlaR 패키지 사용
naive_model <- NaiveBayes(outcome ~ purpose_ + home_ + emp_len_, data=na.omit(loan_data))
naive_model$table## $purpose_
## var
## grouping credit_card debt_consolidation home_improvement major_purchase
## paid off 0.18759649 0.55215915 0.07150104 0.05359270
## default 0.15151515 0.57571347 0.05981209 0.03727229
## var
## grouping medical other small_business
## paid off 0.01424728 0.09990737 0.02099599
## default 0.01433549 0.11561025 0.04574126
##
## $home_
## var
## grouping MORTGAGE OWN RENT
## paid off 0.4894800 0.0808963 0.4296237
## default 0.4313440 0.0832782 0.4853778
##
## $emp_len_
## var
## grouping < 1 Year > 1 Year
## paid off 0.03105289 0.96894711
## default 0.04728508 0.95271492new_loan <- loan_data[147, c('purpose_', 'home_', 'emp_len_')]
row.names(new_loan) <- NULL
new_loan## purpose_ home_ emp_len_
## 1 small_business MORTGAGE > 1 Yearpredict(naive_model, new_loan)## $class
## [1] default
## Levels: paid off default
##
## $posterior
## paid off default
## [1,] 0.3463013 0.6536987이 경우 연체 default를 예상한다.
예측 결과에는 사후추정 posterior도 함께 있다.
수치형 예측변수의 경우에는 binning하여 범주형으로 변환 후, 알고리즘을 적용한다.
나이브 베이즈는 범주형 변수만 사용 가능
판별분석 (Discriminant Analysis)
판별분석에는 여러 가지 기법이 있지만, 그 가운데 가장 일반적인 방법은 선형판별분석 (linear discriminant analysis, LDA)이다.
트리 모델이나 로지스틱 회귀와 같은 더 정교한 기법이 출현하였지만,
여전히 일부 응용 분야에서는 LDA를 사용하고 있다.
예측변수들의 중요성을 측정하거나 효과적으로 특징을 선택하는 방법으로 사용된다.
공분산행렬 개념 이해
판별 분석을 이해하려면 두 개 이상의 변수 사이에 공분산이라는 개념을 먼저 도입해야 한다.
상관관계는 -1에서 1 사이에서 정의되었다면, 공분산은 두 개의 변수에서 사용하는 척도와 동일한 척도에서 정의된다.
피셔의 선형판별
두 개의 연속형 변수를 사용하여 이진 결과변수를 예측하는 분류 문제의 경우,
그룹 안의 편차와 다른 그룹간의 편차를 구분한다.
LDA는 ‘내부’ 제곱합에 대한 ‘사이’ 제곱합의 비율을 최대화 하는 것을 목표로 한다.
이 비율을 최대화 하는 선형결찹을 찾는 것이 두 그룹 사이를 가장 명확하게 나누는 방법이다.
MASS 패키지에서 LDA 함수를 제공한다.
두 예측변수 borrower_score와 payment_inc_ratio를 이용해 선형판별자 가중치 구하기
loan_lda <- lda(outcome ~ borrower_score + payment_inc_ratio, data=loan3000)
loan_lda$scaling## LD1
## borrower_score 7.17583880
## payment_inc_ratio -0.09967559상환과 연체에 대한 확률 계산 가능
pred <- predict(loan_lda)
head(pred$posterior)## default paid off
## 1 0.5535437 0.4464563
## 2 0.5589534 0.4410466
## 3 0.2726962 0.7273038
## 4 0.5062538 0.4937462
## 5 0.6099525 0.3900475
## 6 0.4107406 0.5892594예측 결과 시각화해서 본다면 LDA가 잘 작동하는지 쉽게 알 수 있을 것이다.
predict 함수를 사용하여, 체납에 대한 확률값을 그래프로 시각화해보자.
lda_df <- cbind(loan3000, prob_default=pred$posterior[, 'default'])
x <- seq(from=.33, to=.73, length=100)
y <- seq(from=0, to=20, length=100)
newdata <- data.frame(borrower_score=x, payment_inc_ratio=y)
pred <- predict(loan_lda, newdata=newdata)
lda_df0 <- cbind(newdata, outcome=pred$class)
ggplot(data=lda_df,
aes(x=borrower_score, y=payment_inc_ratio, color=prob_default)) +
geom_point(alpha=.6) +
scale_color_gradient2(low='white', high='blue') +
scale_x_continuous(expand=c(0,0)) +
scale_y_continuous(expand=c(0,0), lim=c(0, 20)) +
geom_line(data=lda_df0, col='green', size=2, alpha=.8) +
theme_bw()## Warning: Removed 18 rows containing missing values (geom_point).
로지스틱 회귀 (Logistic Regression)
로즈 오즈 함수, 또는 로짓 함수는 0과 1 사이의 확률을 무한대의 값으로 매핑해준다.
이후 예측가능한 선형모형을 마들고 컷오프(절사) 기준을 이용하여,
그 값보다 큰 확률이 나오면 1로 분류하는 식의 과정을 통해 클래스 라벨을 구할 수 있다.
로지스틱 회귀는 선형회귀를 확장한 일반화선형모델 (GLM)의 특별한 사례이다.
R에서 로지스틱 회귀를 구하려면 family 인수를 binomial로 지정하고 glm 함수를 사용한다.
logistic_model <- glm(formula = outcome ~ payment_inc_ratio + purpose_ + home_ + emp_len_ + borrower_score, family='binomial', data=loan_data)
logistic_model##
## Call: glm(formula = outcome ~ payment_inc_ratio + purpose_ + home_ +
## emp_len_ + borrower_score, family = "binomial", data = loan_data)
##
## Coefficients:
## (Intercept) payment_inc_ratio
## 1.63809 0.07974
## purpose_debt_consolidation purpose_home_improvement
## 0.24937 0.40774
## purpose_major_purchase purpose_medical
## 0.22963 0.51048
## purpose_other purpose_small_business
## 0.62066 1.21526
## home_OWN home_RENT
## 0.04833 0.15732
## emp_len_ > 1 Year borrower_score
## -0.35673 -4.61264
##
## Degrees of Freedom: 45341 Total (i.e. Null); 45330 Residual
## Null Deviance: 62860
## Residual Deviance: 57510 AIC: 57540outcome이 응답변수. 대출을 모두 갚았다면 0, 연체 상태이면 1의 값을 가진다.
purpose_와 home_으로 시작하는 변수들은 각각 대출목적과 주택소유상태를 나타내는 요인변수.
이러한 요인변수들에 대한 기준 수준은 각각 credit_card와 MORTGAGE
변수 borrower_score는 차용인의 신용도를 나타내는 0~1까지의 점수이다.
일반화선형모델 GLM은 회귀와 함께 두번째로 가장 중요한 모델이다.
이는 두가지 주요 구성 요소로 특징 지어진다.
확률분포 또는 분포군(로지스틱 회귀의 경우 이항분포)
응답을 예측변수에 매핑하는 연결함수(로지스틱 회귀의 경우 로짓)
로지스틱 회귀 모델에서 얻은 예측값과 변환을 통한 확률값
pred <- predict(logistic_model)
summary(pred)## Min. 1st Qu. Median Mean 3rd Qu. Max.
## -2.704774 -0.518825 -0.008539 0.002564 0.505061 3.509606prob <- 1/(1+exp(-pred))
summary(prob)## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 0.06269 0.37313 0.49787 0.50000 0.62365 0.97096해당 값들은 0과 1 사이에 있을뿐 기준값은 아직 명시되지 않았다
logistic_model 평가하기
summary(logistic_model)##
## Call:
## glm(formula = outcome ~ payment_inc_ratio + purpose_ + home_ +
## emp_len_ + borrower_score, family = "binomial", data = loan_data)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -2.51951 -1.06908 -0.05853 1.07421 2.15528
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 1.638092 0.073708 22.224 < 2e-16 ***
## payment_inc_ratio 0.079737 0.002487 32.058 < 2e-16 ***
## purpose_debt_consolidation 0.249373 0.027615 9.030 < 2e-16 ***
## purpose_home_improvement 0.407743 0.046615 8.747 < 2e-16 ***
## purpose_major_purchase 0.229628 0.053683 4.277 1.89e-05 ***
## purpose_medical 0.510479 0.086780 5.882 4.04e-09 ***
## purpose_other 0.620663 0.039436 15.738 < 2e-16 ***
## purpose_small_business 1.215261 0.063320 19.192 < 2e-16 ***
## home_OWN 0.048330 0.038036 1.271 0.204
## home_RENT 0.157320 0.021203 7.420 1.17e-13 ***
## emp_len_ > 1 Year -0.356731 0.052622 -6.779 1.21e-11 ***
## borrower_score -4.612638 0.083558 -55.203 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 62857 on 45341 degrees of freedom
## Residual deviance: 57515 on 45330 degrees of freedom
## AIC: 57539
##
## Number of Fisher Scoring iterations: 4p값은 변수의 중요성을 나타내는 상대적인 지표로 봐야한다.
이진 응답변수가 있는 로지스틱 회귀 모형은 RMSE나 R 제곱이 있을 수 없다.
대신 분류 문제에서 가장 일만적으로 사용되는 측정 지표들을 사용할 수 있다.
(혼동행렬 confusion matrix가 대표적인 예)
선형회귀에 적용되었던 많은 개념이 로지스틱 회귀에도 똑같이 이어진다.
단계적 회귀, 상호작용 항 도입, 스플라인 항 포함 등을 모두 사용할 수 있다.
교란변수가 변수 상관 관련된 문제도 동일하게 고려해야 함. mgcv 패키지 이용하기.
logistic_gam <- gam(outcome ~ s(payment_inc_ratio) + purpose_ +
home_ + emp_len_ + s(borrower_score),
data=loan_data, family='binomial')
logistic_gam##
## Family: binomial
## Link function: logit
##
## Formula:
## outcome ~ s(payment_inc_ratio) + purpose_ + home_ + emp_len_ +
## s(borrower_score)
##
## Estimated degrees of freedom:
## 7.66 4.17 total = 21.83
##
## UBRE score: 0.2681506terms <- predict(logistic_gam, type='terms')
partial_resid <- resid(logistic_model) + terms
df <- data.frame(payment_inc_ratio = loan_data[, 'payment_inc_ratio'],
terms = terms[, 's(payment_inc_ratio)'],
partial_resid = partial_resid[, 's(payment_inc_ratio)'])
ggplot(df, aes(x=payment_inc_ratio, y=partial_resid, solid=FALSE)) +
geom_point(shape=46, alpha=.4) +
geom_line(aes(x=payment_inc_ratio, y=terms),
color='red', alpha=.5, size=1.5) +
labs(y='Partial Residual')
위의 구름은 1의 응답 (연체 default), 아래 구름은 0의 응답 (대출 상환 paid off)이다.
결과변수가 이진형이기 때문에, 로지스틱 회귀의 잔차는 보통 이러한 형태를 띄게 된다.
분류 모델 평가하기
예측 모델링에서, 수많은 모델을 시도해보고 각각에 홀드아웃 표본(시험표본 test sample 또는 타당성검사 표본 validation sample)을 적용하고 성능을 평가하는 것은 아주 일반적이다.
분류 성능을 측정하는 가장 간단한 방법은 정확히 예측한 것들의 비율이 얼마인 지 보는 것이다.
가장 기본적인 컷오프 기준값은 0.5, 즉 50%이다. 확률이 0.5보다 크면 1, 아니면 0으로 분류한다.
혼동행렬
logistic_gam 모델의 혼동행렬(confusion matrix)를 계산한다.
pred <- predict(logistic_gam, newdata = loan_data)
pred_y <- as.numeric(pred > 0)
true_y <- as.numeric(loan_data$outcome == 'default')
true_pos <- (true_y==1) & (pred_y==1)
true_neg <- (true_y==0) & (pred_y==0)
false_pos <- (true_y==0) & (pred_y==1)
false_neg <- (true_y==1) & (pred_y==0)
conf_mat <- matrix(c(sum(true_pos), sum(false_pos),
sum(false_neg), sum(true_neg)), 2, 2)
colnames(conf_mat) <- c('Yhat = 1', 'Yhat = 0')
rownames(conf_mat) <- c('Y = 1', 'Y = 0')
conf_mat## Yhat = 1 Yhat = 0
## Y = 1 14293 8378
## Y = 0 8051 14620결과에서 열은 예측값, 행은 실제 결과를 의미한다.
여기서 가장 중요한 지표 중 하나는 거짓 양성 비율 false negative ratio(정밀도와 대응되는 개념)
예측 결과는 1이지만 실제로는 0일 가능성이 높은 상황이 된다.
회귀 클래스 문제
많은 문제에서 분류해야 할 클래스 간에 불균형이 존재하는 경우가 대부분이다.
보통은 데이터 수가 상대적으로 작은(희귀한) 클래스가 관심의 대상이 되기 때문에,
통상적으로 1로 지정하고, 반대로 수가 많은 경우를 0으로 지정한다.
일반적인 경우 1이 더 중요한 사건을 의미한다.
클래스를 쉽게 분리하기 어려운 경우에는,
가장 정확도가 높은 분류모델은 모든 것을 무조건 0으로 분류하는 모델일 수도 있다.
ex) 쇼핑몰의 경우 방문객의 0.1%만 실제 구매를 하고 나머지는 구매를 하지 않는다.
전반적인 정확도가 비록 떨어지더라도, 실제 구매자를 잘 골라내는 모델이 있다면 그 모델을 선호할 것이다.
정밀도 / 재현율 / 특이도
정밀도란 예측된 양성 결과의 정확도를 의미한다.
(정밀도) = (true positive) / (true positive + false positive)
재현율은 민감도라고 부르기도 하는데, 양성 결과를 예측하는 모델의 능력을 평가한다.
(재현율) = (true positive) / (true positive + false negative)
특이도는 음성 결과를 정확히 예측하는 능력을 측정한다.
(특이도) = (true negative) / (true negative + false positive)
# precision
conf_mat[1,1]/sum(conf_mat[,1])## [1] 0.6396796# recall
conf_mat[1,1]/sum(conf_mat[1,])## [1] 0.630453# specificity
conf_mat[2,2]/sum(conf_mat[2,])## [1] 0.6448767
ROC곡선
재현율과 특이도 사이에는 트레이드오프 관계 (시소 관계)가 있다.
이러한 트레이드오프 관계를 표현하기 위한 지표가
수신자 조작 특성 receiver operating characteristic ROC 곡선이다.
ROC곡선은 x축의 특이도에 대한 y축의 재현율을 표시한다.
idx <- order(-pred)
recall <- cumsum(true_y[idx]==1)/sum(true_y==1)
specificity <- (sum(true_y==0) - cumsum(true_y[idx]==0))/sum(true_y==0)
roc_df <- data.frame(recall = recall, specificity = specificity)
ggplot(roc_df, aes(x=specificity, y=recall)) +
geom_line(color='blue') +
scale_x_reverse(expand=c(0,0)) +
scale_y_continuous(expand=c(0,0)) +
geom_line(data=data.frame(x=(0:100)/100), aes(x=x, y=1-x),
linetype='dotted', color='red')
AUC (Area Underneath the Curve)
AUC는 간단히 말해 ROC 곡선의 아래쪽 면적을 의미함.
AUC 값이 높을 수록 더 좋은 분류기라고 할 수 있다.
ggplot(roc_df, aes(specificity)) +
geom_ribbon(aes(ymin=0, ymax=recall), fill='blue', alpha=.3) +
scale_x_reverse(expand=c(0, 0)) +
scale_y_continuous(expand=c(0, 0)) +
labs(y='recall') +
theme_bw()
sum(roc_df$recall[-1] * diff(1-roc_df$specificity))## [1] 0.6926232모델의 AUC 값은 약 0.7로 상당히 높다고 볼 수 있다.
분류기 성능을 비교하는 지표로 AUC 사용하면 단순히 정확도만을 사용하는 것보다는 나은 결과를 얻을 수 있다.
전체적인 정확도도 높이면서 중요한 1을 더 정확히 분류해야하는 트레이드오프를 얼마나 잘 처리하는지 평가할 수 있기 때문이다.
불균형 데이터 다루기
과소표본추출
앞서 다룬 대출 데이터와 같이 데이터 개수가 충분하다면,
다수의 데이터에 해당하는 클래스에서 과소표본추출(다운샘플링)을 해서
모델링할 때 0과 1의 데이터 개수에 균형을 맞출 수 있다.
다수의 클래스에 속한 데이터들 중에 중복된 데이터들이 많을 것이라는 사실에서 출발함.
소수 클래스의 데이터가 수만 개 정도 있다면 충분히 가능하다.
mean(full_train_set$outcome == 'default')## [1] 0.1889455전체 데이터에서 18% 정도만이 연체 상태이다.
모델을 학습하는 데 전체 데이터를 사용한다면 어떻게 될까?
full_model <- glm(outcome ~ payment_inc_ratio + purpose_ + home_ + emp_len_ + dti + revol_bal + revol_util, data=full_train_set, family='binomial')
pred <- predict(full_model)
mean(pred > 0)## [1] 0.003942094대출의 약 0.39% 정도 만이 연체 상태일 것이라고 예측하므로 기대되는 값보다 1/5 작은 수준이다.
모든 데이터를 동일하게 학습에 사용하다 보니,
대출을 갚는다는 예측이 갚지 않는다는 예측을 압도하는 결과를 보이는 것이다.
균형잡힌 데이터를 이용할 시 이 두개의 값이 비슷하게 나오는 것을 알 수 있다.
과잉표본추출과 상향/하향 가중치
다수 클래스를 과소표본추출하는 대신, 복원추출 방식(부트스트래핑)으로 희귀 클래스의 데이터를 과잉표본추출(업샘플링)해야 한다.
wt <- ifelse(full_train_set$outcome == 'default',
1/mean(full_train_set$outcome == 'default'), 1)
full_model <- glm(outcome ~ payment_inc_ratio + purpose_ + home_ + emp_len_ + dti + revol_bal + revol_util, data=full_train_set, weight=wt, family='binomial')## Warning in eval(family$initialize): non-integer #successes in a binomial glm!pred <- predict(full_model)
mean(pred>0)## [1] 0.5767208연체에 대한 가중치를 1/p로 두었다. (여기서 p는 연체의 확률값)
> 가중치를 적용하는 방식은 회귀 클래스를 업샘플링하거나 다수 클래스를 다운샘플링하는 방법으로 대체할 수 있다.
데이터 생성
기존의 데이터와 다른 데이터를 생성해서 좀 더 로버스트한 분류 규칙을 배울 수 있는 기회를 주고자 하는 것.
합성 소수 과잉표본 기법의 약자인 SMOTE(Synthetic Minority Oversampling Technique) 알고리즘은 업샘플링된 레코드와 비슷한 레코드를 찾고, 원래 라코드와 이웃 레코드의 랜덤 가중평균으로 새로운 합성 레코드를 만든다. 이후 각각의 예측변수에 대해 개별적으로 가중치를 생성한다.
비용 기반 분류
분류 규칙을 결정할 때 정확도나 AUC 만으로는 충분하지 않을 수 있다.
ex)
신규 대출에서 연체로 인해 예상되는 비용이 C라고 하고,
대출 상환을 통해 얻을 수 있는 수익을 R이라고 하자.
기대 수익 = P(Y=0) * R + P(Y=1) * C
여기서 대출 결과를 단순히 연체나 상환, 둘중에 하나로 결정하는 대신에
대출을 통해 얻을 수 있는 기대 수익이 있는지 없는지를 결정하는 것이 더 말이 된다.
연체확률이 더 높더라도 가치가 더 큰 대출을 선호하는 편이 나을 수도 있는 것이다.
예측 결과 분석
대출 데이터의 두 가지 예측 변수 borrower_score와 payment_inc_ratio를 사용해서 구한 4개의 서로 다른 모델에 대한 결정 규칙 알아보기
loan_tree <- rpart(outcome ~ borrower_score + payment_inc_ratio,
data=loan3000,
control = rpart.control(cp=.005))
lda_pred <- lda_df0[, c('borrower_score', 'payment_inc_ratio')]
lda_pred$method = 'LDA'
tree_pred <- data.frame(borrower_score = c(0.375, 0.375, 0.525, 0.525, 0.625, 0.625),
payment_inc_ratio = c(0, 9.732, 9.732, 8.772, 8.772, 20),
method = rep('Tree', 6))
glm0 <- glm(outcome ~ payment_inc_ratio + borrower_score,
data=loan3000, family='binomial')
y <- seq(from=0, to=20, length=100)
x <- (-glm0$coefficients[1] - glm0$coefficients[2]*y)/glm0$coefficients[3]
glm0_pred <- data.frame(borrower_score=x, payment_inc_ratio=y, method='Logistic')
gam1 <- gam(outcome ~ s(payment_inc_ratio) + s(borrower_score),
data=loan3000, family='binomial')
newdata = glm0_pred
gam_fun <- function(x){
rss <- sum(predict(gam1, newdata=data.frame(borrower_score=x, payment_inc_ratio=y))^2)
}
est_x <- nlminb(newdata$borrower_score, gam_fun )
gam1_pred <- data.frame(borrower_score=est_x$par, payment_inc_ratio=y, method="GAM")
loan_fits <- rbind(lda_pred,
tree_pred,
glm0_pred,
gam1_pred)
ggplot(data=loan_fits, aes(x=borrower_score, y=payment_inc_ratio, color=method, linetype=method)) +
geom_line(size=1.5) +
theme(legend.key.width = unit(2,"cm")) +
guides(linetype = guide_legend(override.aes = list(size = 1)))
선형판별분석(LDA), 로지스틱 선형회귀, GAM을 이용한 로지스틱 회귀, 트리 모델 4가지 사용.
최종적으로는 GAM 모델이 트리 모델과 다른 선형모형들을 서로 타협하는 결과를 보여준다.