로지스틱 회귀분석 in R

from http://economia.tistory.com/69

로지스틱 회귀분석이란? (Logistic regression)

로지스턱 회귀분석은 선형 회귀분석과 달리 결과 변수(종속 변수)가 범주형 데이터인 경우에 사용되는 기법입니다.

선형 회귀분석 모델에서는 설명 변수를 입력하면 수치형 결과를 얻게 되는데요, 신장, 시험 성적, 연간 소득 따위를 예측하고 싶을 때는 선형 회귀분석을 사용하면 됩니다.

하지만 내가 예측하고자 하는 것이 수치화 하기 힘든 변수, 예를 들면 어떤 고객이 부도를 낼 것인지의 여부, 타이타닉 호에서 살아남을 것인지의 여부, 어떤 학생의 최종 학력 (초등학고, 중학교, 대학교, 대학원 이상) 일 경우에는 선형 회귀분석 대신 로지스틱 회귀분석을 사용하게 됩니다.

회귀분석이기 때문에 기본적으로 지도 학습 (Supervised Learning) 으로 분류 되며, 다양한 분야에서 분류 및 예측을 위한 기법으로 사용되고 있습니다.  제가 근무하고 있는 라이엇 게임즈에서도 플레이어 이탈을 예측하기 위한 모델링 과정에서 로지스틱 회귀 분석을 비롯한 다양한 지도 학습 기법을 사용한 바 있습니다.

로지스틱 회귀분석 in R

R을 사용하면 아주 간단히 로지스틱 회귀분석 모델을 만들고, 이 모델을 이용해서 신규 데이터가 주어졌을 때 결과 값을 예측하는 작업을 할 수 있습니다.

1. 데이터셋 로드

UCLA에서 제공하는 데이터 셋을 불러와서 데이터의 구조를 확인하면 아래와 같은 결과를 얻을 수 있습니다.

1 2 3

data <- read.csv("http://www.ats.ucla.edu/stat/data/binary.csv") str(data) head(data)

1 2 3 4 5 6 7 8 9 10 11 12 13

'data.frame': 400 obs. of  4 variables: $ admit: int  0 1 1 1 0 1 1 0 1 0 ... $ gre  : int  380 660 800 640 520 760 560 400 540 700 ... $ gpa  : num  3.61 3.67 4 3.19 2.93 3 2.98 3.08 3.39 3.92 ... $ rank : int  3 3 1 4 4 2 1 2 3 2 ... admit gre  gpa rank 1     0 380 3.61    3 2     1 660 3.67    3 3     1 800 4.00    1 4     1 640 3.19    4 5     0 520 2.93    4 6     1 760 3.00    2

admit 변수는 지원자의 입학 여부를 나타내며, 0 혹은 1의 값을 가집니다. gre, gpa 변수는 수치형 변수인 시험 점수이고, rank는 지원자의 학교 등급을 나타내며 1~4의 값을 가지는 범주형 변수입니다.

2. 데이터 클린징 (DATA CLEANSING)

1번 단계에서 str 결과를 보면 rank 변수가 int로 설정되어 있는 것이 보이죠?
분석에 들어가기 전 rank 변수를 factor 변수로 변경해 줍니다.

1 2	data$rank <-as.factor(data$rank) str(data)

1 2 3 4 5

data.frame': 400 obs. of  4 variables: $ admit: int  0 1 1 1 0 1 1 0 1 0 ... $ gre  : int  380 660 800 640 520 760 560 400 540 700 ... $ gpa  : num  3.61 3.67 4 3.19 2.93 3 2.98 3.08 3.39 3.92 ... $ rank : Factor w/ 4 levels "1","2","3","4": 3 3 1 4 4 2 1 2 3 2 ...

rank 변수의 유형이 Factor로 바뀌었습니다. 이제 모델 fitting을 해 볼까요?

3.로지스틱 회귀분석 모델 작성 (MODEL FITTING)

로지스틱 회귀분석 모델에서는 “glm” 함수를 사용하는데요, 선형 회귀분석에서의 “lm”과 유사한 함수라고 보시면 됩니다.

“glm”함수를 사용해서 모델을 만든 후 “summary” 함수로 모델 결과를 확인합니다.

1 2 3 4

train<-data[1:200, ] test<-data[201:400, ] model <- glm(admit ~ gre + gpa + rank, data =train, family = "binomial") summary(model)

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26

Call: glm(formula = admit ~ gre + gpa + rank, family = "binomial",     data = train)   Deviance Residuals:     Min       1Q   Median       3Q      Max   -1.6154  -0.7872  -0.5325   0.8785   2.3242     Coefficients:              Estimate Std. Error z value Pr(>|z|)     (Intercept) -3.790221   1.750305  -2.165 0.030352 *   gre          0.003483   0.001718   2.027 0.042640 *   gpa          0.534007   0.519132   1.029 0.303642     rank2       -0.707194   0.476566  -1.484 0.137826     rank3       -1.803655   0.531352  -3.394 0.000688 *** rank4       -1.959514   0.629189  -3.114 0.001844 ** --- Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1   (Dispersion parameter for binomial family taken to be 1)       Null deviance: 237.18  on 199  degrees of freedom Residual deviance: 205.76  on 194  degrees of freedom AIC: 217.76   Number of Fisher Scoring iterations: 4

