IE 582 - Homework 4
ilaydacelenk
1 Introduction
The aim of this homework is to compare the performance of different machine learning algorithms such as Penalized Regression Approaches, Decision Trees, Support Vector Machines and Tree-Based Ensembles. Four datasets are selected for this purpose.
2 Libraries
We will use library caret as recommended in class with many other libraries.
#importing the libraries
library(caret)
library(dplyr)
library(corrplot)
library(data.table)
library(glmnet)
library(rpart)
library(TunePareto)
library(randomForest)
library(tidyverse)
library(ggplot2)3 Data Sets
Four dataset from UCI and Kaggle are selected. Two of them have separate train and test data, and for the others I will generate train and test by splitting data with 20% test, using createDataPartition function of the package caret. All the test and train samples are larger than 200 and all the data sets have more than 20 features. Characteristics of the selected datasets are given as follows:
- Optical Recognition of Handwritten Digits
- Multi-Class Classification
- More than 50 features
- Only Numerical Features
- Mushroom Types
- Binary Classification
- Only Categorical Features
- Telco Customer Churn
- Binary Classification
- Class Imbalance
- Both Categorical and Numerical Features
- House Prices
- Regression
- More than 50 Features
- Both Categorical and Numerical Features
After reading and pre-processing the datasets, I will implement different approaches in order to make predictions and then I will compare their performances.
4 Approaches
The following methods will be implemented in each dataset and their performances will be compared:
- Penalized Regression Approaches (PRA)
- Decision Trees (DT)
- Random Forests (RF)
- Stochastic gradient Boosting (SGB)
- Support Vector Machines (SVM)
5 Optical Recognition of Handwritten Digits
digits_test_raw <- read.csv("https://raw.githubusercontent.com/BU-IE-582/fall20-ilaydacelenk/master/files/HW4_data/digits_test", header=FALSE)
digits_train_raw <- read.csv("https://raw.githubusercontent.com/BU-IE-582/fall20-ilaydacelenk/master/files/HW4_data/digits_train", header=FALSE)
#read.csv("~/Desktop/projects/fall20-ilaydacelenk/files/HW4_data/digits_test", header=FALSE)
digits_train <- digits_train_raw %>% select(-V1, -V40, -V33)
digits_test <- digits_test_raw %>% select(-V1, -V40, -V33)
X_digits_train <- digits_train %>% select(-V65) %>% scale()
y_digits_train <- digits_train %>% mutate(V65=as.factor(V65)) %>% select(V65)
X_digits_test <- digits_test %>% select(-V65) %>% scale()
y_digits_test <- digits_test %>% mutate(V65=as.factor(V65)) %>% select(V65)
#corrplot(cor(X_digits_train, method = c("pearson")))
#after scaling, we can bind them again
digits_train <- cbind(as.data.table(X_digits_train), y_digits_train)
digits_test <- cbind(as.data.table(X_digits_test), y_digits_test)Optical Recognition of Handwritten Digits dataset is used as a multi-class classification problem with 10 classes and it has more than 50 features. Here, the aim is to distinguish different handwritten digits. For this data set, there are separate train and test samples. The data has 65 columns, where the first 64 columns are the independent variables and the last columns is the dependent variable. Train data has 3823 instances and test data has 1797 instances, both are greater than 200 as desired. Train data is collected from 30 people and test data is collected from 13 different people. All the dependent variables are discrete numeric in the range 0-16 and the dependent variable represents the digit class code which is 0-9. The predictors are created using normalized bitmaps.
When we plot the correlation matrix, we see that the corresponding values to column V1 and V40 are “?” on the matrix. Then we call summary function on X_digits_train and see that V1 and V40 columns are zero for all the instances. Therefore, we remove these columns.
33rd column of the test data is a zero column, therefore when it is scaled, it becomes NA. In order to apply the models, I will also remove it.
5.1 Penalized Regression Approach
set.seed(50)
digits_model_pra<-cv.glmnet(as.matrix(X_digits_train), digits_train$V65, type.measure = "class", family="multinomial", nfolds=10)
digits_pred_pra <- predict(digits_model_pra, newx=X_digits_test ,s=c("lambda.min"), type = "class")
digits_cm_pra_test <- confusionMatrix(factor(digits_pred_pra), digits_test$V65)
digits_cm_pra_train <- confusionMatrix(factor(predict(digits_model_pra, newx=X_digits_train ,s=c("lambda.min"), type = "class")), digits_train$V65)
print(digits_model_pra)FALSE
FALSE Call: cv.glmnet(x = as.matrix(X_digits_train), y = digits_train$V65, type.measure = "class", nfolds = 10, family = "multinomial")
FALSE
FALSE Measure: Misclassification Error
FALSE
FALSE Lambda Measure SE Nonzero
FALSE min 0.0005698 0.02668 0.002095 26
FALSE 1se 0.0008266 0.02877 0.002063 25plot(digits_model_pra)
Lambda value that minimizes the error is 5.697763710^{-4} with error rate 2.67%, the error rate on the test data is 4.4% which is higher.
Accuracy of the model on test data is 95.6%. On the other hand, accuracy on the train data is 98.59%. Those values are very close, therefore there is no significant overfitting.
Error rates on test and train data are 4.4% and 1.41%.
As a result, the model performs better on the train data which may indicate overfitting.
5.2 Decision Tree
set.seed(50)
complexity=c(0.01,0.03,0.05)
min_number_of_observations=c(10,15,25)
#for balanced folds: stratified =TRUE
index <- generateCVRuns(digits_train$V65, ntimes=1, nfold=10, stratified =TRUE)
acc <- 0
best_complexity <- 0
best_min_number_of_observations <- 0
cv_summary=data.table()
for(i in 1:10){
for(c in complexity){
for(m in min_number_of_observations){
test_index <- index$`Run 1`[i]
train_cv <- digits_train[-test_index[[1]],]
test_cv <- digits_train[test_index[[1]],]
fitControl <- rpart.control(cp = c, minbucket = m)
digits_model_DT=rpart(V65~., digits_train, control=fitControl)
pred <- predict(digits_model_DT,data.frame(test_cv),type="class")
cm <- confusionMatrix(pred, test_cv$V65)
acc_temp <- cm$overall[1:1]
cv_summary <- rbind(cv_summary,data.table(Fold=i, Complexity=c, Min_Observations_per_Leaf=m, Accuracy=acc_temp))
if(acc_temp>acc){
acc<-acc_temp
best_complexity <- c
best_min_number_of_observations <- m
}
}
}
}
cv_summary <- cv_summary %>% group_by(Complexity, Min_Observations_per_Leaf) %>% summarise(mean_accuracy = mean(Accuracy))
print(cv_summary)FALSE # A tibble: 9 x 3
FALSE # Groups: Complexity [3]
FALSE Complexity Min_Observations_per_Leaf mean_accuracy
FALSE <dbl> <dbl> <dbl>
FALSE 1 0.01 10 0.769
FALSE 2 0.01 15 0.769
FALSE 3 0.01 25 0.769
FALSE 4 0.03 10 0.679
FALSE 5 0.03 15 0.679
FALSE 6 0.03 25 0.679
FALSE 7 0.05 10 0.612
FALSE 8 0.05 15 0.612
FALSE 9 0.05 25 0.612fitControl <- rpart.control(cp = best_complexity, minbucket = best_min_number_of_observations)
digits_model_DT <- rpart(V65~., digits_train, control=fitControl)
digits_pred_DT <- predict(digits_model_DT, data.frame(X_digits_test), type="class")
digits_cm_DT_test <- confusionMatrix(factor(digits_pred_DT), digits_test$V65)
digits_cm_DT_train <- confusionMatrix(factor(predict(digits_model_DT, data.frame(X_digits_train), type="class")), digits_train$V65)According to the cross validation result, cross validation error rate is 23.12, the error rate on the test data is 25.77% which is higher.
The best complexity parameter is 0.01 and the best minimal number of observations per tree leaf is 10 since it gives the highest accuracy of the model on cross validation test data with 78.53%.
Now I can fit the final model and make predictions. Accuracy of the final model on test data is 74.23%. On the other hand, for the train data the accuracy is 76.88%.
Error rates on test and train data are 25.77% and 23.12%.
As a result, the model performs slightly better on the train data.
5.3 Random Forest
set.seed(50)
m=c(4,6,8)
#for balanced folds: stratified =TRUE
# ntree is J=500 given in hw
J=500
# the minimal number of observations per tree leaf is nodesize 5 given in hw
min_size_of_terminal_nodes=5
n_repeats=2
n_folds=10
fitControl=trainControl(method = "repeatedcv",
number = n_folds,
repeats = n_repeats)
rf_grid=expand.grid(mtry = m)
digits_model_RF=train(V65 ~ ., data = data.frame(digits_train), method = "rf", trControl = fitControl, ntree=J, nodesize = min_size_of_terminal_nodes, tuneGrid = rf_grid)
digits_pred_RF <- predict(digits_model_RF, newdata=X_digits_test)
digits_cm_RF_test <- confusionMatrix(factor(digits_pred_RF), digits_test$V65)
digits_cm_RF_train <- confusionMatrix(factor(predict(digits_model_RF, newdata=X_digits_train)), digits_train$V65)The best mtry parameter is found as 8. I will fit the final model and make predictions setting mtry as 8.
Accuracy of the final model on test data is 97.5%. On the other hand, for the train data the accuracy is 99.95%.
Error rates on test and train data are 2.5% and 0.05%.
As a result, the model accuracy is higher on the train data, almost 100%, which may indicate overfitting.
5.4 Stochastic Gradient Boosting
set.seed(50)
n_folds <- 10
fitControl <- trainControl(method = "cv", number = n_folds)
gbmGrid_forever <- expand.grid(interaction.depth = c(1, 5, 7), n.trees = c(50,100,200), shrinkage = c(0.05, 0.1, 0.2), n.minobsinnode = c(10))
print(gbmGrid_forever)FALSE interaction.depth n.trees shrinkage n.minobsinnode
FALSE 1 1 50 0.05 10
FALSE 2 5 50 0.05 10
FALSE 3 7 50 0.05 10
FALSE 4 1 100 0.05 10
FALSE 5 5 100 0.05 10
FALSE 6 7 100 0.05 10
FALSE 7 1 200 0.05 10
FALSE 8 5 200 0.05 10
FALSE 9 7 200 0.05 10
FALSE 10 1 50 0.10 10
FALSE 11 5 50 0.10 10
FALSE 12 7 50 0.10 10
FALSE 13 1 100 0.10 10
FALSE 14 5 100 0.10 10
FALSE 15 7 100 0.10 10
FALSE 16 1 200 0.10 10
FALSE 17 5 200 0.10 10
FALSE 18 7 200 0.10 10
FALSE 19 1 50 0.20 10
FALSE 20 5 50 0.20 10
FALSE 21 7 50 0.20 10
FALSE 22 1 100 0.20 10
FALSE 23 5 100 0.20 10
FALSE 24 7 100 0.20 10
FALSE 25 1 200 0.20 10
FALSE 26 5 200 0.20 10
FALSE 27 7 200 0.20 10gbmGrid <- expand.grid(interaction.depth = c(5), n.trees = c(100), shrinkage = c(0.1), n.minobsinnode = c(10))
#I got the error "V57 has no variation." Therefore I removed it from the train set
digits_train2 <- digits_train %>% select(-V57)
digits_model_sgb=train(V65 ~ ., data = data.frame(digits_train2), method = "gbm", trControl = fitControl, tuneGrid = gbmGrid_forever, verbose=F)
digits_pred_sgb <- predict(digits_model_sgb, newdata=X_digits_test) %>% as.matrix()
digits_cm_sgb_test <- confusionMatrix(factor(digits_pred_sgb), digits_test$V65)
digits_cm_sgb_train <- confusionMatrix(factor(predict(digits_model_sgb, newdata=X_digits_train)), digits_train$V65)
print(digits_model_sgb[["results"]])FALSE shrinkage interaction.depth n.minobsinnode n.trees Accuracy Kappa
FALSE 1 0.05 1 10 50 0.8624417 0.8471499
FALSE 10 0.10 1 10 50 0.9158068 0.9064481
FALSE 19 0.20 1 10 50 0.9463990 0.9404410
FALSE 4 0.05 5 10 50 0.9529354 0.9477043
FALSE 13 0.10 5 10 50 0.9683538 0.9648364
FALSE 22 0.20 5 10 50 0.9714959 0.9683276
FALSE 7 0.05 7 10 50 0.9620854 0.9578712
FALSE 16 0.10 7 10 50 0.9691440 0.9657144
FALSE 25 0.20 7 10 50 0.9725355 0.9694827
FALSE 2 0.05 1 10 100 0.9131842 0.9035336
FALSE 11 0.10 1 10 100 0.9445720 0.9384113
FALSE 20 0.20 1 10 100 0.9623458 0.9581605
FALSE 5 0.05 5 10 100 0.9680948 0.9645484
FALSE 14 0.10 5 10 100 0.9730584 0.9700638
FALSE 23 0.20 5 10 100 0.9764574 0.9738405
FALSE 8 0.05 7 10 100 0.9717502 0.9686101
FALSE 17 0.10 7 10 100 0.9761977 0.9735521
FALSE 26 0.20 7 10 100 0.9793384 0.9770418
FALSE 3 0.05 1 10 200 0.9471919 0.9413223
FALSE 12 0.10 1 10 200 0.9628619 0.9587339
FALSE 21 0.20 1 10 200 0.9649507 0.9610550
FALSE 6 0.05 5 10 200 0.9759380 0.9732634
FALSE 15 0.10 5 10 200 0.9793364 0.9770395
FALSE 24 0.20 5 10 200 0.9777657 0.9752942
FALSE 9 0.05 7 10 200 0.9751540 0.9723923
FALSE 18 0.10 7 10 200 0.9795968 0.9773288
FALSE 27 0.20 7 10 200 0.9809077 0.9787856
FALSE AccuracySD KappaSD
FALSE 1 0.019612067 0.021788936
FALSE 10 0.019622632 0.021802519
FALSE 19 0.015921354 0.017691203
FALSE 4 0.012320755 0.013690221
FALSE 13 0.005404789 0.006005197
FALSE 22 0.007516921 0.008352553
FALSE 7 0.008586885 0.009541220
FALSE 16 0.006918735 0.007687701
FALSE 25 0.004478137 0.004976086
FALSE 2 0.017833560 0.019814109
FALSE 11 0.014865967 0.016517752
FALSE 20 0.006919845 0.007689247
FALSE 5 0.006581881 0.007313589
FALSE 14 0.004936928 0.005485736
FALSE 23 0.005237421 0.005819778
FALSE 8 0.004743494 0.005270815
FALSE 17 0.004520332 0.005022748
FALSE 26 0.003777682 0.004197459
FALSE 3 0.015748937 0.017498485
FALSE 12 0.006815428 0.007573313
FALSE 21 0.007314194 0.008127080
FALSE 6 0.004892768 0.005436904
FALSE 15 0.004848586 0.005387668
FALSE 24 0.004820809 0.005356488
FALSE 9 0.005536006 0.006151496
FALSE 18 0.005059983 0.005622406
FALSE 27 0.004090623 0.004545276Using gbmGrid_forever and keeping the minimal number of observations per tree leaf fixed as 10, optimal parameters are 5, 0.1 and 100 respective to depth, learning rate and number of trees.
Accuracy of the final model on test data is 95.83%. On the other hand, for the train data the accuracy is 100%.
Error rates on test and train data are 4.17% and 0%.
Again, the model accuracy is higher on the train data, almost 100%, which may indicate overfitting.
6 Mushroom Types
mushroom_raw <- read.csv("https://raw.githubusercontent.com/BU-IE-582/fall20-ilaydacelenk/master/files/HW4_data/mushroom", header=FALSE)
mushroom_colns <- c("class", "cap-shape", "cap-surface", "cap-color", "bruises", "odor", "gill-attachment", "gill-spacing", "gill-size", "gill-color", "stalk-shape", "stalk-root", "stalk-surface-above-ring", "stalk-surface-below-ring", "stalk-color-above-ring", "stalk-color-below-ring", "veil-type", "veil-color", "ring-number", "ring-type", "spore-print-color", "population", "habitat")
colnames(mushroom_raw) <- mushroom_colns
#Balanced binary classification
mushroom_raw %>% group_by(class) %>% summarize(count=n()) FALSE # A tibble: 2 x 2
FALSE class count
FALSE <chr> <int>
FALSE 1 e 4208
FALSE 2 p 3916#veil-type is same for all the rows, therefore I will remove it
mushroom_raw %>% summarise_all(n_distinct)FALSE class cap-shape cap-surface cap-color bruises odor gill-attachment
FALSE 1 2 6 4 10 2 9 2
FALSE gill-spacing gill-size gill-color stalk-shape stalk-root
FALSE 1 2 2 12 2 5
FALSE stalk-surface-above-ring stalk-surface-below-ring stalk-color-above-ring
FALSE 1 4 4 9
FALSE stalk-color-below-ring veil-type veil-color ring-number ring-type
FALSE 1 9 1 4 3 5
FALSE spore-print-color population habitat
FALSE 1 9 6 7mushroom_raw$`veil-type` <- NULL
#no na values
mushroom_raw %>% summarise_all(~ sum(is.na(.))) FALSE class cap-shape cap-surface cap-color bruises odor gill-attachment
FALSE 1 0 0 0 0 0 0 0
FALSE gill-spacing gill-size gill-color stalk-shape stalk-root
FALSE 1 0 0 0 0 0
FALSE stalk-surface-above-ring stalk-surface-below-ring stalk-color-above-ring
FALSE 1 0 0 0
FALSE stalk-color-below-ring veil-color ring-number ring-type spore-print-color
FALSE 1 0 0 0 0 0
FALSE population habitat
FALSE 1 0 0mushroom_raw$class <- factor(mushroom_raw$class)
cols <- names(mushroom_raw[c(2:22)])
mushroom_raw[cols] <- lapply(mushroom_raw[cols], factor)
#method='reference' so that n-1 dummy columns
mushrooms_dummy <- mushroom_raw %>% mlr::createDummyFeatures(., method='reference')
set.seed(50)
train_inds_mushroom = createDataPartition(y = 1:nrow(mushroom_raw), p = 0.8, list = F)
mushroom_train <- mushroom_raw[train_inds_mushroom, ] %>% as.data.frame()
mushroom_test <- mushroom_raw[-train_inds_mushroom, ] %>% as.data.frame()
mushroom_train_dummy <- mushrooms_dummy[train_inds_mushroom, ] %>% as.data.frame()
mushroom_test_dummy <- mushrooms_dummy[-train_inds_mushroom, ] %>% as.data.frame()
X_mushroom_train_dummy <- mushroom_train_dummy %>% select(-p)
X_mushroom_test_dummy <- mushroom_test_dummy %>% select(-p)
y_mushroom_train_dummy <- mushroom_train_dummy %>% select(p)
y_mushroom_test_dummy <- mushroom_test_dummy %>% select(p)
X_mushroom_train <- mushroom_train %>% select(-class)
y_mushroom_train <- mushroom_train %>% select(class)
X_mushroom_test <- mushroom_test %>% select(-class)
y_mushroom_test <- mushroom_test %>% select(class)
ggplot(mushroom_raw, aes(x=class )) + geom_bar(color="blue", fill=rgb(0.1,0.4,0.5,0.7) )
Mushroom Types dataset consists of descriptions of mushrooms from 23 species. This is a binary classification problem where we try to find if a mushroom is edible or poisonous. The first column corresponds to the target(dependent) variable and the rest of the columns are the predictors(independent variable). The feature “veil-type” is same for all the rows and it will not contribute in the model, therefore I will remove it.
Here, there is no separate train and test data, therefore we need to generate them by using createDataPartition function of the package caret. I preferred to use 20% of the total data as the test data. Attribute information can be found here, such as color, size, shape, etc. All of the independent variables are categorical and they are all nominal.
Classes are denoted as “e” and “p”, meaning edible and poisonous, respectively. We also see that classes are balanced.
6.1 Penalized Regression Approach
set.seed(50)
mushroom_model_pra<-cv.glmnet(as.matrix(X_mushroom_train_dummy), mushroom_train_dummy$p, type.measure = "class", family="binomial", nfolds=10)
mushroom_pred_pra <- predict(mushroom_model_pra, newx=as.matrix(X_mushroom_test_dummy),s=c("lambda.min"), type = "class")
mushroom_cm_pra_test <- confusionMatrix(factor(mushroom_pred_pra), factor(y_mushroom_test_dummy$p))
mushroom_cm_pra_train <- confusionMatrix(factor(predict(mushroom_model_pra, newx=as.matrix(X_mushroom_train_dummy), s=c("lambda.min"), type = "class")), factor(y_mushroom_train_dummy$p))
print(mushroom_model_pra)FALSE
FALSE Call: cv.glmnet(x = as.matrix(X_mushroom_train_dummy), y = mushroom_train_dummy$p, type.measure = "class", nfolds = 10, family = "binomial")
FALSE
FALSE Measure: Misclassification Error
FALSE
FALSE Lambda Measure SE Nonzero
FALSE min 0.0002103 0 0 31
FALSE 1se 0.0002103 0 0 31plot(mushroom_model_pra)
Lambda value that minimizes the error is 2.103342810^{-4} with error rate 0%, the error rate on the test data is also 0%.
Accuracy of the model on test data is 100%. On the other hand, accuracy on the train data is 100%. Both of them are 1, which seems like the perfect model but it is very unlikely.
Error rates on test and train data are 0% and 0%.
6.2 Decision Tree
set.seed(50)
complexity=c(0.01,0.03,0.05)
min_number_of_observations=c(10,15,25)
cv_summary=data.table()
#for balanced folds: stratified =TRUE
index <- generateCVRuns(mushroom_train$class, ntimes=1, nfold=10, stratified =TRUE)
acc <- 0
best_complexity <- 0
best_min_number_of_observations <- 0
for(i in 1:10){
for(c in complexity){
for(m in min_number_of_observations){
test_index <- index$`Run 1`[i]
train_cv <- mushroom_train[-test_index[[1]],]
test_cv <- mushroom_train[test_index[[1]],]
fitControl <- rpart.control(cp = c, minbucket = m)
mushroom_model_DT=rpart(class~., mushroom_train, control=fitControl)
pred <- predict(mushroom_model_DT, test_cv, type="class")
cm <- confusionMatrix(pred, factor(test_cv$class))
acc_temp <- cm$overall[1:1]
cv_summary <- rbind(cv_summary,data.table(Fold=i, Complexity=c, Min_Observations_per_Leaf=m, Accuracy=acc_temp))
if(acc_temp>acc){
acc<-acc_temp
best_complexity <- c
best_min_number_of_observations <- m
}
}
}
}
cv_summary <- cv_summary %>% group_by(Complexity, Min_Observations_per_Leaf) %>% summarise(mean_accuracy = mean(Accuracy))
print(cv_summary)FALSE # A tibble: 9 x 3
FALSE # Groups: Complexity [3]
FALSE Complexity Min_Observations_per_Leaf mean_accuracy
FALSE <dbl> <dbl> <dbl>
FALSE 1 0.01 10 0.994
FALSE 2 0.01 15 0.994
FALSE 3 0.01 25 0.994
FALSE 4 0.03 10 0.986
FALSE 5 0.03 15 0.986
FALSE 6 0.03 25 0.986
FALSE 7 0.05 10 0.986
FALSE 8 0.05 15 0.986
FALSE 9 0.05 25 0.986fitControl <- rpart.control(cp = best_complexity, minbucket = best_min_number_of_observations)
mushroom_model_DT <- rpart(class~., mushroom_train, control=fitControl)
mushroom_pred_DT <- predict(mushroom_model_DT, X_mushroom_test, type="class")
mushroom_cm_DT_test <- confusionMatrix(factor(mushroom_pred_DT), factor(mushroom_test$class))
mushroom_cm_DT_train <- confusionMatrix(factor(predict(mushroom_model_DT, X_mushroom_train, type="class")), factor(mushroom_train$class))
print(mushroom_cm_DT_test[["table"]])FALSE Reference
FALSE Prediction e p
FALSE e 838 11
FALSE p 0 775According to the cross validation result, cross validation error rate is 0.57, the error rate on the test data is 0.68% which is slightly higher but still both very low.
The best complexity parameter is 0.01 and the best minimal number of observations per tree leaf is 10 since it gives the highest accuracy of the model on cross validation test data with 99.85%.
Now I can fit the final model and make predictions. Accuracy of the final model on test data is 99.32%. On the other hand, for the train data the accuracy is 99.43%, which are very close and high.
Error rates on test and train data are 0.68% and 0.57%.
This indicates a good model.
6.3 Random Forest
set.seed(50)
m=c(4,6,8)
#for balanced folds: stratified =TRUE
# ntree is J=500 given in hw
J=500
# the minimal number of observations per tree leaf is nodesize 5 given in hw
min_size_of_terminal_nodes=5
n_folds=10
fitControl=trainControl(method = "repeatedcv",
number = n_folds,
repeats = 1, classProbs = TRUE, search="random")
rf_grid=expand.grid(mtry = m)
mushroom_model_RF <- train(class ~ ., data = mushroom_train, method = "rf", trControl = fitControl, ntree=J, nodesize = min_size_of_terminal_nodes, tuneGrid = rf_grid)
print(mushroom_model_RF[["results"]])FALSE mtry Accuracy Kappa AccuracySD KappaSD
FALSE 1 4 0.9998462 0.9996919 0.0004865043 0.0009744477
FALSE 2 6 1.0000000 1.0000000 0.0000000000 0.0000000000
FALSE 3 8 1.0000000 1.0000000 0.0000000000 0.0000000000mushroom_pred_RF <- predict(mushroom_model_RF, newdata=X_mushroom_test)
mushroom_cm_RF_test <- confusionMatrix(factor(mushroom_pred_RF), mushroom_test$class)
mushroom_cm_RF_train <- confusionMatrix(factor(predict(mushroom_model_RF, newdata=X_mushroom_train)), mushroom_train$class)
print(mushroom_cm_RF_test[["table"]])FALSE Reference
FALSE Prediction e p
FALSE e 838 0
FALSE p 0 786According to the repeated cross validation result, the best m parameter is 6 since it gives the highest accuracy 100% of the model on cross validation test data.
Accuracy of the final model on test data is 100%. On the other hand, for the train data the accuracy is 100%. Both are 1, this is not realistic.
Error rates on test and train data are 0% and 0%.
6.4 Stochastic Gradient Boosting
set.seed(50)
n_folds <- 10
fitControl <- trainControl(method = "cv", number = n_folds)
gbmGrid_forever <- expand.grid(interaction.depth = c(1, 5, 7), n.trees = c(50,100,200), shrinkage = c(0.05, 0.1, 0.2), n.minobsinnode = c(10))
print(gbmGrid_forever)FALSE interaction.depth n.trees shrinkage n.minobsinnode
FALSE 1 1 50 0.05 10
FALSE 2 5 50 0.05 10
FALSE 3 7 50 0.05 10
FALSE 4 1 100 0.05 10
FALSE 5 5 100 0.05 10
FALSE 6 7 100 0.05 10
FALSE 7 1 200 0.05 10
FALSE 8 5 200 0.05 10
FALSE 9 7 200 0.05 10
FALSE 10 1 50 0.10 10
FALSE 11 5 50 0.10 10
FALSE 12 7 50 0.10 10
FALSE 13 1 100 0.10 10
FALSE 14 5 100 0.10 10
FALSE 15 7 100 0.10 10
FALSE 16 1 200 0.10 10
FALSE 17 5 200 0.10 10
FALSE 18 7 200 0.10 10
FALSE 19 1 50 0.20 10
FALSE 20 5 50 0.20 10
FALSE 21 7 50 0.20 10
FALSE 22 1 100 0.20 10
FALSE 23 5 100 0.20 10
FALSE 24 7 100 0.20 10
FALSE 25 1 200 0.20 10
FALSE 26 5 200 0.20 10
FALSE 27 7 200 0.20 10gbmGrid <- expand.grid(interaction.depth = c(5), n.trees = c(100), shrinkage = c(0.1), n.minobsinnode = c(10))
mushroom_train_dummy$p <- factor(mushroom_train_dummy$p)
mushroom_model_sgb=train(p ~ ., data = data.frame(mushroom_train_dummy), method = "gbm", trControl = fitControl, tuneGrid = gbmGrid_forever, verbose=F)
mushroom_pred_sgb <- predict(mushroom_model_sgb, newdata=X_mushroom_test_dummy) %>% as.matrix()
mushroom_cm_sgb_test <- confusionMatrix(factor(mushroom_pred_sgb), factor(y_mushroom_test_dummy$p))
mushroom_cm_sgb_train <- confusionMatrix(factor(predict(mushroom_model_sgb, newdata=X_mushroom_train_dummy)), factor(y_mushroom_train_dummy$p))
print(mushroom_model_sgb[["results"]])FALSE shrinkage interaction.depth n.minobsinnode n.trees Accuracy Kappa
FALSE 1 0.05 1 10 50 0.9703077 0.9405516
FALSE 10 0.10 1 10 50 0.9746154 0.9491508
FALSE 19 0.20 1 10 50 0.9918462 0.9836716
FALSE 4 0.05 5 10 50 0.9970769 0.9941441
FALSE 13 0.10 5 10 50 0.9989231 0.9978428
FALSE 22 0.20 5 10 50 1.0000000 1.0000000
FALSE 7 0.05 7 10 50 0.9984615 0.9969181
FALSE 16 0.10 7 10 50 0.9995385 0.9990753
FALSE 25 0.20 7 10 50 1.0000000 1.0000000
FALSE 2 0.05 1 10 100 0.9746154 0.9491508
FALSE 11 0.10 1 10 100 0.9890769 0.9781310
FALSE 20 0.20 1 10 100 0.9949231 0.9898323
FALSE 5 0.05 5 10 100 0.9989231 0.9978428
FALSE 14 0.10 5 10 100 1.0000000 1.0000000
FALSE 23 0.20 5 10 100 1.0000000 1.0000000
FALSE 8 0.05 7 10 100 0.9996923 0.9993836
FALSE 17 0.10 7 10 100 1.0000000 1.0000000
FALSE 26 0.20 7 10 100 1.0000000 1.0000000
FALSE 3 0.05 1 10 200 0.9883077 0.9765924
FALSE 12 0.10 1 10 200 0.9949231 0.9898323
FALSE 21 0.20 1 10 200 0.9989231 0.9978428
FALSE 6 0.05 5 10 200 1.0000000 1.0000000
FALSE 15 0.10 5 10 200 1.0000000 1.0000000
FALSE 24 0.20 5 10 200 1.0000000 1.0000000
FALSE 9 0.05 7 10 200 1.0000000 1.0000000
FALSE 18 0.10 7 10 200 1.0000000 1.0000000
FALSE 27 0.20 7 10 200 1.0000000 1.0000000
FALSE AccuracySD KappaSD
FALSE 1 0.0150397239 0.030082895
FALSE 10 0.0098175332 0.019670301
FALSE 19 0.0039755995 0.007958704
FALSE 4 0.0018418754 0.003690078
FALSE 13 0.0014595128 0.002923910
FALSE 22 0.0000000000 0.000000000
FALSE 7 0.0017764624 0.003558990
FALSE 16 0.0014595128 0.002924008
FALSE 25 0.0000000000 0.000000000
FALSE 2 0.0098175332 0.019670301
FALSE 11 0.0051512304 0.010304375
FALSE 20 0.0029946259 0.005998516
FALSE 5 0.0014595128 0.002923910
FALSE 14 0.0000000000 0.000000000
FALSE 23 0.0000000000 0.000000000
FALSE 8 0.0009730085 0.001949117
FALSE 17 0.0000000000 0.000000000
FALSE 26 0.0000000000 0.000000000
FALSE 3 0.0056270222 0.011258673
FALSE 12 0.0034055298 0.006821198
FALSE 21 0.0014595128 0.002923910
FALSE 6 0.0000000000 0.000000000
FALSE 15 0.0000000000 0.000000000
FALSE 24 0.0000000000 0.000000000
FALSE 9 0.0000000000 0.000000000
FALSE 18 0.0000000000 0.000000000
FALSE 27 0.0000000000 0.000000000Using gbmGrid_forever and keeping the minimal number of observations per tree leaf fixed as 10, optimal parameters are 5, 0.1 and 100 respective to depth, learning rate and number of trees.
Accuracy of the final model on test data is 100%. On the other hand, for the train data the accuracy is 100%. Both accuracies are 1, similar to the other methods.
7 Telco Customer Churn
churn_raw <- read.csv("https://raw.githubusercontent.com/BU-IE-582/fall20-ilaydacelenk/master/files/HW4_data/churn.csv")
#to see the imbalance
ggplot(churn_raw, aes(x=Churn )) + geom_bar(color="blue", fill=rgb(0.1,0.4,0.5,0.7) )
#customer id is unie for every row, so I will remove it
churn_raw$customerID <- NULL
churn_raw %>% summarise_all(~ sum(is.na(.))) #to check na valuesFALSE gender SeniorCitizen Partner Dependents tenure PhoneService MultipleLines
FALSE 1 0 0 0 0 0 0 0
FALSE InternetService OnlineSecurity OnlineBackup DeviceProtection TechSupport
FALSE 1 0 0 0 0 0
FALSE StreamingTV StreamingMovies Contract PaperlessBilling PaymentMethod
FALSE 1 0 0 0 0 0
FALSE MonthlyCharges TotalCharges Churn
FALSE 1 0 11 0#remove na
churn_raw <- churn_raw %>% drop_na()
var_mode <- sapply(churn_raw, mode)
# find character columns
colns <- which(var_mode %in% c("character"))
# character to factor
churn_raw[colns] <- lapply(churn_raw[colns], factor)
glimpse(churn_raw)FALSE Rows: 7,032
FALSE Columns: 20
FALSE $ gender <fct> Female, Male, Male, Male, Female, Female, Male, Fema…
FALSE $ SeniorCitizen <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0…
FALSE $ Partner <fct> Yes, No, No, No, No, No, No, No, Yes, No, Yes, No, Y…
FALSE $ Dependents <fct> No, No, No, No, No, No, Yes, No, No, Yes, Yes, No, N…
FALSE $ tenure <int> 1, 34, 2, 45, 2, 8, 22, 10, 28, 62, 13, 16, 58, 49, …
FALSE $ PhoneService <fct> No, Yes, Yes, No, Yes, Yes, Yes, No, Yes, Yes, Yes, …
FALSE $ MultipleLines <fct> No phone service, No, No, No phone service, No, Yes,…
FALSE $ InternetService <fct> DSL, DSL, DSL, DSL, Fiber optic, Fiber optic, Fiber …
FALSE $ OnlineSecurity <fct> No, Yes, Yes, Yes, No, No, No, Yes, No, Yes, Yes, No…
FALSE $ OnlineBackup <fct> Yes, No, Yes, No, No, No, Yes, No, No, Yes, No, No i…
FALSE $ DeviceProtection <fct> No, Yes, No, Yes, No, Yes, No, No, Yes, No, No, No i…
FALSE $ TechSupport <fct> No, No, No, Yes, No, No, No, No, Yes, No, No, No int…
FALSE $ StreamingTV <fct> No, No, No, No, No, Yes, Yes, No, Yes, No, No, No in…
FALSE $ StreamingMovies <fct> No, No, No, No, No, Yes, No, No, Yes, No, No, No int…
FALSE $ Contract <fct> Month-to-month, One year, Month-to-month, One year, …
FALSE $ PaperlessBilling <fct> Yes, No, Yes, No, Yes, Yes, Yes, No, Yes, No, Yes, N…
FALSE $ PaymentMethod <fct> Electronic check, Mailed check, Mailed check, Bank t…
FALSE $ MonthlyCharges <dbl> 29.85, 56.95, 53.85, 42.30, 70.70, 99.65, 89.10, 29.…
FALSE $ TotalCharges <dbl> 29.85, 1889.50, 108.15, 1840.75, 151.65, 820.50, 194…
FALSE $ Churn <fct> No, No, Yes, No, Yes, Yes, No, No, Yes, No, No, No, …#This time I used caret's function to create dummy variables because mlr's function changed the column names in a very complicated way
create_dummy <- dummyVars("~ .", data=churn_raw, fullRank = T)
churn_dummy = data.frame(predict(create_dummy, newdata = churn_raw))
churn_dummy$Churn.Yes <- as.factor(churn_dummy$Churn.Yes)
set.seed(50)
train_inds_churn = createDataPartition(y = 1:nrow(churn_dummy), p = 0.8, list = F)
churn_train_dummy <- churn_dummy[train_inds_churn, ]
churn_test_dummy <- churn_dummy[-train_inds_churn, ]
churn_train_dummy_smoted <- DMwR::SMOTE(Churn.Yes ~ ., churn_train_dummy, perc.over = 100, perc.under=200)
ggplot(churn_train_dummy_smoted, aes(x=Churn.Yes)) + geom_bar(color="blue", fill=rgb(0.1,0.4,0.5,0.7))
levels(churn_train_dummy_smoted$Churn.Yes) <- c("NO", "YES")
levels(churn_test_dummy$Churn.Yes) <- c("NO", "YES")
levels(churn_train_dummy$Churn.Yes) <- c("NO", "YES")
#model using the smote data and make predictions on non-smoted test data
# also validate using non-smote train data
X_churn_train_dummy_smoted <- churn_train_dummy_smoted %>% select(-Churn.Yes)
y_churn_train_dummy_smoted <- churn_train_dummy_smoted %>% select(Churn.Yes)
X_churn_train_dummy <- churn_train_dummy %>% select(-Churn.Yes)
y_churn_train_dummy <- churn_train_dummy %>% select(Churn.Yes)
X_churn_test_dummy <- churn_test_dummy %>% select(-Churn.Yes)
y_churn_test_dummy <- churn_test_dummy %>% select(Churn.Yes)Telco Customer Churn data contains information on customers churned and not churned. Therefore, it is a binary classification problem. There are 7032 instances with 1869 customers have churned and 5163 customers have not churned. Therefore, the data is imbalanced since only 27% of the data is class of churned.
There are 21 columns and some features are numerical, some are categorical. Again, there is no separate train and test data, I use 20% of the data for test.
We see that classes are imbalanced. Therefore I use smote function in order to balance the classes. This function does oversampling in other words it creates more samples from the minority class. For the models below, I used smote data for the fitting the model and then performed on the non-smote train and test data.
7.1 Penalized Regression Approach
set.seed(50)
churn_model_pra<-cv.glmnet(as.matrix(X_churn_train_dummy_smoted), churn_train_dummy_smoted$Churn.Yes, type.measure = "class", family="binomial", nfolds=10)
print(churn_model_pra)FALSE
FALSE Call: cv.glmnet(x = as.matrix(X_churn_train_dummy_smoted), y = churn_train_dummy_smoted$Churn.Yes, type.measure = "class", nfolds = 10, family = "binomial")
FALSE
FALSE Measure: Misclassification Error
FALSE
FALSE Lambda Measure SE Nonzero
FALSE min 0.00039 0.2269 0.004010 28
FALSE 1se 0.07072 0.2304 0.002871 7churn_pred_pra <- predict(churn_model_pra, newx=as.matrix(X_churn_test_dummy),s=c("lambda.min"), type = "class")
churn_cm_pra_test <- confusionMatrix(factor(churn_pred_pra), factor(y_churn_test_dummy$Churn.Yes))
churn_cm_pra_train <- confusionMatrix(factor(predict(churn_model_pra, newx=as.matrix(X_churn_train_dummy), s=c("lambda.min"), type = "class")), factor(y_churn_train_dummy$Churn.Yes))
print(churn_model_pra)FALSE
FALSE Call: cv.glmnet(x = as.matrix(X_churn_train_dummy_smoted), y = churn_train_dummy_smoted$Churn.Yes, type.measure = "class", nfolds = 10, family = "binomial")
FALSE
FALSE Measure: Misclassification Error
FALSE
FALSE Lambda Measure SE Nonzero
FALSE min 0.00039 0.2269 0.004010 28
FALSE 1se 0.07072 0.2304 0.002871 7plot(churn_model_pra)
Lambda value that minimizes the error is 3.86317310^{-4} with error rate 22.69% on the cross validation data, the error rate on the test data is 24.57% which is higher.
Accuracy of the model on non-smote test data is 75.43%. On the other hand, accuracy on the non-smote train data is 74.93%. Accuracy on the test data is larger, this means that the model is not memorizing i.e there is no overfitting.
Error rates on test and train data are 24.57% and 25.07%.
Since the accuracy of test and train data very close, there is no overfitting. But still the accuracy is not good enough.
7.2 Decision Tree
set.seed(50)
complexity=c(0.01,0.03,0.05)
min_number_of_observations=c(10,15,25)
cv_summary=data.table()
#for balanced folds: stratified =TRUE
index <- generateCVRuns(churn_train_dummy_smoted$Churn.Yes, ntimes=1, nfold=10, stratified =TRUE)
acc <- 0
best_complexity <- 0
best_min_number_of_observations <- 0
for(i in 1:10){
for(c in complexity){
for(m in min_number_of_observations){
test_index <- index$`Run 1`[i]
train_cv <- churn_train_dummy_smoted[-test_index[[1]],]
test_cv <- churn_train_dummy_smoted[test_index[[1]],]
fitControl <- rpart.control(cp = c, minbucket = m)
churn_model_DT=rpart(Churn.Yes~., churn_train_dummy_smoted, control=fitControl)
pred <- predict(churn_model_DT, test_cv, type="class")
cm <- confusionMatrix(pred, factor(test_cv$Churn.Yes))
acc_temp <- cm$overall[1:1]
cv_summary <- rbind(cv_summary,data.table(Fold=i, Complexity=c, Min_Observations_per_Leaf=m, Accuracy=acc_temp))
if(acc_temp>acc){
acc<-acc_temp
best_complexity <- c
best_min_number_of_observations <- m
}
}
}
}
cv_summary <- cv_summary %>% group_by(Complexity, Min_Observations_per_Leaf) %>% summarise(mean_accuracy = mean(Accuracy))
print(cv_summary)FALSE # A tibble: 9 x 3
FALSE # Groups: Complexity [3]
FALSE Complexity Min_Observations_per_Leaf mean_accuracy
FALSE <dbl> <dbl> <dbl>
FALSE 1 0.01 10 0.774
FALSE 2 0.01 15 0.774
FALSE 3 0.01 25 0.774
FALSE 4 0.03 10 0.774
FALSE 5 0.03 15 0.774
FALSE 6 0.03 25 0.774
FALSE 7 0.05 10 0.744
FALSE 8 0.05 15 0.744
FALSE 9 0.05 25 0.744fitControl <- rpart.control(cp = best_complexity, minbucket = best_min_number_of_observations)
churn_model_DT <- rpart(Churn.Yes~., churn_train_dummy_smoted, control=fitControl)
churn_pred_DT <- predict(churn_model_DT, X_churn_test_dummy, type="class")
churn_cm_DT_test <- confusionMatrix(factor(churn_pred_DT), factor(y_churn_test_dummy$Churn.Yes))
churn_cm_DT_train <- confusionMatrix(factor(predict(churn_model_DT, X_churn_train_dummy, type="class")), factor(y_churn_train_dummy$Churn.Yes))
print(churn_cm_DT_test[["table"]])FALSE Reference
FALSE Prediction NO YES
FALSE NO 715 87
FALSE YES 290 312According to the cross validation result, the best complexity parameter is 0.01 and the best minimal number of observations per tree leaf is 10 since it gives the highest accuracy of the model on cross validation test data with 81.29%.
Now I can fit the final model and make predictions. Accuracy of the final model on test data is 73.15%. On the other hand, for the train data the accuracy is 73.67%. Accuracy values for the non-smote test and train data are very close to each other. So there is no overfitting. But the cross validation results suggested 81.29% accuracy on the smote train data. So we can say that the model is not learning enough and there is no overfitting.
Error rates on test and train data are 26.85% and 26.33%.
7.3 Random Forest
set.seed(50)
m=c(4,6,8)
#for balanced folds: stratified =TRUE
# ntree is J=500 given in hw
J=500
# the minimal number of observations per tree leaf is nodesize 5 given in hw
min_size_of_terminal_nodes=5
n_folds=10
fitControl=trainControl(method = "repeatedcv", number = n_folds, repeats = 1, classProbs = TRUE, search="random", savePred=TRUE)
rf_grid=expand.grid(mtry = m)
churn_model_RF <- train(Churn.Yes ~ ., data = churn_train_dummy_smoted, method = "rf", trControl = fitControl, ntree=J, nodesize = min_size_of_terminal_nodes, tuneGrid = rf_grid)
print(churn_model_RF[["results"]])FALSE mtry Accuracy Kappa AccuracySD KappaSD
FALSE 1 4 0.8370748 0.6741497 0.01413373 0.02826747
FALSE 2 6 0.8481293 0.6962585 0.01418480 0.02836960
FALSE 3 8 0.8522109 0.7044218 0.01496534 0.02993068churn_pred_RF <- predict(churn_model_RF, newdata=X_churn_test_dummy)
churn_cm_RF_test <- confusionMatrix(factor(churn_pred_RF), churn_test_dummy$Churn.Yes)
churn_cm_RF_train <- confusionMatrix(factor(predict(churn_model_RF, newdata=X_churn_train_dummy)), churn_train_dummy$Churn.Yes)
print(churn_cm_RF_test[["table"]])FALSE Reference
FALSE Prediction NO YES
FALSE NO 779 120
FALSE YES 226 279According to the repeated cross validation result, the best m parameter is 8 since it gives the highest accuracy of the model on cross validation test data.
Accuracy of the final model on test data is 75.36%. On the other hand, for the train data the accuracy is 89.36%.
Error rates on test and train data are 24.64% and 10.64%.
The model is performing better on the train set. Therefore, there is significant overfitting in the model.
7.4 Stochastic Gradient Boosting
set.seed(50)
n_folds <- 10
fitControl <- trainControl(method = "cv", number = n_folds)
gbmGrid_forever <- expand.grid(interaction.depth = c(1, 5, 7), n.trees = c(50,100,200), shrinkage = c(0.05, 0.1, 0.2), n.minobsinnode = c(10))
gbmGrid <- expand.grid(interaction.depth = c(5), n.trees = c(100), shrinkage = c(0.1), n.minobsinnode = c(10))
churn_model_sgb=train(Churn.Yes ~ ., data = data.frame(churn_train_dummy_smoted), method = "gbm", trControl = fitControl, tuneGrid = gbmGrid_forever, verbose=F)
churn_pred_sgb <- predict(churn_model_sgb, newdata=X_churn_test_dummy) %>% as.matrix()
churn_cm_sgb_test <- confusionMatrix(factor(churn_pred_sgb), factor(y_churn_test_dummy$Churn.Yes))
churn_cm_sgb_train <- confusionMatrix(factor(predict(churn_model_sgb, newdata=X_churn_train_dummy)), factor(y_churn_train_dummy$Churn.Yes))
print(churn_model_sgb[["results"]])FALSE shrinkage interaction.depth n.minobsinnode n.trees Accuracy Kappa
FALSE 1 0.05 1 10 50 0.7678571 0.5357143
FALSE 10 0.10 1 10 50 0.7806122 0.5612245
FALSE 19 0.20 1 10 50 0.7846939 0.5693878
FALSE 4 0.05 5 10 50 0.7959184 0.5918367
FALSE 13 0.10 5 10 50 0.8112245 0.6224490
FALSE 22 0.20 5 10 50 0.8193878 0.6387755
FALSE 7 0.05 7 10 50 0.7996599 0.5993197
FALSE 16 0.10 7 10 50 0.8141156 0.6282313
FALSE 25 0.20 7 10 50 0.8260204 0.6520408
FALSE 2 0.05 1 10 100 0.7806122 0.5612245
FALSE 11 0.10 1 10 100 0.7863946 0.5727891
FALSE 20 0.20 1 10 100 0.7926871 0.5853741
FALSE 5 0.05 5 10 100 0.8090136 0.6180272
FALSE 14 0.10 5 10 100 0.8217687 0.6435374
FALSE 23 0.20 5 10 100 0.8246599 0.6493197
FALSE 8 0.05 7 10 100 0.8120748 0.6241497
FALSE 17 0.10 7 10 100 0.8238095 0.6476190
FALSE 26 0.20 7 10 100 0.8312925 0.6625850
FALSE 3 0.05 1 10 200 0.7875850 0.5751701
FALSE 12 0.10 1 10 200 0.7925170 0.5850340
FALSE 21 0.20 1 10 200 0.7991497 0.5982993
FALSE 6 0.05 5 10 200 0.8209184 0.6418367
FALSE 15 0.10 5 10 200 0.8294218 0.6588435
FALSE 24 0.20 5 10 200 0.8251701 0.6503401
FALSE 9 0.05 7 10 200 0.8263605 0.6527211
FALSE 18 0.10 7 10 200 0.8333333 0.6666667
FALSE 27 0.20 7 10 200 0.8392857 0.6785714
FALSE AccuracySD KappaSD
FALSE 1 0.01517433 0.03034866
FALSE 10 0.01213205 0.02426410
FALSE 19 0.01565283 0.03130566
FALSE 4 0.01397826 0.02795653
FALSE 13 0.01254872 0.02509744
FALSE 22 0.01522824 0.03045648
FALSE 7 0.01477842 0.02955685
FALSE 16 0.01346417 0.02692833
FALSE 25 0.01334429 0.02668858
FALSE 2 0.01465176 0.02930351
FALSE 11 0.01647311 0.03294621
FALSE 20 0.01412577 0.02825155
FALSE 5 0.01260111 0.02520222
FALSE 14 0.01395065 0.02790129
FALSE 23 0.01705304 0.03410607
FALSE 8 0.01128816 0.02257631
FALSE 17 0.01295696 0.02591391
FALSE 26 0.01687020 0.03374040
FALSE 3 0.01530507 0.03061015
FALSE 12 0.01454167 0.02908335
FALSE 21 0.01184929 0.02369859
FALSE 6 0.01458691 0.02917381
FALSE 15 0.01794556 0.03589111
FALSE 24 0.01438167 0.02876335
FALSE 9 0.01227557 0.02455113
FALSE 18 0.01639097 0.03278193
FALSE 27 0.01809181 0.03618361Using gbmGrid_forever and keeping the minimal number of observations per tree leaf fixed as 10, optimal parameters are 5, 0.1 and 100 respective to depth, learning rate and number of trees.
Accuracy of the final model on test data is 75.57%. On the other hand, for the train data the accuracy is 82.96%. The model performs better on the train data, which may indicate overfitting.
8 House Prices
houseprice_raw <- read.csv("https://raw.githubusercontent.com/BU-IE-582/fall20-ilaydacelenk/master/files/HW4_data/house_price_train.csv")
houseprice <- houseprice_raw
houseprice$Alley[is.na(houseprice$Alley)] <- "none"
houseprice$FireplaceQu[is.na(houseprice$FireplaceQu)] <- "none"
houseprice$PoolQC[is.na(houseprice$PoolQC)] <- "none"
houseprice$Fence[is.na(houseprice$Fence)] <- "none"
houseprice$MiscFeature[is.na(houseprice$MiscFeature)] <- "none"
houseprice <- houseprice %>% drop_na()
glimpse(houseprice)FALSE Rows: 1,094
FALSE Columns: 81
FALSE $ Id <int> 1, 2, 3, 4, 5, 6, 7, 9, 10, 11, 12, 14, 16, 19, 20, 21,…
FALSE $ MSSubClass <int> 60, 20, 60, 70, 60, 50, 20, 50, 190, 20, 60, 20, 45, 20…
FALSE $ MSZoning <chr> "RL", "RL", "RL", "RL", "RL", "RL", "RL", "RM", "RL", "…
FALSE $ LotFrontage <int> 65, 80, 68, 60, 84, 85, 75, 51, 50, 70, 85, 91, 51, 66,…
FALSE $ LotArea <int> 8450, 9600, 11250, 9550, 14260, 14115, 10084, 6120, 742…
FALSE $ Street <chr> "Pave", "Pave", "Pave", "Pave", "Pave", "Pave", "Pave",…
FALSE $ Alley <chr> "none", "none", "none", "none", "none", "none", "none",…
FALSE $ LotShape <chr> "Reg", "Reg", "IR1", "IR1", "IR1", "IR1", "Reg", "Reg",…
FALSE $ LandContour <chr> "Lvl", "Lvl", "Lvl", "Lvl", "Lvl", "Lvl", "Lvl", "Lvl",…
FALSE $ Utilities <chr> "AllPub", "AllPub", "AllPub", "AllPub", "AllPub", "AllP…
FALSE $ LotConfig <chr> "Inside", "FR2", "Inside", "Corner", "FR2", "Inside", "…
FALSE $ LandSlope <chr> "Gtl", "Gtl", "Gtl", "Gtl", "Gtl", "Gtl", "Gtl", "Gtl",…
FALSE $ Neighborhood <chr> "CollgCr", "Veenker", "CollgCr", "Crawfor", "NoRidge", …
FALSE $ Condition1 <chr> "Norm", "Feedr", "Norm", "Norm", "Norm", "Norm", "Norm"…
FALSE $ Condition2 <chr> "Norm", "Norm", "Norm", "Norm", "Norm", "Norm", "Norm",…
FALSE $ BldgType <chr> "1Fam", "1Fam", "1Fam", "1Fam", "1Fam", "1Fam", "1Fam",…
FALSE $ HouseStyle <chr> "2Story", "1Story", "2Story", "2Story", "2Story", "1.5F…
FALSE $ OverallQual <int> 7, 6, 7, 7, 8, 5, 8, 7, 5, 5, 9, 7, 7, 5, 5, 8, 7, 8, 5…
FALSE $ OverallCond <int> 5, 8, 5, 5, 5, 5, 5, 5, 6, 5, 5, 5, 8, 5, 6, 5, 7, 5, 7…
FALSE $ YearBuilt <int> 2003, 1976, 2001, 1915, 2000, 1993, 2004, 1931, 1939, 1…
FALSE $ YearRemodAdd <int> 2003, 1976, 2002, 1970, 2000, 1995, 2005, 1950, 1950, 1…
FALSE $ RoofStyle <chr> "Gable", "Gable", "Gable", "Gable", "Gable", "Gable", "…
FALSE $ RoofMatl <chr> "CompShg", "CompShg", "CompShg", "CompShg", "CompShg", …
FALSE $ Exterior1st <chr> "VinylSd", "MetalSd", "VinylSd", "Wd Sdng", "VinylSd", …
FALSE $ Exterior2nd <chr> "VinylSd", "MetalSd", "VinylSd", "Wd Shng", "VinylSd", …
FALSE $ MasVnrType <chr> "BrkFace", "None", "BrkFace", "None", "BrkFace", "None"…
FALSE $ MasVnrArea <int> 196, 0, 162, 0, 350, 0, 186, 0, 0, 0, 286, 306, 0, 0, 0…
FALSE $ ExterQual <chr> "Gd", "TA", "Gd", "TA", "Gd", "TA", "Gd", "TA", "TA", "…
FALSE $ ExterCond <chr> "TA", "TA", "TA", "TA", "TA", "TA", "TA", "TA", "TA", "…
FALSE $ Foundation <chr> "PConc", "CBlock", "PConc", "BrkTil", "PConc", "Wood", …
FALSE $ BsmtQual <chr> "Gd", "Gd", "Gd", "TA", "Gd", "Gd", "Ex", "TA", "TA", "…
FALSE $ BsmtCond <chr> "TA", "TA", "TA", "Gd", "TA", "TA", "TA", "TA", "TA", "…
FALSE $ BsmtExposure <chr> "No", "Gd", "Mn", "No", "Av", "No", "Av", "No", "No", "…
FALSE $ BsmtFinType1 <chr> "GLQ", "ALQ", "GLQ", "ALQ", "GLQ", "GLQ", "GLQ", "Unf",…
FALSE $ BsmtFinSF1 <int> 706, 978, 486, 216, 655, 732, 1369, 0, 851, 906, 998, 0…
FALSE $ BsmtFinType2 <chr> "Unf", "Unf", "Unf", "Unf", "Unf", "Unf", "Unf", "Unf",…
FALSE $ BsmtFinSF2 <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0…
FALSE $ BsmtUnfSF <int> 150, 284, 434, 540, 490, 64, 317, 952, 140, 134, 177, 1…
FALSE $ TotalBsmtSF <int> 856, 1262, 920, 756, 1145, 796, 1686, 952, 991, 1040, 1…
FALSE $ Heating <chr> "GasA", "GasA", "GasA", "GasA", "GasA", "GasA", "GasA",…
FALSE $ HeatingQC <chr> "Ex", "Ex", "Ex", "Gd", "Ex", "Ex", "Ex", "Gd", "Ex", "…
FALSE $ CentralAir <chr> "Y", "Y", "Y", "Y", "Y", "Y", "Y", "Y", "Y", "Y", "Y", …
FALSE $ Electrical <chr> "SBrkr", "SBrkr", "SBrkr", "SBrkr", "SBrkr", "SBrkr", "…
FALSE $ X1stFlrSF <int> 856, 1262, 920, 961, 1145, 796, 1694, 1022, 1077, 1040,…
FALSE $ X2ndFlrSF <int> 854, 0, 866, 756, 1053, 566, 0, 752, 0, 0, 1142, 0, 0, …
FALSE $ LowQualFinSF <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0…
FALSE $ GrLivArea <int> 1710, 1262, 1786, 1717, 2198, 1362, 1694, 1774, 1077, 1…
FALSE $ BsmtFullBath <int> 1, 0, 1, 1, 1, 1, 1, 0, 1, 1, 1, 0, 0, 1, 0, 0, 0, 0, 1…
FALSE $ BsmtHalfBath <int> 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0…
FALSE $ FullBath <int> 2, 2, 2, 1, 2, 1, 2, 2, 1, 1, 3, 2, 1, 1, 1, 3, 1, 2, 1…
FALSE $ HalfBath <int> 1, 0, 1, 0, 1, 1, 0, 0, 0, 0, 0, 0, 0, 1, 0, 1, 0, 0, 0…
FALSE $ BedroomAbvGr <int> 3, 3, 3, 3, 4, 1, 3, 2, 2, 3, 4, 3, 2, 3, 3, 4, 3, 3, 3…
FALSE $ KitchenAbvGr <int> 1, 1, 1, 1, 1, 1, 1, 2, 2, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1…
FALSE $ KitchenQual <chr> "Gd", "TA", "Gd", "Gd", "Gd", "TA", "Gd", "TA", "TA", "…
FALSE $ TotRmsAbvGrd <int> 8, 6, 6, 7, 9, 5, 7, 8, 5, 5, 11, 7, 5, 6, 6, 9, 6, 7, …
FALSE $ Functional <chr> "Typ", "Typ", "Typ", "Typ", "Typ", "Typ", "Typ", "Min1"…
FALSE $ Fireplaces <int> 0, 1, 1, 1, 1, 0, 1, 2, 2, 0, 2, 1, 0, 0, 0, 1, 1, 1, 1…
FALSE $ FireplaceQu <chr> "none", "TA", "TA", "Gd", "TA", "none", "Gd", "TA", "TA…
FALSE $ GarageType <chr> "Attchd", "Attchd", "Attchd", "Detchd", "Attchd", "Attc…
FALSE $ GarageYrBlt <int> 2003, 1976, 2001, 1998, 2000, 1993, 2004, 1931, 1939, 1…
FALSE $ GarageFinish <chr> "RFn", "RFn", "RFn", "Unf", "RFn", "Unf", "RFn", "Unf",…
FALSE $ GarageCars <int> 2, 2, 2, 3, 3, 2, 2, 2, 1, 1, 3, 3, 2, 2, 1, 3, 1, 2, 2…
FALSE $ GarageArea <int> 548, 460, 608, 642, 836, 480, 636, 468, 205, 384, 736, …
FALSE $ GarageQual <chr> "TA", "TA", "TA", "TA", "TA", "TA", "TA", "Fa", "Gd", "…
FALSE $ GarageCond <chr> "TA", "TA", "TA", "TA", "TA", "TA", "TA", "TA", "TA", "…
FALSE $ PavedDrive <chr> "Y", "Y", "Y", "Y", "Y", "Y", "Y", "Y", "Y", "Y", "Y", …
FALSE $ WoodDeckSF <int> 0, 298, 0, 0, 192, 40, 255, 90, 0, 0, 147, 160, 48, 0, …
FALSE $ OpenPorchSF <int> 61, 0, 42, 35, 84, 30, 57, 0, 4, 0, 21, 33, 112, 102, 0…
FALSE $ EnclosedPorch <int> 0, 0, 0, 272, 0, 0, 0, 205, 0, 0, 0, 0, 0, 0, 0, 0, 205…
FALSE $ X3SsnPorch <int> 0, 0, 0, 0, 0, 320, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,…
FALSE $ ScreenPorch <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0…
FALSE $ PoolArea <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0…
FALSE $ PoolQC <chr> "none", "none", "none", "none", "none", "none", "none",…
FALSE $ Fence <chr> "none", "none", "none", "none", "none", "MnPrv", "none"…
FALSE $ MiscFeature <chr> "none", "none", "none", "none", "none", "Shed", "none",…
FALSE $ MiscVal <int> 0, 0, 0, 0, 0, 700, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,…
FALSE $ MoSold <int> 2, 5, 9, 2, 12, 10, 8, 4, 1, 2, 7, 8, 7, 6, 5, 11, 6, 9…
FALSE $ YrSold <int> 2008, 2007, 2008, 2006, 2008, 2009, 2007, 2008, 2008, 2…
FALSE $ SaleType <chr> "WD", "WD", "WD", "WD", "WD", "WD", "WD", "WD", "WD", "…
FALSE $ SaleCondition <chr> "Normal", "Normal", "Normal", "Abnorml", "Normal", "Nor…
FALSE $ SalePrice <int> 208500, 181500, 223500, 140000, 250000, 143000, 307000,…houseprice %>% summarise_all(~ sum(is.na(.))) #to check no na valuesFALSE Id MSSubClass MSZoning LotFrontage LotArea Street Alley LotShape LandContour
FALSE 1 0 0 0 0 0 0 0 0 0
FALSE Utilities LotConfig LandSlope Neighborhood Condition1 Condition2 BldgType
FALSE 1 0 0 0 0 0 0 0
FALSE HouseStyle OverallQual OverallCond YearBuilt YearRemodAdd RoofStyle RoofMatl
FALSE 1 0 0 0 0 0 0 0
FALSE Exterior1st Exterior2nd MasVnrType MasVnrArea ExterQual ExterCond Foundation
FALSE 1 0 0 0 0 0 0 0
FALSE BsmtQual BsmtCond BsmtExposure BsmtFinType1 BsmtFinSF1 BsmtFinType2
FALSE 1 0 0 0 0 0 0
FALSE BsmtFinSF2 BsmtUnfSF TotalBsmtSF Heating HeatingQC CentralAir Electrical
FALSE 1 0 0 0 0 0 0 0
FALSE X1stFlrSF X2ndFlrSF LowQualFinSF GrLivArea BsmtFullBath BsmtHalfBath FullBath
FALSE 1 0 0 0 0 0 0 0
FALSE HalfBath BedroomAbvGr KitchenAbvGr KitchenQual TotRmsAbvGrd Functional
FALSE 1 0 0 0 0 0 0
FALSE Fireplaces FireplaceQu GarageType GarageYrBlt GarageFinish GarageCars
FALSE 1 0 0 0 0 0 0
FALSE GarageArea GarageQual GarageCond PavedDrive WoodDeckSF OpenPorchSF
FALSE 1 0 0 0 0 0 0
FALSE EnclosedPorch X3SsnPorch ScreenPorch PoolArea PoolQC Fence MiscFeature
FALSE 1 0 0 0 0 0 0 0
FALSE MiscVal MoSold YrSold SaleType SaleCondition SalePrice
FALSE 1 0 0 0 0 0 0# Utilities column gives no information since it is same for all the instances
unique(houseprice$Utilities)FALSE [1] "AllPub"#So we drop it
houseprice$Utilities <- NULL
var_mode <- sapply(houseprice, mode)
var_class <- sapply(houseprice, class)
# find character columns
ind1 <- which(var_mode %in% c("character"))
# character to factor
houseprice[ind1] <- lapply(houseprice[ind1], as.factor)
set.seed(5)
train_indices_house = createDataPartition(y = 1:nrow(houseprice), p = 0.8, list = F)
house_train <- houseprice[train_indices_house, ] %>% as.data.frame()
house_test <- houseprice[-train_indices_house, ] %>% as.data.frame()
X_house_train <- house_train %>% select(-SalePrice)
y_house_train <- house_train %>% select(SalePrice)
X_house_test <- house_test %>% select(-SalePrice)
y_house_test <- house_test %>% select(SalePrice)House prices data set is taken from an ongoing Kaggle competition where the participants try to predict sale prices according to the given information about the houses. It is a regression problem with more than 50 features, some are numerical and some are categorical. It has separate train and test data. Since it is an ongoing competition, the test data does not have the actual values of sale prices. Therefore, I will use the train data as my only data and extract test data from it.
Train data has 1460 instances and using summarise_all function on the train data, we see that the NA values are as follows: 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0.
Since I do not want to lose information, I vill create a new category called “None” for the categorical missing values Alley, FireplaceQu, PoolQC, Fence and MiscFeature. For the other na values, I will remove the instances.
8.1 Penalized Regression Approach
set.seed(50)
house_model_pra<-cv.glmnet(data.matrix(X_house_train), house_train$SalePrice, type.measure = "mse", family="gaussian", nfolds=10)
print(house_model_pra)FALSE
FALSE Call: cv.glmnet(x = data.matrix(X_house_train), y = house_train$SalePrice, type.measure = "mse", nfolds = 10, family = "gaussian")
FALSE
FALSE Measure: Mean-Squared Error
FALSE
FALSE Lambda Measure SE Nonzero
FALSE min 208 1.704e+09 388558177 71
FALSE 1se 10331 2.073e+09 380458145 9house_pred_pra <- predict(house_model_pra, newx=data.matrix(X_house_test) ,s=c("lambda.min"), type = "response")
mape_test_pra <- MLmetrics::MAPE(house_pred_pra, house_test$SalePrice)
mape_train_pra <- MLmetrics::MAPE(predict(house_model_pra, newx=data.matrix(X_house_train) ,s=c("lambda.min"), type = "response"), house_train$SalePrice)
print(house_model_pra)FALSE
FALSE Call: cv.glmnet(x = data.matrix(X_house_train), y = house_train$SalePrice, type.measure = "mse", nfolds = 10, family = "gaussian")
FALSE
FALSE Measure: Mean-Squared Error
FALSE
FALSE Lambda Measure SE Nonzero
FALSE min 208 1.704e+09 388558177 71
FALSE 1se 10331 2.073e+09 380458145 9plot(house_model_pra)
Lambda value that minimizes the error is 207.5690562. Mean absolute percentage error (MAPE) of the model on test data is 12.82%. On the other hand, MAPE on the train data is 10.7%. MAPE is higher on the test data.
8.2 Decision Tree
set.seed(50)
complexity=c(0.01,0.03,0.05)
min_number_of_observations=c(10,15,25)
cv_summary=data.table()
index <- generateCVRuns(house_train$SalePrice, ntimes=1, nfold=10, stratified =TRUE)
best_mape <- 1000
best_complexity <- 0
best_min_number_of_observations <- 0
for(i in 1:10){
for(c in complexity){
for(m in min_number_of_observations){
test_index <- index$`Run 1`[i]
train_cv <- house_train[-test_index[[1]],]
test_cv <- house_train[test_index[[1]],]
fitControl <- rpart.control(cp = c, minbucket = m)
house_model_DT=rpart(SalePrice~., house_train, control=fitControl)
pred <- predict(house_model_DT,test_cv)
mape_DT <- MLmetrics::MAPE(pred, test_cv$SalePrice)
cv_summary <- rbind(cv_summary,data.table(Fold=i, Complexity=c, Min_Observations_per_Leaf=m, MAPE=mape_DT))
if(best_mape>mape_DT){
best_mape<-mape_DT
best_complexity <- c
best_min_number_of_observations <- m
best_pred <- pred
}
}
}
}
cv_summary <- cv_summary %>% group_by(Complexity, Min_Observations_per_Leaf) %>% summarise(mean_MAPE = mean(MAPE))
print(cv_summary)FALSE # A tibble: 9 x 3
FALSE # Groups: Complexity [3]
FALSE Complexity Min_Observations_per_Leaf mean_MAPE
FALSE <dbl> <dbl> <dbl>
FALSE 1 0.01 10 0.147
FALSE 2 0.01 15 0.147
FALSE 3 0.01 25 0.148
FALSE 4 0.03 10 0.196
FALSE 5 0.03 15 0.196
FALSE 6 0.03 25 0.196
FALSE 7 0.05 10 0.196
FALSE 8 0.05 15 0.196
FALSE 9 0.05 25 0.196fitControl <- rpart.control(cp = best_complexity, minbucket = best_min_number_of_observations)
house_model_DT=rpart(SalePrice~., house_train, control=fitControl)
house_pred_DT <- predict(house_model_DT, X_house_test)
train_pred <- predict(house_model_DT, X_house_train)
mape_test_DT <- MLmetrics::MAPE(house_pred_DT, house_test$SalePrice)
mape_train_DT <- MLmetrics::MAPE(predict(house_model_DT,X_house_train), house_train$SalePrice)According to the cross validation result, the best complexity parameter is 0.01 and the best minimal number of observations per tree leaf is 15 since it gives the lowest mape of the model on cross validation test data with 0.111149.
Now I can fit the final model and make predictions on the actual test data. MAPE of the final model on test data is 14.78% and on the train data it is 15.08%. They are very close to each other, which indicates a good model.
8.3 Random Forest
set.seed(50)
m=c(4,6,8)
#for balanced folds: stratified =TRUE
# ntree is J=500 given in hw
J=500
# the minimal number of observations per tree leaf is nodesize 5 given in hw
min_size_of_terminal_nodes=5
n_folds=10
fitControl=trainControl(method = "repeatedcv", number = n_folds, repeats = 1, verboseIter = TRUE, returnResamp = "all")
rf_grid=expand.grid(mtry = m)
house_model_RF=train(SalePrice ~ ., data = house_train, method = "rf", trControl = fitControl, ntree=J, nodesize = min_size_of_terminal_nodes, tuneGrid = rf_grid)FALSE + Fold01.Rep1: mtry=4
FALSE - Fold01.Rep1: mtry=4
FALSE + Fold01.Rep1: mtry=6
FALSE - Fold01.Rep1: mtry=6
FALSE + Fold01.Rep1: mtry=8
FALSE - Fold01.Rep1: mtry=8
FALSE + Fold02.Rep1: mtry=4
FALSE - Fold02.Rep1: mtry=4
FALSE + Fold02.Rep1: mtry=6
FALSE - Fold02.Rep1: mtry=6
FALSE + Fold02.Rep1: mtry=8
FALSE - Fold02.Rep1: mtry=8
FALSE + Fold03.Rep1: mtry=4
FALSE - Fold03.Rep1: mtry=4
FALSE + Fold03.Rep1: mtry=6
FALSE - Fold03.Rep1: mtry=6
FALSE + Fold03.Rep1: mtry=8
FALSE - Fold03.Rep1: mtry=8
FALSE + Fold04.Rep1: mtry=4
FALSE - Fold04.Rep1: mtry=4
FALSE + Fold04.Rep1: mtry=6
FALSE - Fold04.Rep1: mtry=6
FALSE + Fold04.Rep1: mtry=8
FALSE - Fold04.Rep1: mtry=8
FALSE + Fold05.Rep1: mtry=4
FALSE - Fold05.Rep1: mtry=4
FALSE + Fold05.Rep1: mtry=6
FALSE - Fold05.Rep1: mtry=6
FALSE + Fold05.Rep1: mtry=8
FALSE - Fold05.Rep1: mtry=8
FALSE + Fold06.Rep1: mtry=4
FALSE - Fold06.Rep1: mtry=4
FALSE + Fold06.Rep1: mtry=6
FALSE - Fold06.Rep1: mtry=6
FALSE + Fold06.Rep1: mtry=8
FALSE - Fold06.Rep1: mtry=8
FALSE + Fold07.Rep1: mtry=4
FALSE - Fold07.Rep1: mtry=4
FALSE + Fold07.Rep1: mtry=6
FALSE - Fold07.Rep1: mtry=6
FALSE + Fold07.Rep1: mtry=8
FALSE - Fold07.Rep1: mtry=8
FALSE + Fold08.Rep1: mtry=4
FALSE - Fold08.Rep1: mtry=4
FALSE + Fold08.Rep1: mtry=6
FALSE - Fold08.Rep1: mtry=6
FALSE + Fold08.Rep1: mtry=8
FALSE - Fold08.Rep1: mtry=8
FALSE + Fold09.Rep1: mtry=4
FALSE - Fold09.Rep1: mtry=4
FALSE + Fold09.Rep1: mtry=6
FALSE - Fold09.Rep1: mtry=6
FALSE + Fold09.Rep1: mtry=8
FALSE - Fold09.Rep1: mtry=8
FALSE + Fold10.Rep1: mtry=4
FALSE - Fold10.Rep1: mtry=4
FALSE + Fold10.Rep1: mtry=6
FALSE - Fold10.Rep1: mtry=6
FALSE + Fold10.Rep1: mtry=8
FALSE - Fold10.Rep1: mtry=8
FALSE Aggregating results
FALSE Selecting tuning parameters
FALSE Fitting mtry = 8 on full training setprint(house_model_RF[["results"]])FALSE mtry RMSE Rsquared MAE RMSESD RsquaredSD MAESD
FALSE 1 4 40624.42 0.8172214 24374.87 15772.15 0.07984813 5519.292
FALSE 2 6 37703.41 0.8318682 22159.44 14986.48 0.07889018 5043.560
FALSE 3 8 36225.18 0.8409230 21039.76 14620.85 0.07647966 4878.379house_pred_RF <- predict(house_model_RF, newdata=X_house_test)
mape_test_RF <- MLmetrics::MAPE(house_pred_RF, house_test$SalePrice)
mape_train_RF <- MLmetrics::MAPE(predict(house_model_RF,X_house_train), house_train$SalePrice)The best m parameter is found as 8 and the final model is created using it. MAPE of the final model on test data is 8.88% and on the train data it is 6.26%. They are close to each other, which indicates a good model.
8.4 Stochastic Gradient Boosting
set.seed(50)
n_folds <- 10
fitControl <- trainControl(method = "cv", number = n_folds)
gbmGrid_forever <- expand.grid(interaction.depth = c(1, 5, 7), n.trees = c(50,100,200), shrinkage = c(0.05, 0.1, 0.2), n.minobsinnode = c(10))
print(gbmGrid_forever)FALSE interaction.depth n.trees shrinkage n.minobsinnode
FALSE 1 1 50 0.05 10
FALSE 2 5 50 0.05 10
FALSE 3 7 50 0.05 10
FALSE 4 1 100 0.05 10
FALSE 5 5 100 0.05 10
FALSE 6 7 100 0.05 10
FALSE 7 1 200 0.05 10
FALSE 8 5 200 0.05 10
FALSE 9 7 200 0.05 10
FALSE 10 1 50 0.10 10
FALSE 11 5 50 0.10 10
FALSE 12 7 50 0.10 10
FALSE 13 1 100 0.10 10
FALSE 14 5 100 0.10 10
FALSE 15 7 100 0.10 10
FALSE 16 1 200 0.10 10
FALSE 17 5 200 0.10 10
FALSE 18 7 200 0.10 10
FALSE 19 1 50 0.20 10
FALSE 20 5 50 0.20 10
FALSE 21 7 50 0.20 10
FALSE 22 1 100 0.20 10
FALSE 23 5 100 0.20 10
FALSE 24 7 100 0.20 10
FALSE 25 1 200 0.20 10
FALSE 26 5 200 0.20 10
FALSE 27 7 200 0.20 10gbmGrid <- expand.grid(interaction.depth = c(5), n.trees = c(100), shrinkage = c(0.1), n.minobsinnode = c(10))
house_model_sgb=train(SalePrice ~ ., data = data.frame(house_train), method = "gbm", trControl = fitControl, tuneGrid = gbmGrid_forever, verbose=F, preProc = "zv")
house_pred_sgb <- predict(house_model_sgb, newdata=X_house_test) %>% as.matrix()
mape_test_sgb <- MLmetrics::MAPE(house_pred_sgb, house_test$SalePrice)
mape_train_sgb <- MLmetrics::MAPE(predict(house_model_sgb,X_house_train), house_train$SalePrice)
print(house_model_sgb[["results"]])FALSE shrinkage interaction.depth n.minobsinnode n.trees RMSE Rsquared
FALSE 1 0.05 1 10 50 44521.42 0.7826962
FALSE 10 0.10 1 10 50 38408.08 0.8032589
FALSE 19 0.20 1 10 50 37445.01 0.8072553
FALSE 4 0.05 5 10 50 36099.63 0.8351669
FALSE 13 0.10 5 10 50 32962.68 0.8510153
FALSE 22 0.20 5 10 50 34409.15 0.8378721
FALSE 7 0.05 7 10 50 34990.46 0.8439563
FALSE 16 0.10 7 10 50 32806.52 0.8517628
FALSE 25 0.20 7 10 50 34706.08 0.8330642
FALSE 2 0.05 1 10 100 37708.59 0.8117490
FALSE 11 0.10 1 10 100 35963.62 0.8204663
FALSE 20 0.20 1 10 100 35910.31 0.8223209
FALSE 5 0.05 5 10 100 33289.27 0.8488720
FALSE 14 0.10 5 10 100 32976.26 0.8512610
FALSE 23 0.20 5 10 100 33402.05 0.8479358
FALSE 8 0.05 7 10 100 32578.09 0.8554044
FALSE 17 0.10 7 10 100 33108.73 0.8492098
FALSE 26 0.20 7 10 100 34865.35 0.8324573
FALSE 3 0.05 1 10 200 35357.52 0.8283620
FALSE 12 0.10 1 10 200 35152.83 0.8308586
FALSE 21 0.20 1 10 200 36569.91 0.8154462
FALSE 6 0.05 5 10 200 32973.84 0.8521257
FALSE 15 0.10 5 10 200 33089.51 0.8515491
FALSE 24 0.20 5 10 200 33480.09 0.8489250
FALSE 9 0.05 7 10 200 32212.29 0.8573041
FALSE 18 0.10 7 10 200 33447.93 0.8470034
FALSE 27 0.20 7 10 200 35357.16 0.8285757
FALSE MAE RMSESD RsquaredSD MAESD
FALSE 1 29711.56 13121.65 0.06813512 4834.704
FALSE 10 24497.17 12446.17 0.08192922 4570.988
FALSE 19 23371.05 11535.03 0.07317592 3898.590
FALSE 4 22069.63 12809.76 0.07876319 4648.553
FALSE 13 19179.76 12349.68 0.08287159 3968.324
FALSE 22 19726.90 12427.48 0.07677772 3564.959
FALSE 7 20649.17 13209.61 0.08137265 4297.298
FALSE 16 19074.60 12043.05 0.07299341 3773.001
FALSE 25 19903.76 13460.45 0.07692152 3171.009
FALSE 2 24113.06 12197.80 0.07975538 4525.645
FALSE 11 22268.44 11702.81 0.08218408 3924.894
FALSE 20 21940.38 10450.86 0.06738165 3175.232
FALSE 5 19064.55 12693.97 0.08113335 4006.936
FALSE 14 18528.62 12668.10 0.08311735 3640.495
FALSE 23 19589.87 10974.88 0.06234621 3040.750
FALSE 8 18162.49 13137.52 0.08425215 3658.773
FALSE 17 18659.77 12162.98 0.07683851 3614.880
FALSE 26 19619.13 13157.76 0.07689512 3023.074
FALSE 3 22016.43 11723.06 0.07566280 3878.104
FALSE 12 21522.30 11508.42 0.07385661 3754.427
FALSE 21 21993.53 10850.82 0.06774750 3074.557
FALSE 6 18257.35 13101.04 0.07898225 3930.895
FALSE 15 18329.00 12495.87 0.07646508 3530.505
FALSE 24 19598.82 10922.22 0.05935772 3196.555
FALSE 9 17484.25 12831.47 0.08418385 3399.600
FALSE 18 18781.07 11909.28 0.07562016 3494.744
FALSE 27 20213.15 12970.35 0.07481806 3159.616Using gbmGrid_forever and keeping the minimal number of observations per tree leaf fixed as 10, optimal parameters are 5, 0.1 and 100 respective to depth, learning rate and number of trees.
MAPE of the final model on test data is 7.8% and on the train data it is 6.38%. They are close to each other, which indicates a good model. These values are the smallest ones compared to other methods.
9 Comparison
For the observations below, I will use the model order PRA, DT, RF and SGB.
For Optical Recognition of Handwritten Digits dataset, the accuracy of the models on the test data are: 95.6%, 74.23%, 97.5% and 95.83%, respective to the model order. RF and SGB methods indicates overfitting since they both have 100% accuracy on train data. DT does not indicate overfitting but since accuracy is near 75% for both test and train data, it may ve underfitting. On the other hand, PRA seem to work better on this data since it’s test and train accuracy values are close and high, also train accuracy is not 100%.
For Mushroom Types dataset, the accuracy of the models on the test data are: 100%, 99.32%, 100% and 100%, respective to the model order. These values are very high, so there might be a data leakage problem.
For Customer Churn dataset, the accuracy of the models on the test data are: 75.43%, 73.15%, 75.36% and 75.57%, respective to the model order.
For House Prices dataset, the MAPE of the models on the train data are: 10.7%, 15.08%, 6.26% and 6.38%, respective to the model order. For House Prices dataset, the MAPE of the models on the test data are: 12.82%, 14.78%, 8.88% and 7.8%, respective to the model order. Here, DT model performed better on the test data. On the other hand SGB model gave the smallest MAPE for both train and test data, since the errors are very close, this is a good model.
10 References
[3] Churn Data