[R] 기계학습 과제...
반응형
학교 다닐 때.. 기계학습 수업 때 했던 과제 포스팅해봅니다. 지금와서 보면 많이 부족했네요..
필요패키지
set.seed(2017) # fix seed number
library(knitr)
library(cvTools)
library(MASS)
library(class)
library(mgcv)
library(rpart)
library(partykit)
library(randomForest)
library(e1071)
library(ggplot2)Abstract
ISLR 패키지에 내장되어 있는 OJ 데이터를 가지고 최적의 분류모형을 찾아보았다. Logistic regression, Linear discriminant analysis 등 모수적 가정을 하는 분류모형과 비모수적 방법인 k-Nearest Neighbor, Generalized additive model 등을 적용해보았으며, 더 나아가 이 밖에 실생활에서 널리 쓰일 수 있는 Tree 기반 모형인 Random forest도 적용하였고 1990년 중반에 분류모형으로써 센세이션을 일으킨 Support vector machine도 적용해보았다.
데이터의 분포 및 성질에 따라 각각의 방법을 적용함으로써 존재하는 장단점이 달라질 수 있지만, OJ 데이터에서는 본 내용에서 소개하는 분류모형 중 비교적 가장 최신의 기법인 Support vector machine이 가장 예측력이 좋게 나왔다.
1. Introduction
- 분류 문제는 관심이 대상이 되는 반응변수와 관련된 설명변수를 기반으로 사전 정의 된 그룹 또는 클래스에 객체를 할당해야할 때 발생한다.
- Statistical classification 이란 기계학습에서 지도학습(Supervised learning)의 일종으로 분류하고자 하는 반응변수와 관련된 설명변수를 가지고 있는 데이터로부터 분류 모형을 설정 및 식별하고 새로운 데이터가 주어졌을 때 반응변수의 클래스를 예측하는 통계적 방법이다.
- 일반적으로 전체 데이터의 일부분을 가지고 Training data로 정의하여 모형을 만들고 나머지를 Test data로 정의하여 모형의 예측력을 평가한다.
- 여기서는
OJ데이터를 가지고 다양한 분류 모델링을 적용하여 최적의 분류 예측모형을 찾는데 중점을 두겠다.
1. 1. Orange Juice(OJ) data
OJ데이터는 18개의 변수와 1070개의 관측치로 이루어진 데이터이다.- 이 데이터에는 특정 고객이 Citrus Hill 또는 Minute Maid Orange Juice를 구매한 1070건의 구매기록이 포함되어 있는 자료이며
ISLR패키지안에 내장되어 있다.## Purchase WeekofPurchase StoreID PriceCH PriceMM DiscCH DiscMM SpecialCH ## 1 CH 237 1 1.75 1.99 0.00 0.0 0 ## 2 CH 239 1 1.75 1.99 0.00 0.3 0 ## 3 CH 245 1 1.86 2.09 0.17 0.0 0 ## 4 MM 227 1 1.69 1.69 0.00 0.0 0 ## 5 CH 228 7 1.69 1.69 0.00 0.0 0 ## 6 CH 230 7 1.69 1.99 0.00 0.0 0 ## SpecialMM LoyalCH SalePriceMM SalePriceCH PriceDiff Store7 PctDiscMM ## 1 0 0.500000 1.99 1.75 0.24 No 0.000000 ## 2 1 0.600000 1.69 1.75 -0.06 No 0.150754 ## 3 0 0.680000 2.09 1.69 0.40 No 0.000000 ## 4 0 0.400000 1.69 1.69 0.00 No 0.000000 ## 5 0 0.956535 1.69 1.69 0.00 Yes 0.000000 ## 6 1 0.965228 1.99 1.69 0.30 Yes 0.000000 ## PctDiscCH ListPriceDiff STORE ## 1 0.000000 0.24 1 ## 2 0.000000 0.24 1 ## 3 0.091398 0.23 1 ## 4 0.000000 0.00 1 ## 5 0.000000 0.00 0 ## 6 0.000000 0.30 0- 출처 : Stine, Robert A., Foster, Dean P., Waterman, Richard P. Business Analysis Using Regression (1998). Published by Springer.
-
library(ISLR) data(OJ) head(OJ)
- 데이터의 구조 및 변수에 대한 설명은 다음과 같다.
## 'data.frame': 1070 obs. of 18 variables: ## $ Purchase : Factor w/ 2 levels "CH","MM": 1 1 1 2 1 1 1 1 1 1 ... ## $ WeekofPurchase: num 237 239 245 227 228 230 232 234 235 238 ... ## $ StoreID : num 1 1 1 1 7 7 7 7 7 7 ... ## $ PriceCH : num 1.75 1.75 1.86 1.69 1.69 1.69 1.69 1.75 1.75 1.75 ... ## $ PriceMM : num 1.99 1.99 2.09 1.69 1.69 1.99 1.99 1.99 1.99 1.99 ... ## $ DiscCH : num 0 0 0.17 0 0 0 0 0 0 0 ... ## $ DiscMM : num 0 0.3 0 0 0 0 0.4 0.4 0.4 0.4 ... ## $ SpecialCH : num 0 0 0 0 0 0 1 1 0 0 ... ## $ SpecialMM : num 0 1 0 0 0 1 1 0 0 0 ... ## $ LoyalCH : num 0.5 0.6 0.68 0.4 0.957 ... ## $ SalePriceMM : num 1.99 1.69 2.09 1.69 1.69 1.99 1.59 1.59 1.59 1.59 ... ## $ SalePriceCH : num 1.75 1.75 1.69 1.69 1.69 1.69 1.69 1.75 1.75 1.75 ... ## $ PriceDiff : num 0.24 -0.06 0.4 0 0 0.3 -0.1 -0.16 -0.16 -0.16 ... ## $ Store7 : Factor w/ 2 levels "No","Yes": 1 1 1 1 2 2 2 2 2 2 ... ## $ PctDiscMM : num 0 0.151 0 0 0 ... ## $ PctDiscCH : num 0 0 0.0914 0 0 ... ## $ ListPriceDiff : num 0.24 0.24 0.23 0 0 0.3 0.3 0.24 0.24 0.24 ... ## $ STORE : num 1 1 1 1 0 0 0 0 0 0 ...Variables Description Purchase A factor with levels CH and MM indicating whether the customer purchased Citrus Hill or Minute Maid Orange Juice WeekofPurchase Week of purchase StoreID Store ID PriceHH Price charged for CH PriceMM Price charged for MM DiscCH Discount offered for CH DiscMM Discount offered for MM SpecialCH Indicator of special on CH SpecialMM Indicator of special on MM LoyalCH Customer brand loyalty for CH SalePriceMM Sale price for MM : PriceMM - DiscMM SalePriceCH Sale price for CH : PriceCH - DiscHH PriceDiff Sale price of MM less sale price of CH : SalePriceMM - SalePriceCH Store7 A factor with levels No and Yes indication whether the sale is at Store 7 PctDiscMM Percentage discount for MM : DiscMM / PriceMM PctDiscCH Percentage discount for CH : DiscCH / PriceCH ListPriceDiff List price of MM less list price of CH : PriceMM - PriceCH STORE Which of 5 possible stores the sale occured at -
str(OJ) - 변수의 정의에 맞게끔 STORE 변수의 속성을 factor 형식으로 변경할 필요가 보이며, StoreID 변수와 STORE 변수는 똑같이 5개의 수준을 갖는 변수이므로 차후 모델링하는데 있어서 둘 중 STORE 변수 하나만을 이용할 것 이다.
- SpeicalCH, SpecialMM 또한 0 아니면 1의 값을 갖는 indicator로서 변수의 속성을 factor 형식으로 변경할 필요가 있다.
-
OJ$StoreID <- as.factor(OJ$StoreID) OJ$STORE <- as.factor(OJ$STORE) OJ$SpecialCH <- as.factor(OJ$SpecialCH) OJ$SpecialMM <- as.factor(OJ$SpecialMM) summary()함수를 이용하여 데이터의 요약치를 살펴본 결과는 다음과 같다.## Purchase WeekofPurchase StoreID PriceCH PriceMM ## CH:653 Min. :227.0 1:157 Min. :1.690 Min. :1.690 ## MM:417 1st Qu.:240.0 2:222 1st Qu.:1.790 1st Qu.:1.990 ## Median :257.0 3:196 Median :1.860 Median :2.090 ## Mean :254.4 4:139 Mean :1.867 Mean :2.085 ## 3rd Qu.:268.0 7:356 3rd Qu.:1.990 3rd Qu.:2.180 ## Max. :278.0 Max. :2.090 Max. :2.290 ## DiscCH DiscMM SpecialCH SpecialMM LoyalCH ## Min. :0.00000 Min. :0.0000 0:912 0:897 Min. :0.000011 ## 1st Qu.:0.00000 1st Qu.:0.0000 1:158 1:173 1st Qu.:0.325257 ## Median :0.00000 Median :0.0000 Median :0.600000 ## Mean :0.05186 Mean :0.1234 Mean :0.565782 ## 3rd Qu.:0.00000 3rd Qu.:0.2300 3rd Qu.:0.850873 ## Max. :0.50000 Max. :0.8000 Max. :0.999947 ## SalePriceMM SalePriceCH PriceDiff Store7 ## Min. :1.190 Min. :1.390 Min. :-0.6700 No :714 ## 1st Qu.:1.690 1st Qu.:1.750 1st Qu.: 0.0000 Yes:356 ## Median :2.090 Median :1.860 Median : 0.2300 ## Mean :1.962 Mean :1.816 Mean : 0.1465 ## 3rd Qu.:2.130 3rd Qu.:1.890 3rd Qu.: 0.3200 ## Max. :2.290 Max. :2.090 Max. : 0.6400 ## PctDiscMM PctDiscCH ListPriceDiff STORE ## Min. :0.0000 Min. :0.00000 Min. :0.000 0:356 ## 1st Qu.:0.0000 1st Qu.:0.00000 1st Qu.:0.140 1:157 ## Median :0.0000 Median :0.00000 Median :0.240 2:222 ## Mean :0.0593 Mean :0.02731 Mean :0.218 3:196 ## 3rd Qu.:0.1127 3rd Qu.:0.00000 3rd Qu.:0.300 4:139 ## Max. :0.4020 Max. :0.25269 Max. :0.440-
summary(OJ) - 변수를 자세히 살펴보면 각 설명변수가 서로 연관되어 있음을 알 수 있다.
- 따라서 앞으로 분석하는데 있어서 다음의 변수들을 사용할 것 이다.
- Purchase : Class of response
- WeekofPurchase : Week of purchase
- SpecialCH : Indicator of special on CH
- SpecialMM : Indicator of special on MM
- LoyalCH : Customer brand loyalty for CH, 구매자들의 CH에 대한 브랜드 충성도를 반영
- SalePriceMM : PriceMM - DiscMM, 실제 구매자들이 구매한 금액을 반영함으로써 해석상 용이 도모
- SalePriceCH : PriceCH - DiscHH, 실제 구매자들이 구매한 금액을 반영함으로써 해석상 용이 도모
- STORE : Which of 5 possible stores the sale occured at
- ListPriceDiff : PriceMM - PriceCH, 할인되지 않은 실제 가격을 반영하지 못했지만 그 가격 차이를 반영함으로써 실제 가격에 의한 원인요소를 어느정도 반영할 있을 것이라고 예상
1. 2. Methods and functions
Logistic regression
- Logistic regression은 다음과 같은 Logistic function을 이용하여 특정 클래스로 분류될 확률을 추정하는 회귀 모형이다.\[p(X) = \frac{e^{\beta_0 + \beta_1 X + \beta_2 X + ... + \beta_p X}}{1 + e^{\beta_0 + \beta_1 X + \beta_2 X + ... + \beta_p X}}\]
- R에서 기본으로 제공하는
stats패키지 안에 내장되어 있는glm()함수를 이용하여 모형을 적합시킬 수 있다. glm(formula, data, family = binomial(), ...)
Linear discriminant analysis (LDA)
- Logistic regression과 달리 분류하고자 하는 클래스가 3개 이상일 때도 적용시킬 수 있는 방법이다.
- 오히려 클래스의 구분이 명확한 경우 Logistic regression의 회귀계수 추정치의 불안정성을 보완해주는 대안으로서 널리 쓰이고 있는 방법이다.
- Linear discriminant function의 기본 형태는 다음과 같은 선형결합의 형태로 주어진다.
- \[\mathbf{d}=\mathbf{x}^{\text{T}}\mathbf{b} = b_1 x_1 + b_2 x_2 + ... + b_p x_p\]
- 각 클래스가 공통분산을 가지는 정규분포로부터 나온 것으로 가정하기 때문에 표본의 크기가 비교적 작으면서 반응변수의 분포가 근사적으로 정규분포를 따르게 되면 Logistic regression 보다 더 안정적인 방법이다.
MASS패키지 안에 내장되어 있는lda()함수를 이용하여 모형을 적합시킬 수 있다.lda(formulam, data, ...)DiscriMiner패키지 안에 내장되어 있는linDA()함수 역시 LDA 모형을 적합시킬 수 있는 함수이다.
Quadratic discriminant analysis (QDA)
- LDA와 달리 각 클래스에 대해 서로 같거나 다른 분산을 가정하고 추정하는 방법이다.
- 클래스를 결정하게 되는 경계는 1차식이 아닌 2차식의 형태를 띄게 된다.
- 클래스별 분산에 큰 차이가 있으면서 관측치가 충분할 때 LDA보다 예측력이 좋을 것이라고 예상할 수 있다.
- 이 또한
MASS패키지 안에 내장되어 있는qda()함수를 이용하여 모형을 적합시킬 수 있다. qda(formulam, data, ...)DiscriMiner패키지 안에 내장되어 있는quaDA()함수 역시 QDA 모형을 적합시킬 수 있는 함수이다.
K-nearest neighbors (kNN)
- 모수를 가정하지 않은 비모수적인 방법으로 기준으로부터 가까운 k개의 관측치들을 가지고 클래스를 분류를 하는 방법이다.
- Decision boundary가 비선형인 경우 LDA 및 Logistic regression 보다 예측력이 좋을 것이라고 예상할 수 있다.
- 하지만 어떤 설명변수가 중요한 것인지 판단할 수 있는 유의성 판단의 부재가 단점이다.
- 여기서 k는 Cross-Validation 을 통하여 오분류율을 가장 최소로 하는 k로 결정하게 된다.
class패키지 안에 내장되어 있는knn()함수를 이용하여 모형을 적합시킬 수 있다.knn(train, test, cl, k = 1, ...)
Generalized additive models (GAM)
- 각 설명변수들을 알려지지 않는(unknown) 비선형 함수 \(f_j\)를 가지고 변수들간 가법성(additivity)을 유지한 상태로 다중 선형회귀 분석을 하는 것과 비슷하다.
- GAM 역시도 Logistic regression(GLM)과 유사하게 link function을 이용할 수 있다.
- \[\log \bigg(\frac{p(X)}{1-p(X)}\bigg) = \beta_0 + f_1(x_{i1}) + f_2(x_{i2}) + ... + f_p(x_{ip})\]
- 비선형 함수를 고려하기 때문에 때때로 조금 더 정확한 예측력을 보일 수 있다.
mgcv패키지 안에 내장되어 있는gam()함수를 이용하여 모형을 적합시킬 수 있다.gam(formula, data, family = binomial(), ...)
Classification trees
- 결과에 대해 비교적 쉬운 해석을 제공하고 의사결정이 편리하다는 장점이 있다.
- Dummay variables이 존재해도 다루기가 쉽다.
- 너무 많은 의사결정나무 가지가 나오면 pruning을 거쳐 조금 단순화할 수 있다.
- 모형의 예측력을 평가하는 척도로 Gini Index, Cross-entropy 등을 이용할 수 있다.
- Gini index가 0에 가까울수록 순도가 높다고 판단
rpart패키지 안에 내장되어 있는rpart()함수를 이용하여 모형을 적합시킬 수 있다.rpart(formula, data, ...)
Bagging & Random forest
- Tree 모형의 high variance 및 예측력이 떨어지는 단점을 보완한 방법으로써 대표적으로 앙상블 기법인 Bagging과 Random forest를 이용할 수 있다.
- Bagging은 데이터로부터 복원추출을 통하여 \(B\)개의 Bootstrap sample을 가지고 \(B\)개의 Tree 모형을 만들어 모형의 예측력을 평가하는 방법이다.
- \(B\)개의 모형으로부터 나온 값의 평균 또는 최빈값을 가지고 예측력을 평가한다. 이 경우 일반 Tree 모형보다 분산을 작게 해준다.
- 또한, re-sampling 하는 과정에서 추출된 데이터를 가지고 modeling을 하고 추출되지 않은 데이터는 자동적으로 모형을 평가하는데 쓰이게 된다.(Out of Bag, 자동적으로 cross-validation)
- Random forest는 Bagging과 유사하나, Bootstrap sampling 할 때 과정이 조금 상이하다.
- 관측치 뿐만 아니라 변수 또한 임의로 뽑아서 modeling을 한다. (\(p\)개의 변수가 있을 때 일반적으로 \(\sqrt{p}\))
- 변수를 임의로 뽑아서 하는 이유는 서로 상관성이 높은 변수들이 있으면 모형의 분산이 커지기 때문이다.(Decorrelation)
- 따라서 Bagging에 비해서 예측력이 더 좋다.
- Random forest는 변수의 중요도(Variable Importance)를 확인할 수 있다.
randomForest패키지 안에 내장되어 있는randomForest()함수를 이용하여 모형을 적합시킬 수 있다.randomForest(formula, data, subset, mtry, importance = TRUE, ...)
Support vector machine
- 커널(Kenrel) 함수에 기반하고 초평면(Hyperplane) \(f(x)\)을 고려하여 초평면으로부터 support vector와의 간격 \(m\)값이 가장 큰 초평면을 decision boundary로 판단하는 방법이다. (Maximum margin hyperplane)
- 초평면의 선형결합 형태는 다음과 같다.
- \[f(x) = \beta_0 + \beta_1 x_1 + ... + \beta_p x_p = 0\]
- \(m\) 값은 제약조건이 있는 최적화 문제(Optimization problem)로 Linear한 형태의 경우 다음과 같이 정의할 수 있다.
- Maximize \(m\) such that \(\sum^{p}_{j=1} \beta_j^2 = 1\)
\[y_i(\beta_0 + \beta_1 x_{i1} + ... + \beta_p x_{ip}) \ge m(1-\epsilon_i)\]slack variable \(\epsilon_i \ge 0\), and \(\sum \epsilon_i \le C\)
- Maximize \(m\) such that \(\sum^{p}_{j=1} \beta_j^2 = 1\)
- 상당히 robust한 장점도 있으며 Linear한 형태 이 외에도 커널 함수를 접목하여 Nonlinear decision boundary 에도 적용시킬 수 있다.
e1071패키지 안에 내장되어 있는svm()함수를 이용하여 모형을 적합시킬 수 있다.svm(formula, data, kernel, cost, scale, ...)
1. 3. Evaluate performance
- 위에서 언급했듯 전체 데이터를 가지고 k-Folds Cross-Validation을 통해 Training data와 Test data로 분할하여 모델링을 진행할 것 이다. \((k=5)\)
- k-Folds Cross-Validation을 통해 과적합(overfitting)을 피하고 결과에 대한 일반화를 도모할 수 있다.
- Training data와 Test data를 분류하는 데이터의 index는
cvTools패키지 안에 내장되어 있는cvFolds()함수를 이용하여 할당 할 수 있다. - 각 Fold에 해당하는 값을 다음과 같이
index객체에 저장한다. -
index <- cvFolds( nrow(OJ), K = 5 ) - 따라서 Training data와 Test data를 다음과 같이 할당한다. (\(k = 1, 2, 3, 4, 5\))
-
# example OJ_train <- OJ[index$which == 1, ] OJ_test <- OJ[index$which != 1, ] - Training data를 가지고 분류 모델링을 한 후 Test data를 적용시켜 각 모형마다 오분류율을 구하고 이를 가지고 예측력을 평가한다.
2. Modeling
- 데이터에서 주어진 변수들 중 Purchase, WeekofPurchase, SpecialCH, SpecialMM, LoyalCH, SalePriceMM, SalepriceCH, STORE, ListPriceDiff 변수들을 활용하여 모형에 적용시킬 것 이다.
- 여러 가지 모형을 적용시킨 후 예측력이 가장 좋은 분석모형(오분류율이 가장 낮은 모형)을 최적의 모형으로 판단하겠다.
- Test MSE를 구하는 과정에서 5-folds Cross-Validation을 근거로 하였다.
- 분석 방법 안에서 최적의 모형을 찾는 것이 아닌 각 모형들간의 비교를 주된 목적으로 하므로 분석하는데 있어 사용된 함수들의 구체적인 옵션 argument들은 기본값으로 하였다.
- k-folds Cross-Validation을 활용하였기에 빠른 계산과 코드의 간결성을 고려하여
for문을 이용하여 각 함수를 적용하였다.
2. 1. Logistic regression
- k-Fold Cross-Validation에 의해 오분류율이 할당되는 NULL 벡터
logistic_error를 생성한다. -
# generate null vector assigned error rate logistic_error <- NULL - 모형식은 다음과 같이 정의할 수 있다.
\[\Pr(\mathbf X) = \frac{e^{\beta_0 + \beta_1 \text{WeekofPurchase} + \beta_2 \text{SpecialCH} + ... + \beta_8 \text{ListPriceDiff}}}{1 + e^{\beta_0 + \beta_1 \text{WeekofPurchase} + \beta_2 \text{SpecialCH} + ... + \beta_8 \text{ListPriceDiff}}}\]
for()문을 통해 Cross-Validation을 적용하면서glm()함수를 이용하여 로지스틱 회귀분석에 적합시키고 오분류율을 계산한다.- 반응변수가 어느 class에 할당될 지를 결정하는 확률을 구하기 위해
predict()함수를 이용하였고ifelse()함수를 통하여 예측된 class를 정의하였다.
for (k in 1:5) { # k-Fold STEP OJ_train <- OJ[index$which == k, ] OJ_test <- OJ[index$which != k, ] # Fit logistic regression logistic_fit <- glm(Purchase ~ WeekofPurchase + SpecialCH + SpecialMM + LoyalCH + SalePriceMM + SalePriceCH + STORE + ListPriceDiff, data = OJ_train, family = binomial()) # Classify response class logistic_pred <- ifelse( predict(logistic_fit, newdata = OJ_test, type = "response") >= 0.5, "MM", "CH") # Calculate mis-classification rate logistic_error[k] = mean( OJ_test$Purchase != logistic_pred ) }- 반응변수가 어느 class에 할당될 지를 결정하는 확률을 구하기 위해
- 마지막 fold를 이용하여 적합시킨 모형을 가지고
summary()함수를 통해 변수의 유의성을 살펴본 결과 브랜드의 충성도(LoyalCH)가 가장 중요한 요인으로 예상할 수 있다.## ## Call: ## glm(formula = Purchase ~ WeekofPurchase + SpecialCH + SpecialMM + ## LoyalCH + SalePriceMM + SalePriceCH + STORE + ListPriceDiff, ## family = binomial(), data = OJ_train) ## ## Deviance Residuals: ## Min 1Q Median 3Q Max ## -2.4806 -0.5556 -0.1855 0.4805 3.0382 ## ## Coefficients: ## Estimate Std. Error z value Pr(>|z|) ## (Intercept) 5.47181 4.14369 1.321 0.1867 ## WeekofPurchase -0.01622 0.01517 -1.069 0.2850 ## SpecialCH1 0.54619 0.77959 0.701 0.4835 ## SpecialMM1 0.07300 0.62429 0.117 0.9069 ## LoyalCH -6.52606 0.97330 -6.705 2.01e-11 *** ## SalePriceMM -2.50900 1.11507 -2.250 0.0244 * ## SalePriceCH 3.93520 2.22765 1.767 0.0773 . ## STORE1 0.25298 0.70270 0.360 0.7188 ## STORE2 0.58416 0.61092 0.956 0.3390 ## STORE3 0.40467 0.69164 0.585 0.5585 ## STORE4 -0.45714 0.82703 -0.553 0.5804 ## ListPriceDiff -4.06142 2.59399 -1.566 0.1174 ## --- ## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 ## ## (Dispersion parameter for binomial family taken to be 1) ## ## Null deviance: 287.56 on 213 degrees of freedom ## Residual deviance: 153.77 on 202 degrees of freedom ## AIC: 177.77 ## ## Number of Fisher Scoring iterations: 6 -
summary(logistic_fit) - 로지스틱 회귀분석을 적용했을 경우 Test MSE는 다음과 같다.
## [1] 0.1817757 -
( logistic_MSE <- mean(logistic_error) )
2. 2. Linear discriminant analysis
- k-Fold Cross-Validation에 의해 오분류율이 할당되는 NULL 벡터
lda_error를 생성한다. -
# generate null vector assigned error rate lda_error <- NULL for()문을 통해 Cross-Validation을 적용하면서lda()함수를 이용하여 LDA 모형에 적합시키고 오분류율을 계산한다.- 반응변수의 예측 class를 구하기 위해
predict()함수를 이용하였다.
for (k in 1:5) { # k-Fold STEP OJ_train <- OJ[index$which == k, ] OJ_test <- OJ[index$which != k, ] # Fit LDA model lda_fit <- lda(Purchase ~ WeekofPurchase + SpecialCH + SpecialMM + LoyalCH + SalePriceMM + SalePriceCH + STORE + ListPriceDiff, data = OJ_train) # Predicted response class lda_pred <- predict(lda_fit, newdata = OJ_test, type = "response")$class # Calculate mis-classification rate lda_error[k] = mean( OJ_test$Purchase != lda_pred ) }- 반응변수의 예측 class를 구하기 위해
- 마지막 fold를 이용하여 적합시킨 모형을 가지고 class별 각 변수의 평균과 계수추정치를 살펴보면 다음과 같다.
## Call: ## lda(Purchase ~ WeekofPurchase + SpecialCH + SpecialMM + LoyalCH + ## SalePriceMM + SalePriceCH + STORE + ListPriceDiff, data = OJ_train) ## ## Prior probabilities of groups: ## CH MM ## 0.6028037 0.3971963 ## ## Group means: ## WeekofPurchase SpecialCH1 SpecialMM1 LoyalCH SalePriceMM SalePriceCH ## CH 257.4961 0.20155039 0.1395349 0.7342251 2.028372 1.799225 ## MM 252.2118 0.09411765 0.2235294 0.3234680 1.898000 1.845765 ## STORE1 STORE2 STORE3 STORE4 ListPriceDiff ## CH 0.1085271 0.1782946 0.1085271 0.17054264 0.2376744 ## MM 0.1411765 0.2823529 0.2941176 0.05882353 0.1943529 ## ## Coefficients of linear discriminants: ## LD1 ## WeekofPurchase -0.00658402 ## SpecialCH1 0.17684793 ## SpecialMM1 0.15029155 ## LoyalCH -3.79436566 ## SalePriceMM -1.30844756 ## SalePriceCH 1.28178918 ## STORE1 0.10136868 ## STORE2 0.25239586 ## STORE3 0.30360433 ## STORE4 -0.18613585 ## ListPriceDiff -2.42348020 -
lda_fit - LDA를 적용했을 경우 Test MSE는 다음과 같다.
## [1] 0.1810748 -
( lda_MSE <- mean(lda_error) )
2. 3. Quadratic discriminant analysis
- k-Fold Cross-Validation에 의해 오분류율이 할당되는 NULL 벡터
qda_error를 생성한다. -
# generate null vector assigned error rate qda_error <- NULL for()문을 통해 Cross-Validation을 적용하면서lda()함수를 이용하여 LDA 모형에 적합시키고 오분류율을 계산한다.- 반응변수의 예측 class를 구하기 위해
predict()함수를 이용하였다.
for (k in 1:5) { # k-Fold STEP OJ_train <- OJ[index$which == k, ] OJ_test <- OJ[index$which != k, ] # Fit QDA model qda_fit <- qda(Purchase ~ WeekofPurchase + SpecialCH + SpecialMM + LoyalCH + SalePriceMM + SalePriceCH + STORE + ListPriceDiff, data = OJ_train) # Predicted response class qda_pred <- predict(qda_fit, newdata = OJ_test, type = "response")$class # Calculate mis-classification rate qda_error[k] = mean( OJ_test$Purchase != lda_pred ) }- 반응변수의 예측 class를 구하기 위해
- QDA를 적용했을 경우 Test MSE는 다음과 같다.
## [1] 0.2030374 -
( qda_MSE <- mean(qda_error) )
2. 4. k-Nearest neighbor
- k-Fold Cross-Validation에 의해 오분류율이 할당되는 행렬
knn_error를 생성한다. -
# generate null matrix assigned error rate knn_error <- matrix(NA, nrow = 150, ncol = 5) for()문을 통해 Cross-Validation을 적용하면서knn.cv()함수를 이용해 kNN 방법을 적용시킨다.- 여기서 k 값을 1~150까지 적용하여 Test MSE를 최소로 하는 최적의 k 값을 구한다.
-
for (i in 1:150) { for (k in 1:5) { # k-Fold STEP OJ_train <- OJ[index$which == k, ] OJ_test <- OJ[index$which != k, ] # Fit kNN & Predicted response class knn_fit <- knn.cv(train = OJ_train[, c(2, 8:12, 17:18)], cl = OJ_train$Purchase, k = i) # Calculate mis-classification rate knn_error[i, k] = mean( OJ_test$Purchase != knn_fit ) } } - 그 결과 최적의 k값은 101이다.
- 최적의 k값을 가지고 kNN을 적용했을 경우 Test MSE는 다음과 같다.
( knncv_MSE = mean( OJ_test$Purchase != knn_fit.opt ) )## [1] 0.3878505 -
knn_fit.opt <- knn.cv(train = OJ_train[, c(2, 8:12, 17:18)], cl = OJ_train$Purchase, k = which.min( rowMeans(knn_error) ))
2. 5. Generalized additive models
- k-Fold Cross-Validation에 의해 오분류율이 할당되는 NULL 벡터
gam_error를 생성한다. -
# generate null vector assigned error rate gam_error <- NULL for()문을 통해 Cross-Validation을 적용하면서gam()함수를 이용하여 GAM에 적합시키고 오분류율을 계산한다.- \[\log \bigg(\frac{p(\textbf{X})}{1-p(\textbf{X})}\bigg) = \beta_0 + f_1(\text{WeekofPurchase}_{i}) + f_2(\text{LoyalCH}_{i}) + ... + f_5(\text{ListPriceDiff}_{i}) + \text{SpecialCH}_i + \text{SpeicalMM}_i + \text{STORE}_i\]
- 여기서 factor형인 변수들을 제외한 나머지 변수들은
s()함수를 이용하여 해당 함수를 smooth function \(f_p(x)\)화 시키겠다.s()함수의 나머지 arguments는 기본값으로 하였다.- Logistic regression과 마찬가지로 반응변수가 어느 class에 할당될 지를 결정하는 확률을 구하기 위해
predict()함수를 이용하였고ifelse()함수를 통하여 예측된 class를 정의하였다.
for (k in 1:5) { # k-Fold STEP OJ_train <- OJ[index$which == k, ] OJ_test <- OJ[index$which != k, ] # Fit generalized additive model gam_fit <- gam(Purchase ~ s(WeekofPurchase) + SpecialCH + SpecialMM + s(LoyalCH) + s(SalePriceMM) + s(SalePriceCH) + STORE + s(ListPriceDiff), data = OJ_train, family = binomial()) # Classify response class gam_pred <- ifelse( predict(gam_fit, newdata = OJ_test, type = "response") >= 0.5, "MM", "CH") # Calculate mis-classification rate gam_error[k] = mean( OJ_test$Purchase != gam_pred ) } - Logistic regression과 마찬가지로 반응변수가 어느 class에 할당될 지를 결정하는 확률을 구하기 위해
- 마지막 fold를 이용하여 적합시킨 모형을 가지고
summary()함수를 통해 parametric term과 smooth function term의 유의성을 확인할 수 있다.## ## Family: binomial ## Link function: logit ## ## Formula: ## Purchase ~ s(WeekofPurchase) + SpecialCH + SpecialMM + s(LoyalCH) + ## s(SalePriceMM) + s(SalePriceCH) + STORE + s(ListPriceDiff) ## ## Parametric coefficients: ## Estimate Std. Error z value Pr(>|z|) ## (Intercept) -1.2700 0.4845 -2.621 0.00876 ** ## SpecialCH1 0.9424 0.8669 1.087 0.27695 ## SpecialMM1 0.1002 0.6896 0.145 0.88449 ## STORE1 0.5238 0.7688 0.681 0.49562 ## STORE2 0.8531 0.6610 1.291 0.19681 ## STORE3 0.4079 0.7626 0.535 0.59278 ## STORE4 -0.9537 0.9123 -1.045 0.29583 ## --- ## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 ## ## Approximate significance of smooth terms: ## edf Ref.df Chi.sq p-value ## s(WeekofPurchase) 1.000 1.000 3.842 0.0500 * ## s(LoyalCH) 4.707 5.721 48.770 9.18e-09 *** ## s(SalePriceMM) 2.761 3.471 7.664 0.0796 . ## s(SalePriceCH) 2.974 3.694 7.470 0.0917 . ## s(ListPriceDiff) 1.002 1.004 2.184 0.1396 ## --- ## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 ## ## R-sq.(adj) = 0.56 Deviance explained = 52.9% ## UBRE = -0.18561 Scale est. = 1 n = 214 -
summary(gam_fit) - GAM을 적용했을 경우 Test MSE는 다음과 같다.
## [1] 0.1880841 -
( gam_MSE <- mean(gam_error) )
2. 6. Classification trees
- k-Fold Cross-Validation에 의해 오분류율이 할당되는 NULL 벡터
tree_error를 생성한다. -
# generate null vector assigned error rate tree_error <- NULL for()문을 통해 Cross-Validation을 적용하면서rpart()함수를 이용하여 Tree 모형에 적합시키고 오분류율을 계산한다.- 반응변수의 예측 class를 구하기 위해
predict()함수를 이용하였다.
for (k in 1:5) { # k-Fold STEP OJ_train <- OJ[index$which == k, ] OJ_test <- OJ[index$which != k, ] # Fit tree model tree_fit <- rpart(Purchase ~ WeekofPurchase + SpecialCH + SpecialMM + LoyalCH + SalePriceMM + SalePriceCH + STORE + ListPriceDiff, data = OJ_train) # Classify response class tree_pred <- predict(tree_fit, newdata = OJ_test, type = "class") # Calculate mis-classification rate tree_error[k] = mean( OJ_test$Purchase != tree_pred ) }- 반응변수의 예측 class를 구하기 위해
- 더 나은 시각화를 위해
partykit패키지 안에 내장되어 있는as.party()함수를 통하여 적합된 모형 객체의 class를 적절히 변형시킨 후plot()함수를 이용하여 Tree를 출력해 본 결과는 다음과 같다.
-
plot( as.party(tree_fit) ) - Tree 모형을 적용했을 경우 Test MSE는 다음과 같다.
## [1] 0.203271 -
( tree_MSE <- mean(tree_error) )
2. 7. Random forest
- k-Fold Cross-Validation에 의해 오분류율이 할당되는 NULL 벡터
rf_error를 생성한다. -
# generate null vector assigned error rate rf_error <- NULL for()문을 통해 Cross-Validation을 적용하면서randomForest()함수를 이용하여 모형에 적합시키고 오분류율을 계산한다.- 여기서 뽑히는 변수의 갯수를 설정하는
mtryargument는 \(\sqrt8 \approx 3\), random forest를 구성하는 나무의 갯수를 설정하는ntreeargument는 기본값으로 설정하였다.- 반응변수의 예측 class를 구하기 위해
predict()함수를 이용하였다.
for (k in 1:5) { # k-Fold STEP OJ_train <- OJ[index$which == k, ] OJ_test <- OJ[index$which != k, ] # Fit random forest model rf_fit <- randomForest(Purchase ~ WeekofPurchase + SpecialCH + SpecialMM + LoyalCH + SalePriceMM + SalePriceCH + STORE + ListPriceDiff, data = OJ_train, mtry = sqrt(8), importance = TRUE) # Classify response class rf_pred <- predict(rf_fit, newdata = OJ_test) # Calculate mis-classification rate rf_error[k] = mean( OJ_test$Purchase != rf_pred ) } - 반응변수의 예측 class를 구하기 위해
varImpPlot함수를 이용하여 변수의 중요도를 살펴보면 다음과 같다.
-
varImpPlot(rf_fit, main = "Variable importance as measured by Random forest") - Random forest 모형을 적용했을 경우 Test MSE는 다음과 같다.
## [1] 0.1850467 -
( rf_MSE <- mean(rf_error) )
2. 8. Support vector machine
- k-Fold Cross-Validation에 의해 오분류율이 할당되는 NULL 벡터
svm_error를 생성한다. -
# generate null vector assigned error rate svm_error <- NULL for()문을 통해 Cross-Validation을 적용하면서svm()함수를 이용하여 모형에 적합시키고 오분류율을 계산한다.- 여기서 커널 함수는 radial basis kernel을 이용하였으며 (
kernel = "radial"), 그에 따른cost나gammaargument는 기본값으로 설정하였다.for (k in 1:5) { # k-Fold STEP OJ_train <- OJ[index$which == k, ] OJ_test <- OJ[index$which != k, ] # Fit SVM model svm_fit <- svm(Purchase ~ WeekofPurchase + SpecialCH + SpecialMM + LoyalCH + SalePriceMM + SalePriceCH + STORE + ListPriceDiff, data = OJ_train, kernel = "radial") # Classify response class svm_pred <- predict(svm_fit, newdata = OJ_test) # Calculate mis-classification rate svm_error[k] = mean( OJ_test$Purchase != svm_pred ) } - \[K(\mathbf u, \mathbf v) = \exp \bigg(-\frac{||\mathbf u - \mathbf v||^2}{2\sigma^2} \bigg)\]
- Support vector machine 모형을 적용했을 경우 Test MSE는 다음과 같다.
## [1] 0.1778037 -
( svm_MSE <- mean(svm_error) )
3. Conclusion
- Logistic regression 이나 Random forest 모형에서 보였듯이 LoyalCH 가 구매하고자 하는 브랜드에 가장 영향을 많이 미쳤으며, 판매가격 (SalePriceMM) 또한 적지 않은 영향을 미쳤다.
- k-folds Cross-Validation를 기반으로 위에서 언급한 각 분류모형들을 적용하여 Test MSE를 비교하여 보았다.
3. 1. Comparison Test MSE
- 각 분류모형에서 구한 MSE는 다음과 같다.
## Method MSE ## 1 Logistic 0.1818 ## 2 LDA 0.1811 ## 3 QDA 0.2030 ## 4 kNNCV 0.3879 ## 5 GAM 0.1881 ## 6 Tree 0.2033 ## 7 Random Forest 0.1850 ## 8 SVM 0.1778 -
mse.value <- round( c(logistic_MSE, lda_MSE, qda_MSE, knncv_MSE, gam_MSE, tree_MSE, rf_MSE, svm_MSE), 4) MSE <- data.frame( Method = c("Logistic", "LDA", "QDA", "kNNCV", "GAM", "Tree", "Random Forest", "SVM"), MSE = mse.value ) MSE - MSE를 비교한 그림은 다음과 같다. 값을 시각화하는데 쓰인 함수는
ggplot2패키지 안에 내장되어 있는ggplot()함수를 이용하였다.p + geom_text(aes( x = 4, y = mse.value[1]-0.1, label = mse.value[1], size = 5, vjust = 0.5)) + geom_text(aes( x = 3, y = mse.value[2]-0.1, label = mse.value[2], size = 5, vjust = 0.5)) + geom_text(aes( x = 5, y = mse.value[3]-0.1, label = mse.value[3], size = 5, vjust = 0.5)) + geom_text(aes( x = 2, y = mse.value[4]-0.1, label = mse.value[4], size = 5, vjust = 0.5)) + geom_text(aes( x = 1, y = mse.value[5]-0.1, label = mse.value[5], size = 5, vjust = 0.5)) + geom_text(aes( x = 8, y = mse.value[6]-0.1, label = mse.value[6], size = 5, vjust = 0.5)) + geom_text(aes( x = 6, y = mse.value[7]-0.1, label = mse.value[7], size = 5, vjust = 0.5)) + geom_text(aes( x = 7, y = mse.value[8]-0.1, label = mse.value[8], size = 5, vjust = 0.5))
-
theme_set(theme_bw()) p <- ggplot(MSE, aes(Method, MSE, fill = Method)) + geom_bar(stat = "identity") + labs(title = "Comparison") p <- p + theme(legend.position = "none", plot.title = element_text( hjust = 0.5 ), axis.text.x = element_text(size = 13, angle = 30), axis.text.y = element_text(size = 10)) - k-Nearest Neighbor 방법이 Test MSE가 가장 높게 나왔으며 Support vector machine을 이용한 방법이 Test MSE가 가장 낮게 출력되었다.
3. 2. Summary
- 가격에 따른 구매차이보다는 브랜드 충성도에 따른 구매 차이가 더 큰 것으로 보인다.
- 몇 가지 가정을 전제로하는 로지스틱 회귀모형이나 LDA 모형의 경우에는 비슷한 performance를 보이며, 가장 높은 오분류율을 나타낸 kNN의 경우 이 데이터의 반응변수를 예측하는 데 있어서 다소 부적절하다고 생각된다.
- 또한 예상되어진대로 Tree 모형보다 보완점을 개선한 Random forest 모형이 더 좋은 예측력을 나타내고 있다.
- Smooth function을 기반으로 하는 비모수적인 관점인 GAM 모형 역시 Random forest와 비슷한 performance를 나타내고 있다.
- 하지만 결과적으로 Test MSE 관점에서 보면 Support vector machine을 이용한 분류모형이 가장 최적의 모형으로 판단되어진다.
반응형
'ETC' 카테고리의 다른 글
| [R] R을 활용한 테일러 전개 (0) | 2017.07.04 |
|---|---|
| [R] R을 활용한 미분과 적분 계산 (0) | 2017.07.04 |
TAGS.