German Credit: Linear Regression Analysis
Sangamesh K S
October 8, 2017
Introduction
we will try to understand the data using linear regression and perform feature selection and other test on the data and try to optimize the model.Lets import the data with the name mydata using choose.file and read.csv function
mydata<-read.csv("C:/Users/Sangmesh/Google Drive/Big Data using R/Kaggle/germancredit.csv",header = TRUE)
Now Lets view the structure of the data
str(mydata)
## 'data.frame': 1000 obs. of 21 variables:
## $ Default : int 0 1 0 0 1 0 0 0 0 1 ...
## $ checkingstatus1: Factor w/ 4 levels "A11","A12","A13",..: 1 2 4 1 1 4 4 2 4 2 ...
## $ duration : int 6 48 12 42 24 36 24 36 12 30 ...
## $ history : Factor w/ 5 levels "A30","A31","A32",..: 5 3 5 3 4 3 3 3 3 5 ...
## $ purpose : Factor w/ 10 levels "A40","A41","A410",..: 5 5 8 4 1 8 4 2 5 1 ...
## $ amount : int 1169 5951 2096 7882 4870 9055 2835 6948 3059 5234 ...
## $ savings : Factor w/ 5 levels "A61","A62","A63",..: 5 1 1 1 1 5 3 1 4 1 ...
## $ employ : Factor w/ 5 levels "A71","A72","A73",..: 5 3 4 4 3 3 5 3 4 1 ...
## $ installment : int 4 2 2 2 3 2 3 2 2 4 ...
## $ status : Factor w/ 4 levels "A91","A92","A93",..: 3 2 3 3 3 3 3 3 1 4 ...
## $ others : Factor w/ 3 levels "A101","A102",..: 1 1 1 3 1 1 1 1 1 1 ...
## $ residence : int 4 2 3 4 4 4 4 2 4 2 ...
## $ property : Factor w/ 4 levels "A121","A122",..: 1 1 1 2 4 4 2 3 1 3 ...
## $ age : int 67 22 49 45 53 35 53 35 61 28 ...
## $ otherplans : Factor w/ 3 levels "A141","A142",..: 3 3 3 3 3 3 3 3 3 3 ...
## $ housing : Factor w/ 3 levels "A151","A152",..: 2 2 2 3 3 3 2 1 2 2 ...
## $ cards : int 2 1 1 1 2 1 1 1 1 2 ...
## $ job : Factor w/ 4 levels "A171","A172",..: 3 3 2 3 3 2 3 4 2 4 ...
## $ liable : int 1 1 2 2 2 2 1 1 1 1 ...
## $ tele : Factor w/ 2 levels "A191","A192": 2 1 1 1 1 2 1 2 1 1 ...
## $ foreign : Factor w/ 2 levels "A201","A202": 1 1 1 1 1 1 1 1 1 1 ...
the structure of the data reveals many things most of the list are integer or character. Now my next quest is there any missing values? Let’s check it using is.na function and View function. As we know the data is a treated data free from missing values and all the data in the lists are converted into coded characters.which help us in computation
Now let me show you a small sample how the data is distributed with a example plotting it on a barograph.
barplot(table(mydata$checkingstatus1),main ="Plot of checking Status")
you can see the data is distributed in 4 coded characters i.e A11 to A14. and with the help of data visualization I have shown you the ##frequency distribution of character
barplot(table(mydata$history))
now let’s see for history in mydata and plot a bar plot if I want to have a check the frequency distribution in tabular format. Please check the same below
table(mydata$history)
##
## A30 A31 A32 A33 A34
## 40 49 530 88 293
Now we came to know about the data. Its important to know what to do with this data? what are the various algorithms we can apply.
I want to start from basic lm function and let me try with 2 or 3 algorithms and test which fit the best among these 2 or 3 algos
1st I am interested applying linear models to understand the data and its significance and co-linearity
Data preparation
There is a small problem applying the linear model. the data are not numerical so we will convert the data into numerical format
So we will re code the vectors into numerical format
Now I will remove “A” form mydata and convert everything into numerical. So we will come up with a new dataset and keep the original mydata undisturbed
newdata<-mydata
After creating a new data
We will look at data wangling of newdata that we have recently
table(newdata$checkingstatus1)
##
## A11 A12 A13 A14
## 274 269 63 394
Now we will remove all A and convert it into a numaric format
newdata$checkingstatus1<-gsub("A11","11",newdata$checkingstatus1)
newdata$checkingstatus1<-gsub("A12","12",newdata$checkingstatus1)
newdata$checkingstatus1<-gsub("A13","13",newdata$checkingstatus1)
newdata$checkingstatus1<-gsub("A14","14",newdata$checkingstatus1)
now lets check it it done?
table(newdata$checkingstatus1)
##
## 11 12 13 14
## 274 269 63 394
now lets look at history
table(newdata$history)
##
## A30 A31 A32 A33 A34
## 40 49 530 88 293
and perfrom the same
newdata$history<-gsub("A30","30",newdata$history)
newdata$history<-gsub("A31","31",newdata$history)
newdata$history<-gsub("A32","32",newdata$history)
newdata$history<-gsub("A33","33",newdata$history)
newdata$history<-gsub("A34","34",newdata$history)
table(newdata$history)
##
## 30 31 32 33 34
## 40 49 530 88 293
repete the same process for purpose
table(newdata$purpose)
##
## A40 A41 A410 A42 A43 A44 A45 A46 A48 A49
## 234 103 12 181 280 12 22 50 9 97
newdata$purpose<-gsub("A40","40",newdata$purpose)
newdata$purpose<-gsub("A41","41",newdata$purpose)
newdata$purpose<-gsub("A410","410",newdata$purpose)
newdata$purpose<-gsub("A42","42",newdata$purpose)
newdata$purpose<-gsub("A44","44",newdata$purpose)
newdata$purpose<-gsub("A45","45",newdata$purpose)
newdata$purpose<-gsub("A46","46",newdata$purpose)
newdata$purpose<-gsub("A48","48",newdata$purpose)
newdata$purpose<-gsub("A49","49",newdata$purpose)
newdata$purpose<-gsub("A43","43",newdata$purpose)
table(newdata$purpose)
##
## 40 41 410 42 43 44 45 46 48 49
## 234 103 12 181 280 12 22 50 9 97
Lets look at saving
table(newdata$savings)
##
## A61 A62 A63 A64 A65
## 603 103 63 48 183
newdata$savings<-gsub("A61","61",newdata$savings)
newdata$savings<-gsub("A62","62",newdata$savings)
newdata$savings<-gsub("A63","63",newdata$savings)
newdata$savings<-gsub("A64","64",newdata$savings)
newdata$savings<-gsub("A65","65",newdata$savings)
table(newdata$savings)
##
## 61 62 63 64 65
## 603 103 63 48 183
Lets look ar employ
table(newdata$employ)
##
## A71 A72 A73 A74 A75
## 62 172 339 174 253
newdata$employ<-gsub("A71","71",newdata$employ)
newdata$employ<-gsub("A72","72",newdata$employ)
newdata$employ<-gsub("A73","73",newdata$employ)
newdata$employ<-gsub("A74","74",newdata$employ)
newdata$employ<-gsub("A75","75",newdata$employ)
table(newdata$employ)
##
## 71 72 73 74 75
## 62 172 339 174 253
Perform the same with all the data list
table(newdata$status)
##
## A91 A92 A93 A94
## 50 310 548 92
newdata$status<-gsub("A91","91",newdata$status)
newdata$status<-gsub("A92","92",newdata$status)
newdata$status<-gsub("A93","93",newdata$status)
newdata$status<-gsub("A94","94",newdata$status)
table(newdata$status)
##
## 91 92 93 94
## 50 310 548 92
table(newdata$others)
##
## A101 A102 A103
## 907 41 52
newdata$others<-gsub("A101","101",newdata$others)
newdata$others<-gsub("A102","102",newdata$others)
newdata$others<-gsub("A103","103",newdata$others)
table(newdata$others)
##
## 101 102 103
## 907 41 52
table(newdata$property)
##
## A121 A122 A123 A124
## 282 232 332 154
newdata$property<-gsub("A121","121",newdata$property)
newdata$property<-gsub("A122","122",newdata$property)
newdata$property<-gsub("A123","123",newdata$property)
newdata$property<-gsub("A124","124",newdata$property)
table(newdata$property)
##
## 121 122 123 124
## 282 232 332 154
table(newdata$otherplans)
##
## A141 A142 A143
## 139 47 814
newdata$otherplans<-gsub("A141","141",newdata$otherplans)
newdata$otherplans<-gsub("A142","142",newdata$otherplans)
newdata$otherplans<-gsub("A143","143",newdata$otherplans)
table(newdata$otherplans)
##
## 141 142 143
## 139 47 814
table(newdata$housing)
##
## A151 A152 A153
## 179 713 108
newdata$housing<-gsub("A151","151",newdata$housing)
newdata$housing<-gsub("A152","152",newdata$housing)
newdata$housing<-gsub("A153","153",newdata$housing)
table(newdata$housing)
##
## 151 152 153
## 179 713 108
table(newdata$job)
##
## A171 A172 A173 A174
## 22 200 630 148
newdata$job<-gsub("A171","171",newdata$job)
newdata$job<-gsub("A172","172",newdata$job)
newdata$job<-gsub("A173","173",newdata$job)
newdata$job<-gsub("A174","174",newdata$job)
table(newdata$job)
##
## 171 172 173 174
## 22 200 630 148
table(newdata$tele)
##
## A191 A192
## 596 404
newdata$tele<-gsub("A191","191",newdata$tele)
newdata$tele<-gsub("A192","192",newdata$tele)
table(newdata$tele)
##
## 191 192
## 596 404
table(newdata$foreign)
##
## A201 A202
## 963 37
newdata$foreign<-gsub("A201","201",newdata$foreign)
newdata$foreign<-gsub("A202","202",newdata$foreign)
table(newdata$foreign)
##
## 201 202
## 963 37
After converting character to numeric the r is not going to read it as numeric so we will convert it to numeric
newdata$foreign<-as.numeric(newdata$foreign)
newdata$checkingstatus1<-as.numeric(newdata$checkingstatus1)
newdata$duration<-as.numeric(newdata$duration)
newdata$Default<-as.numeric(newdata$Default)
newdata$history<-as.numeric(newdata$history)
newdata$purpose<-as.numeric(newdata$purpose)
newdata$amount<-as.numeric(newdata$amount)
newdata$savings<-as.numeric(newdata$savings)
newdata$employ<-as.numeric(newdata$employ)
newdata$installment<-as.numeric(newdata$installment)
newdata$status<-as.numeric(newdata$status)
newdata$others<-as.numeric(newdata$others)
newdata$residence<-as.numeric(newdata$residence)
newdata$property<-as.numeric(newdata$property)
newdata$age<-as.numeric(newdata$age)
newdata$otherplans<-as.numeric(newdata$otherplans)
newdata$housing<-as.numeric(newdata$housing)
newdata$cards<-as.numeric(newdata$cards)
newdata$job<-as.numeric(newdata$job)
newdata$liable<-as.numeric(newdata$liable)
newdata$tele<-as.numeric(newdata$tele)
newdata$foreign<-as.numeric(newdata$foreign)
str(newdata)
## 'data.frame': 1000 obs. of 21 variables:
## $ Default : num 0 1 0 0 1 0 0 0 0 1 ...
## $ checkingstatus1: num 11 12 14 11 11 14 14 12 14 12 ...
## $ duration : num 6 48 12 42 24 36 24 36 12 30 ...
## $ history : num 34 32 34 32 33 32 32 32 32 34 ...
## $ purpose : num 43 43 46 42 40 46 42 41 43 40 ...
## $ amount : num 1169 5951 2096 7882 4870 ...
## $ savings : num 65 61 61 61 61 65 63 61 64 61 ...
## $ employ : num 75 73 74 74 73 73 75 73 74 71 ...
## $ installment : num 4 2 2 2 3 2 3 2 2 4 ...
## $ status : num 93 92 93 93 93 93 93 93 91 94 ...
## $ others : num 101 101 101 103 101 101 101 101 101 101 ...
## $ residence : num 4 2 3 4 4 4 4 2 4 2 ...
## $ property : num 121 121 121 122 124 124 122 123 121 123 ...
## $ age : num 67 22 49 45 53 35 53 35 61 28 ...
## $ otherplans : num 143 143 143 143 143 143 143 143 143 143 ...
## $ housing : num 152 152 152 153 153 153 152 151 152 152 ...
## $ cards : num 2 1 1 1 2 1 1 1 1 2 ...
## $ job : num 173 173 172 173 173 172 173 174 172 174 ...
## $ liable : num 1 1 2 2 2 2 1 1 1 1 ...
## $ tele : num 192 191 191 191 191 192 191 192 191 191 ...
## $ foreign : num 201 201 201 201 201 201 201 201 201 201 ...
Let me apply linear model on to the data
fit_lm_1<-lm(Default~.,data = newdata)
summary(fit_lm_1)
##
## Call:
## lm(formula = Default ~ ., data = newdata)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.8203 -0.2974 -0.1150 0.3458 1.0573
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 6.078e+01 1.651e+01 3.680 0.000246 ***
## checkingstatus1 -9.972e-02 1.082e-02 -9.216 < 2e-16 ***
## duration 4.187e-03 1.471e-03 2.846 0.004518 **
## history -6.426e-02 1.366e-02 -4.705 2.91e-06 ***
## purpose -3.074e-04 3.321e-04 -0.926 0.354880
## amount 1.631e-05 6.850e-06 2.381 0.017455 *
## savings -3.405e-02 8.473e-03 -4.019 6.29e-05 ***
## employ -2.649e-02 1.157e-02 -2.290 0.022216 *
## installment 4.704e-02 1.307e-02 3.600 0.000335 ***
## status -4.465e-02 1.864e-02 -2.395 0.016824 *
## others -5.962e-02 2.777e-02 -2.147 0.032025 *
## residence 3.507e-03 1.258e-02 0.279 0.780390
## property 3.285e-02 1.434e-02 2.291 0.022169 *
## age -1.094e-03 1.283e-03 -0.853 0.393976
## otherplans -4.883e-02 1.873e-02 -2.607 0.009286 **
## housing -4.866e-02 2.753e-02 -1.767 0.077463 .
## cards 4.076e-02 2.518e-02 1.618 0.105894
## job -3.702e-03 2.257e-02 -0.164 0.869744
## liable 2.895e-02 3.667e-02 0.789 0.430048
## tele -5.265e-02 2.941e-02 -1.790 0.073747 .
## foreign -1.103e-01 6.998e-02 -1.576 0.115267
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.4038 on 979 degrees of freedom
## Multiple R-squared: 0.24, Adjusted R-squared: 0.2244
## F-statistic: 15.45 on 20 and 979 DF, p-value: < 2.2e-16
Now lets check the VIF colinearity
library(car)
## Warning: package 'car' was built under R version 3.4.1
vif(fit_lm_1)
## checkingstatus1 duration history purpose
## 1.134732 1.928471 1.341004 1.086295
## amount savings employ installment
## 2.290836 1.098126 1.197026 1.309721
## status others residence property
## 1.067936 1.078030 1.180578 1.389551
## age otherplans housing cards
## 1.305363 1.070639 1.311005 1.296714
## job liable tele foreign
## 1.333466 1.080431 1.277645 1.070195
There is no co-linearity in the model
By applying linear model we are able to identify that checkingstatus,history,savings and installment are highly significant
Duration and otherplans have moderate significance
Others are not hiving significance or have least significance
Lets load leaps packages to carry further analysis
library(leaps)
## Warning: package 'leaps' was built under R version 3.4.1
Now we will do subset and later we will do feature selection and try to find out the min and max features for the model.
We will use regsubset to create a subset
sub.fit<-regsubsets(Default~.,data = newdata)
now from the subset we will create summary to slelect the best fit
best.summary<-summary(sub.fit)
Now we will use which.min fuction to understand the min features to be involved in the model with the help of RSS
which.min(best.summary$rss)
## [1] 8
Now with the help of data data visuvilization let me show performence of Mallow’s Cp
which.min(best.summary$cp)
## [1] 8
Lets do some data Viz stuff
par(mfrow=c(1,2))
plot(best.summary$cp,xlab = "Number of Features",ylab = "Mallow's CP")
plot(sub.fit,scale = "Cp")
Let me try the same with BIC
which.min(best.summary$bic)
## [1] 5
Data visuvalisation of BIC
par(mfrow=c(1,2))
plot(best.summary$bic,xlab = "Number of Features",ylab = "BIC")
plot(sub.fit,scale = "bic")
let me plot adjusted R^2
plot(best.summary$adjr2)
Now my main objective for the best model feature selection is to select the model with highest adjt R^2 and low Cp
There is no need for stepwise i.e forward and backward test for optimizing feature selection we are clearly able to see checking status, history,duration,savings,installment, status, others and otherplans are significantly influencing the model
thus we have done possible feature selection. now we will create a model based on the analysis of feature selection
fit.lm.2<-lm(Default~checkingstatus1+duration+history+savings+installment+status+others+otherplans,data=newdata)
summary(fit.lm.2)
##
## Call:
## lm(formula = Default ~ checkingstatus1 + duration + history +
## savings + installment + status + others + otherplans, data = newdata)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.7322 -0.3043 -0.1158 0.3621 1.1285
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 24.770992 4.302522 5.757 1.14e-08 ***
## checkingstatus1 -0.103332 0.010801 -9.567 < 2e-16 ***
## duration 0.006787 0.001081 6.281 5.04e-10 ***
## history -0.063290 0.012250 -5.166 2.89e-07 ***
## savings -0.036733 0.008411 -4.367 1.39e-05 ***
## installment 0.030879 0.011649 2.651 0.00816 **
## status -0.050493 0.018404 -2.744 0.00619 **
## others -0.072948 0.027386 -2.664 0.00785 **
## otherplans -0.049104 0.018465 -2.659 0.00796 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.4073 on 991 degrees of freedom
## Multiple R-squared: 0.2172, Adjusted R-squared: 0.2109
## F-statistic: 34.37 on 8 and 991 DF, p-value: < 2.2e-16
Lets again check for the colineraity
vif(fit.lm.2)
## checkingstatus1 duration history savings
## 1.111266 1.022862 1.060301 1.063681
## installment status others otherplans
## 1.022836 1.022688 1.030762 1.022322
The model looks good without significant colinerity now lets do Breush Pagan test on the new linear model that we have selected
library(lmtest)
## Warning: package 'lmtest' was built under R version 3.4.1
## Loading required package: zoo
##
## Attaching package: 'zoo'
## The following objects are masked from 'package:base':
##
## as.Date, as.Date.numeric
bptest(fit.lm.2)
##
## studentized Breusch-Pagan test
##
## data: fit.lm.2
## BP = 96.231, df = 8, p-value < 2.2e-16
we got a more extreme P value we have non constant variance homoscedasticity is one of the assumptions for linear regression so we will drop this and let apply other model.
If you are interested applying linear regression. Kindly do the transformation of variables with log normal and carry further analysis. Based on these facts I am able to know what are the feature affecting the output.
Logistic regression will work fine for the same.So in our next article I will be covering the same.
No comments:
Post a Comment