Logistic Regression Using Caret Package

library(data.table)
autoF.dt <- fread("AutoFinanaceData.csv")
attach(autoF.dt)

# dimension
dim(autoF.dt)

[1] 28906    21

• Defaulter Flag

“1” if customer has delayed even once ,“0” otherwise

addmargins(table(autoF.dt$DefaulterFlag),1)

    0     1   Sum 
 8332 20574 28906 

round(prop.table(table(autoF.dt$DefaulterFlag))*100,2)

    0     1 
28.82 71.18 

71.18 % of 28,906 customers have defaulted on their Auto Loan

• TENORYR

• DWNPMFR

• FULLPDC

# data columns
colnames(autoF.dt)

 [1] "Agmt No"        "ContractStatus" "StartDate"      "AGE"           
 [5] "NOOFDEPE"       "MTHINCTH"       "SALDATFR"       "TENORYR"       
 [9] "DWNPMFR"        "PROFBUS"        "QUALHSC"        "QUAL_PG"       
[13] "SEXCODE"        "FULLPDC"        "FRICODE"        "WASHCODE"      
[17] "Region"         "Branch"         "DefaulterFlag"  "DefaulterType" 
[21] "DATASET"       

# structure of the data table
str(autoF.dt)

# descriptive statistics of the dataframe
library(psych)
describe(autoF.dt)[, c(1:5)]

Reserve 80% for training and 20% of test

# loading the package
library(caTools)
# fixing the observations 
set.seed(123)
# splitting the data 
split = sample.split(autoF.dt$DefaulterFlag, SplitRatio = 0.75)
# creating the training set
trainingSet = subset(autoF.dt, split == TRUE)
# creating the test set
testSet = subset(autoF.dt, split == FALSE)

# dimensions of the full data
dim(autoF.dt)

[1] 28906    21

# dimensions of the training data
dim(trainingSet)

[1] 21679    21

# dimensions of the Testing data
dim(testSet)

[1] 7227   21

library(caret)
# control parameters
objControl <- trainControl(method = "boot", 
                           number = 2, 
                           returnResamp = 'none', 
                           summaryFunction = twoClassSummary, 
                           classProbs = TRUE,
                           savePredictions = TRUE)

set.seed(766)
# model building using caret package
caretLogitModel <- train(trainingSet[,c(4:17)],
                      trainingSet$DefaulterFlag,
                      method = 'glmStepAIC',
                      trControl = objControl,
                      metric = "ROC",verbose = FALSE)

# list of all variables
colnames(trainingSet)

 [1] "Agmt No"        "ContractStatus" "StartDate"      "AGE"           
 [5] "NOOFDEPE"       "MTHINCTH"       "SALDATFR"       "TENORYR"       
 [9] "DWNPMFR"        "PROFBUS"        "QUALHSC"        "QUAL_PG"       
[13] "SEXCODE"        "FULLPDC"        "FRICODE"        "WASHCODE"      
[17] "Region"         "Branch"         "DefaulterFlag"  "DefaulterType" 
[21] "DATASET"       

# summary of the model
caretLogitModelAIC <- caretLogitModel$finalModel
summary(caretLogitModelAIC)

• MTHINCTH (monthly Income in Thousand)

# important variables
varImp(caretLogitModelAIC) 

                Overall
AGE            8.732406
NOOFDEPE       5.252202
SALDATFR       9.256242
TENORYR       17.025141
DWNPMFR       10.310381
PROFBUS1       4.029424
QUALHSC1       4.642604
QUAL_PG1       3.805121
SEXCODE1       3.900838
FULLPDC1      33.721322
FRICODE1       4.737867
WASHCODE1      5.555757
RegionAP2      3.221608
RegionChennai 10.042751
RegionKA1      4.623264
RegionKE2      3.966970
RegionTN1      5.932138
RegionTN2      4.211713
RegionVellore  4.117365

# creating single value data frame

newdata <- data.frame(
AGE = mean(trainingSet$AGE),
NOOFDEPE = mean(trainingSet$NOOFDEPE),
SALDATFR = mean(trainingSet$SALDATFR),
TENORYR  = mean(trainingSet$TENORYR),
DWNPMFR  = mean(trainingSet$DWNPMFR),
PROFBUS  = "0",
QUALHSC  = "0",
QUAL_PG  = "0",
SEXCODE  = "0",
FULLPDC  = "0",
FRICODE  = "1",
WASHCODE = "1",
Region   = "TN1")

# predicting probability
pred1 <- predict(caretLogitModelAIC, newdata, type = "response")
pred1

        1 
0.7195374 

The probability of Mr. A. Kumar defaulting = 71.95%

Did not submit FULLPDC

# creating single value data frame

newdata <- data.frame(
AGE = mean(trainingSet$AGE),
NOOFDEPE = mean(trainingSet$NOOFDEPE),
SALDATFR = mean(trainingSet$SALDATFR),
TENORYR  = mean(trainingSet$TENORYR),
DWNPMFR  = mean(trainingSet$DWNPMFR),
PROFBUS  = "0",
QUALHSC  = "0",
QUAL_PG  = "0",
SEXCODE  = "0",
FULLPDC  = "1",
FRICODE  = "1",
WASHCODE = "1",
Region   = "TN1")

# predicting probability
pred2 <- predict(caretLogitModelAIC, newdata, type = "response")
pred2

      1 
0.42682 

The probability of Mr. B. Kumar defaulting = 42.68%

Mr. B. Kumar is similar to A. Kumar, except that he submitted FULLPDC

# table for mean and sd
library(psych)
describe(trainingSet)[c(4:5,7:9),c(3,4)]

          mean   sd
AGE      36.45 9.82
NOOFDEPE  2.85 1.62
SALDATFR  0.44 0.46
TENORYR   1.29 0.52
DWNPMFR   0.38 0.16

