Split dataset

At first, we’re splitting the dataset into 2 parts, the training & test dataset. We’ll use the training dataset to train our model & the testing dataset to check the performance of the model.

# Split the dataset in two parts
set.seed(1)
inTrain = createDataPartition(bank$Churn, p = .70)[[1]]

# Assign the 70% of observations to training data
training <- bank[inTrain,]

# Assign the remaining 30 % of observations to testing data
testing <- bank[-inTrain,]

Prepare the recipe of data for modelling

Here we’re using the recipes package to apply the same pre-processing steps to training & test data.

# Create the recipe object
recipe_obj <- 
  recipe(Churn ~ ., data = training) %>% 
  step_zv(all_predictors()) %>%       # check any zero variance features
  step_center(all_numeric()) %>% 
  step_scale(all_numeric()) %>%         # scale the numeric features
  prep()

Processing data according the recipe

train_data <- bake(recipe_obj, training)

test_data  <- bake(recipe_obj, testing)

Setting the train controls for modelling

We’re using TrainControl() function to specify some parameters during model training, e.g. the resampling method, number of folds in K-Fold etc.

train_ctr <- trainControl(method = 'cv', 
                          number = 10,
                          classProbs = TRUE,
                          summaryFunction = twoClassSummary)

Develop multiple machine learning models

Develop a logistic regression model

Logistic_model <- train(Churn ~ .,
                        data = train_data,
                        method = 'glm',
                        family = 'binomial',
                        trControl = train_ctr,
                        metric = 'ROC')

Develop a random forest model

rf_model <- train(Churn ~ .,
                  data = train_data,
                  method = 'rf',
                  trControl = train_ctr,
                  tuneLength = 5,
                  metric = 'ROC')

Develop an XGBoost model

xgb_model <- train(Churn ~ ., data = train_data,
                        method = 'xgbTree',
                        trControl = train_ctr,
                        tuneLength = 5,
                        metric = 'ROC')

Now we can save all models in a single file

save(Logistic_model, rf_model, xgb_model, file = "/Users/manos/OneDrive/Projects/R/All_Projects/Churn/data/ml_models_bank.RDA")

We can load all trained models

load(file = "/Users/manos/OneDrive/Projects/R/All_Projects/Churn/data/ml_models_bank.RDA")

Model Comparison

In this step, we’ll compare the models’ accuracy.

model_list <- resamples(list(Logistic = Logistic_model,
                             Random_forest = rf_model,
                             XgBoost = xgb_model))

summary(model_list)

## 
## Call:
## summary.resamples(object = model_list)
## 
## Models: Logistic, Random_forest, XgBoost 
## Number of resamples: 10 
## 
## ROC 
##                    Min.   1st Qu.    Median      Mean   3rd Qu.      Max. NA's
## Logistic      0.7310692 0.7629073 0.7713738 0.7685092 0.7793850 0.7893184    0
## Random_forest 0.8257989 0.8573118 0.8669290 0.8653269 0.8735644 0.8900139    0
## XgBoost       0.8445059 0.8651163 0.8716812 0.8717819 0.8818719 0.8984382    0
## 
## Sens 
##                    Min.   1st Qu.    Median      Mean   3rd Qu.      Max. NA's
## Logistic      0.9408602 0.9573794 0.9605027 0.9612578 0.9654857 0.9820789    0
## Random_forest 0.9425494 0.9546864 0.9631957 0.9608965 0.9681900 0.9802513    0
## XgBoost       0.9498208 0.9609515 0.9641577 0.9648430 0.9677275 0.9784946    0
## 
## Spec 
##                    Min.   1st Qu.    Median      Mean   3rd Qu.      Max. NA's
## Logistic      0.1818182 0.2126465 0.2245395 0.2307003 0.2376761 0.3076923    0
## Random_forest 0.3732394 0.4534128 0.4825175 0.4803211 0.5131119 0.5594406    0
## XgBoost       0.4265734 0.4797843 0.4964789 0.4902098 0.5034965 0.5524476    0

• Based on ROC (AUC value) the best model is XgBoost. The AUC value of the best model (mean of logistic regression) is 0.87. In general, models with an AUC value > 0.7 are considered as useful, depending of course on the context of the problem.

Model evaluation

# Predictions on test data
pred_logistic <- predict(Logistic_model, newdata = test_data, type = 'prob')

pred_rf <- predict(rf_model, newdata = test_data, type = 'prob')

pred_xgb <- predict(xgb_model, newdata = test_data, type = 'prob')


evaluation_tbl <- tibble(true_class     = test_data$Churn,
                         logistic_churn = pred_logistic$Yes,
                         rf_churn       = pred_rf$Yes,
                         xgb_churn      = pred_xgb$Yes)

evaluation_tbl