4. 모델 해석

표의 가장 우측 열의 p-value가 0.05보다 작은 변수는 gre, rank3, rank4 이네요.
이 세가지가 모두 통계적으로 유의하다고 해석 할 수 있습니다. rank3 의 p-value가 가장 낮은 것으로 보아 예측력이 좀더 강한 것으로 보이네요. Estimate 열의 부호가 마이너스(-) 이므로 랭크가 낮은 학교의 지원자의 합격 확률은 낮아진다고 볼 수 있겠네요.

(보다 정확히는 랭크 4 학교의 지원자는 그렇지 않은 지원자에 비해 합격 확률의 Log가 1.95만큼 감소하게 된다는 의미입니다.)

anova 분석을 해 볼까요?

1	anova(model, test="Chisq")

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16

Analysis of Deviance Table   Model: binomial, link: logit   Response: admit   Terms added sequentially (first to last)          Df Deviance Resid. Df Resid. Dev  Pr(>Chi)     NULL                   199     237.18               gre   1  12.8706       198     224.31 0.0003338 *** gpa   1   0.5607       197     223.75 0.4539719     rank  3  17.9894       194     205.76 0.0004421 *** --- Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Null deviance와 residual deviance의 차이는 Null 모델과 비교하여 해당 변수가 모델에 포함 되었을 때 모델의 성능이 얼마나 나아지는지를 보여주는 수치입니다. rank 가 포함되었을 때 resid. Dev 의 값이 역시 가장 많이 하락하네요.

로지스틱 회귀분석에도 선형 회귀 분석에서의   R2 와 유사한 개념이 존재하는데요, Mcfadden R2으로 모델 fit을 확인할 수 있습니다.
“pscl” 패키지의 pR2 함수를 사용하여 Mcfadden R2 를 간단히 알아 볼 수 있는데요.
R2값이 0.2인 것으로 보아 모델이 train 데이터 셋의 분산의 약 20% 정도 밖에 설명해 주지 못한다고 봐야 겠네요.

1 2 3 4

library(pscl) pR2(model)          llh      llhNull           G2     McFadden         r2ML         r2CU -102.8802728 -118.5906635   31.4207814    0.1324758    0.1453847    0.2093280

 

5. 모델 평가

로지스틱 회귀 모델의 평가를 하기 위해 ROC 커브를 그려볼까요?

ROC 커브 아래의 영역 AUC가 1에 가까울 수록 모델의 예측 정확도가 높다고 할 수 있있습니다. “ROCR” 패키지를 사용하면 로지스틱 회귀 모델의 ROC 곡선을 그릴 수 있습니다.

1 2 3 4 5

library(ROCR) p <- predict(model, newdata=test, type="response") pr <- prediction(p, test$admit) prf <- performance(pr, measure = "tpr", x.measure = "fpr") plot(prf)

아래 그래프는 train 데이터 셋으로 만든 모델을 test 데이터 셋에 적용했을 경우의 ROC 곡선입니다.  ROC 곡선이 reference line인 1번 라인에서 2번 라인으로 가까워 질 수록 모델의 정확도는 향상되는데요, 그래프가 reference line에 가까운 것으로 보아 모델의 성능이 그다지 뛰어나지는 않은 것 같네요.

 

ROC 그래프 아래 영역의 넓이인  AUC를 구해볼까요?

ROC 곡선이 2번 방향을 향할 수록 AUC는 넓어지는데요,
모든 case의 결과를 완벽하게 예측하는 모델인 경우 AUC는 1,  무작위로 예측한 모델과 다를 바가 없는 경우 (=reference line) AUC는 0.5의 값을 가지게 됩니다.

AUC를 판단하는 대략적인 기준은 아래와 같습니다.

excellent =  0.9~1
good = 0.8~0.9
fair = 0.7~0.8
poor = 0.6~0.7
fail = 0.5~0.6

예시 모델의 AUC를 볼까요? 0.632로 poor 한 수준으로 보입니다.

1 2 3 4

auc <- performance(pr, measure = "auc") auc <- auc@y.values[[1]] auc [1] 0.6323834

 

결과를 종합해 보면 모델은 유의하지만, 예측력은 그다지 높지는 않은 것으로 보이네요.
설명 변수를 추가하고, 유의하지 않은 변수를 제거하고 여러 번 모델을 다시 fitting 하는 작업이 필요할 것 같습니다.

작업 하고 계신 모델의 성능이 좋지 않다고 너무 낙담하지는 마세요.
모델의 성능을 높이기 위해서는 아래 세 가지 작업을 반복 수행해야 합니다.
1. 적당한 설명 변수를 입력 하는 것
2. 적절한 모델/알고리즘을 선택하는 것
3. 위 1~2번을 반복

대부분의 머신 러닝을 사용한 예측 모델은 위 과정을 반복하면서 더욱 정교해지니까요!