Mr. C. Kumar is similar to A. Kumar, expect his loan Tenure is (Mean + SD of Loan Tenure) Years

# creating single value data frame
newdataX <- data.frame(
AGE = mean(trainingSet$AGE),
NOOFDEPE = mean(trainingSet$NOOFDEPE),
SALDATFR = mean(trainingSet$SALDATFR),
TENORYR  = mean(trainingSet$TENORYR) + sd(trainingSet$TENORYR),
DWNPMFR  = mean(trainingSet$DWNPMFR),
PROFBUS  = "0",
QUALHSC  = "0",
QUAL_PG  = "0",
SEXCODE  = "0",
FULLPDC  =  "0",
FRICODE  = "1",
WASHCODE = "1",
Region   = "TN1")

# predicting probability
probX <- predict(caretLogitModelAIC, newdataX, type = "response")
probX

        1 
0.7934194 

            TENORYR PROB  
1 "mean-sd" "0.76"  "0.63"
1 "mean"    "1.28"  "0.72"
1 "mean+sd" "1.8"   "0.79"

Prob(Default) when TENORYR is low (e.g. mean - SD) = 63%

Prob(Default) when TENORYR is average = 72%

Prob(Default) when TENORYR is high (e.g. mean + SD) = 79%

Thus, an increase in TENORYR Increases Prob(Default)

• Classification is done based on the predicted probability of default.

• We could assume the threshold probability above which customers will default to be 50%.

• But that leads to a meaningless answer.

# predicted probabilities
predProbTest <- predict(caretLogitModelAIC, testSet, type = "response")

# plot of probabilities
plot(predProbTest, 
     main = "Scatterplot of Probabilities of default (test data)", 
     xlab = "Customer ID", 
     ylab = "Predicted Probability of default")

# minimum value of predicted probability
min(predProbTest)

[1] 0.1672974

# maximum value of predicted probability
max(predProbTest)

[1] 0.977689

# mean value of predicted probability
mean(predProbTest)

[1] 0.7101482

library(caret)
# confusion matrix using caret package
yPred <- ifelse(predProbTest > 0.5, "Yes", "No")
predY <- as.factor(yPred)
levels(testSet$DefaulterFlag) <- c("No", "Yes")
confusionMatrix(data = predY, reference = testSet$DefaulterFlag, positive = "Yes")

library(caret)
# confusion matrix using caret package
yPred <- ifelse(predProbTest > 0.7101482, "Yes", "No")
predY <- as.factor(yPred)
levels(testSet$DefaulterFlag) <- c("No", "Yes")
confusionMatrix(data = predY, reference = testSet$DefaulterFlag, positive = "Yes")

In some applications of ROC curves, you want the point closest to the TPR of (1) and FPR of (0). This cut point is “optimal” in the sense it weighs both sensitivity and specificity equally. To deterimine this cutoff, you can use the code below. The code takes in BOTH the performance object and prediction object and gives the optimal cutoff value of your predictions.

library(ROCR)
# predicted prob on test data
predProbTest <- predict(caretLogitModelAIC, testSet, type = "response")
# prediction
lgPredObj <- prediction(predProbTest, testSet$DefaulterFlag)
# performance
lgPerfObj <- performance(lgPredObj, "tpr",x.measure = "fpr")
# function for optimal cut point
opt.cut = function(lgPerfObj, lgPredObj){
    cut.ind = mapply(FUN=function(x, y, p){
        d = (x - 0)^2 + (y-1)^2
        ind = which(d == min(d))
        c(sensitivity = y[[ind]], specificity = 1-x[[ind]], 
            cutoff = p[[ind]])
    }, lgPerfObj@x.values, lgPerfObj@y.values, lgPredObj@cutoffs)
}
print(opt.cut(lgPerfObj, lgPredObj))

                 [,1]
sensitivity 0.6940124
specificity 0.6615458
cutoff      0.7121607

library(caret)
# confusion matrix using caret package
yPred <- ifelse(predProbTest > 0.7121607, "Yes", "No")
predY <- as.factor(yPred)
levels(testSet$DefaulterFlag) <- c("No", "Yes")
confusionMatrix(data = predY, reference = testSet$DefaulterFlag, positive = "Yes")

library(ROCR)
# predicted probabilities
predProbTest <- predict(caretLogitModelAIC, testSet, type = "response")
lgPredObj <- prediction(predProbTest,testSet$DefaulterFlag)
lgPerfObj <- performance(lgPredObj, "tpr","fpr")
plot(lgPerfObj,main = "ROC Curve",col = 2,lwd = 2)
abline(a = 0,b = 1,lwd = 2,lty = 3,col = "black")

# predicted probability on whole data
predProbW <- predict(caretLogitModelAIC, autoF.dt, type = "response")

# creating a table with consumer ID and pred probability
tab1 <- cbind(ID = autoF.dt$`Agmt No`,predProbW)
# converting as data frame
tab2 <- as.data.frame(tab1)
predProbW1<- round(as.numeric(predProbW)*100,2)
# ordering the table
tab3 <- tab2[with(tab2,  order(- predProbW1)),]
# few rows of the table 
head(tab3,13)

              ID         predProbW
12641 KA02201217 0.977689003767456
12700 KA02201218  0.97737022872726
2972  TN17100419  0.97702167957542
7779  KA02201297 0.976724012868744
5577  KA02200459 0.976324628602939
11587 AP18100009 0.973341881883995
4816  KA02200700 0.972678056191383
11406 KA02200099  0.97266595289672
10235 KA02201299 0.972491525598801
10744 KA02201231 0.972381407526044
11701 KA02201298 0.971704299085467
9089  KA02201116 0.971427460239132
12320 KA02201114 0.971426280754097