## # A tibble: 2,999 x 4
##    true_class logistic_churn rf_churn xgb_churn
##    <fct>               <dbl>    <dbl>     <dbl>
##  1 Yes                0.125     0.276    0.297 
##  2 No                 0.234     0.178    0.0503
##  3 Yes                0.289     0.758    0.997 
##  4 No                 0.112     0.158    0.0631
##  5 No                 0.0600    0.002    0.0139
##  6 No                 0.160     0.006    0.0324
##  7 No                 0.0991    0.002    0.0349
##  8 No                 0.0669    0.042    0.0267
##  9 No                 0.232     0.164    0.246 
## 10 No                 0.0311    0.056    0.0181
## # … with 2,989 more rows

Roc curve

ROC curve or receiver operating characteristic curve is a plot that illustrates the diagnostic ability of a binary classifier system as its discrimination threshold is varied.

i.e. how well the model predicts at different thresholds

# set the second level as the positive class
options(yardstick.event_first = FALSE)

# creating data for ploting ROC curve
roc_curve_logistic <- roc_curve(evaluation_tbl, true_class, logistic_churn) %>% 
  mutate(model = 'logistic')

roc_curve_rf <- roc_curve(evaluation_tbl, true_class, rf_churn) %>% 
  mutate(model = 'RF')

roc_curve_xgb <- roc_curve(evaluation_tbl, true_class, xgb_churn) %>% 
  mutate(model = 'XGB')

# combine all the roc curve data
roc_curve_combine_tbl <- Reduce(rbind, list(roc_curve_logistic, roc_curve_rf,
                                            roc_curve_xgb))

head(roc_curve_combine_tbl,10)

## # A tibble: 10 x 4
##    .threshold specificity sensitivity model   
##         <dbl>       <dbl>       <dbl> <chr>   
##  1  -Inf         0              1     logistic
##  2     0.0115    0              1     logistic
##  3     0.0118    0.000419       1     logistic
##  4     0.0134    0.000838       1     logistic
##  5     0.0137    0.00126        1     logistic
##  6     0.0141    0.00168        1     logistic
##  7     0.0144    0.00209        1     logistic
##  8     0.0149    0.00251        1     logistic
##  9     0.0152    0.00293        1     logistic
## 10     0.0154    0.00293        0.998 logistic

# Plot ROC curves

roc_curve_combine_tbl %>% 
  ggplot(aes(x = 1- specificity, y = sensitivity, color = model))+
  geom_line(size = 1)+
  geom_abline(linetype = 'dashed')+
  scale_color_tableau()+
  labs(title = 'ROC curve Comparison',
       x = '1 - Specificity',
       y = 'Sensitivity')

The largest the AUC value, the better is the model accuracy.

Feature Importance

It is important to understand how a model works, especially in a case of a churn model. A common measure that is used is the global feature importance. It measures the contribution of each variable towards the churn prediction.

This importance is calculated explicitly for each attribute in the dataset, allowing attributes to be ranked and compared to each other. Importance is calculated for a single decision tree by the amount that each attribute split point improves the performance measure, weighted by the number of observations the node is responsible for. The feature importances are then averaged across all of the decision trees within the model.

ggplot(varImp(xgb_model)) +
  labs(title = "Feature Importance for XGBoost churn model") 

It seems that the top 2 most important feature are:
- Age
- Number of products

Predict churn for new customers

When we want to run predictions on customers we can follow the process below.

new_data <- read_csv(file = "/Users/manos/OneDrive/Projects/R/All_Projects/Churn/data/new_data_bank.csv")

new_data_recipe  <- bake(recipe_obj, new_data)

new_dat_pred <- 
  predict(Logistic_model, newdata = new_data_recipe, type = 'prob') %>% 
  select(Yes) %>% 
  rename(churn_prob = Yes) %>% 
  bind_cols(new_data) %>% 
  mutate(churn_group = ntile(churn_prob, n = 10)) %>% 
  select(churn_prob, churn_group, everything()) %>% 
  mutate(churn_prob = round(churn_prob, 2))

new_dat_pred %>% 
  select(churn_prob, churn_group, Age, NumOfProducts, Geography, HasCrCard, Tenure) %>% 
  head(10) %>% 
  kable() %>%
  kableExtra::kable_styling()

Create a churn risk ranking

Although we developed a model that can predict pretty well if a customer will churn, the model output probabilities are not sufficient in the business context. We need some metric that will be understood & easily used by all stakeholders and remove the complexities of e.g. explaining a threshold to non-technical stakeholder.

So instead of an actual probability, a churn risk ranking would be more useful. So we break up the probabilities variable into 10 churn risk buckets. Now a customer has a churn risk from 1 (lowest probability) to 10 (highest probability) .

Tactics for churn prevention

An initial tactic is to develop different sales offers (or marketing campaigns) for the different churn risk groups.

For example, customers that belong in churn risk groups 10 & 9 have a significantly higher churn risk than for example 1 & 2. So it will be crucial to offer them something more to retain them.