Setup
Primero cargamos todos lo paquetes que vamos a utilizar
library(Epi)
library(DescTools)
library(plyr)
library(dplyr)
##
## Attaching package: 'dplyr'
## The following objects are masked from 'package:plyr':
##
## arrange, count, desc, failwith, id, mutate, rename, summarise,
## summarize
## The following objects are masked from 'package:stats':
##
## filter, lag
## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, union
library(Hmisc)
## Loading required package: lattice
## Loading required package: survival
## Loading required package: Formula
## Loading required package: ggplot2
##
## Attaching package: 'Hmisc'
## The following objects are masked from 'package:dplyr':
##
## src, summarize
## The following objects are masked from 'package:plyr':
##
## is.discrete, summarize
## The following objects are masked from 'package:DescTools':
##
## %nin%, Label, Mean, Quantile
## The following objects are masked from 'package:base':
##
## format.pval, units
library(Formula)
library(rms)
## Loading required package: SparseM
##
## Attaching package: 'SparseM'
## The following object is masked from 'package:base':
##
## backsolve
library(ggplot2)
library(carData)
library(effects)
## Use the command
## lattice::trellis.par.set(effectsTheme())
## to customize lattice options for effects plots.
## See ?effectsTheme for details.
Importing data
dat1 <- read.csv("C:\\Users\\ASPIRE\\Documents\\Trabajos\\usmp\\UCI_Credit_Card.csv")
Data exploration and preparation A continuacion vamos familiarizandonos con la data e identificando inconsistencias, usaremos las funciones head para revisar las 6 primeras files de todas las columnas, str para ver los labels, summary nos da una idea general de la estadistica descriptiva de los valores de cada variable. Para una evaluacion mas detallada de cada variable usamos la funcion Desc.
Observamos principalmente que: las variables no estan correctamente codificadas por ejemplo EDUCATION es interger y deberia estar codificado como factor. Esto es principalmente importante en la elaboracion de los modelos de regression.
head(dat1)
str(dat1)
## 'data.frame': 30000 obs. of 25 variables:
## $ ID : int 1 2 3 4 5 6 7 8 9 10 ...
## $ LIMIT_BAL : num 20000 120000 90000 50000 50000 50000 500000 100000 140000 20000 ...
## $ SEX : int 2 2 2 2 1 1 1 2 2 1 ...
## $ EDUCATION : int 2 2 2 2 2 1 1 2 3 3 ...
## $ MARRIAGE : int 1 2 2 1 1 2 2 2 1 2 ...
## $ AGE : int 24 26 34 37 57 37 29 23 28 35 ...
## $ PAY_0 : int 2 -1 0 0 -1 0 0 0 0 -2 ...
## $ PAY_2 : int 2 2 0 0 0 0 0 -1 0 -2 ...
## $ PAY_3 : int -1 0 0 0 -1 0 0 -1 2 -2 ...
## $ PAY_4 : int -1 0 0 0 0 0 0 0 0 -2 ...
## $ PAY_5 : int -2 0 0 0 0 0 0 0 0 -1 ...
## $ PAY_6 : int -2 2 0 0 0 0 0 -1 0 -1 ...
## $ BILL_AMT1 : num 3913 2682 29239 46990 8617 ...
## $ BILL_AMT2 : num 3102 1725 14027 48233 5670 ...
## $ BILL_AMT3 : num 689 2682 13559 49291 35835 ...
## $ BILL_AMT4 : num 0 3272 14331 28314 20940 ...
## $ BILL_AMT5 : num 0 3455 14948 28959 19146 ...
## $ BILL_AMT6 : num 0 3261 15549 29547 19131 ...
## $ PAY_AMT1 : num 0 0 1518 2000 2000 ...
## $ PAY_AMT2 : num 689 1000 1500 2019 36681 ...
## $ PAY_AMT3 : num 0 1000 1000 1200 10000 657 38000 0 432 0 ...
## $ PAY_AMT4 : num 0 1000 1000 1100 9000 ...
## $ PAY_AMT5 : num 0 0 1000 1069 689 ...
## $ PAY_AMT6 : num 0 2000 5000 1000 679 ...
## $ default.payment.next.month: int 1 1 0 0 0 0 0 0 0 0 ...
summary(dat1)
## ID LIMIT_BAL SEX EDUCATION
## Min. : 1 Min. : 10000 Min. :1.000 Min. :0.000
## 1st Qu.: 7501 1st Qu.: 50000 1st Qu.:1.000 1st Qu.:1.000
## Median :15000 Median : 140000 Median :2.000 Median :2.000
## Mean :15000 Mean : 167484 Mean :1.604 Mean :1.853
## 3rd Qu.:22500 3rd Qu.: 240000 3rd Qu.:2.000 3rd Qu.:2.000
## Max. :30000 Max. :1000000 Max. :2.000 Max. :6.000
## MARRIAGE AGE PAY_0 PAY_2
## Min. :0.000 Min. :21.00 Min. :-2.0000 Min. :-2.0000
## 1st Qu.:1.000 1st Qu.:28.00 1st Qu.:-1.0000 1st Qu.:-1.0000
## Median :2.000 Median :34.00 Median : 0.0000 Median : 0.0000
## Mean :1.552 Mean :35.49 Mean :-0.0167 Mean :-0.1338
## 3rd Qu.:2.000 3rd Qu.:41.00 3rd Qu.: 0.0000 3rd Qu.: 0.0000
## Max. :3.000 Max. :79.00 Max. : 8.0000 Max. : 8.0000
## PAY_3 PAY_4 PAY_5 PAY_6
## Min. :-2.0000 Min. :-2.0000 Min. :-2.0000 Min. :-2.0000
## 1st Qu.:-1.0000 1st Qu.:-1.0000 1st Qu.:-1.0000 1st Qu.:-1.0000
## Median : 0.0000 Median : 0.0000 Median : 0.0000 Median : 0.0000
## Mean :-0.1662 Mean :-0.2207 Mean :-0.2662 Mean :-0.2911
## 3rd Qu.: 0.0000 3rd Qu.: 0.0000 3rd Qu.: 0.0000 3rd Qu.: 0.0000
## Max. : 8.0000 Max. : 8.0000 Max. : 8.0000 Max. : 8.0000
## BILL_AMT1 BILL_AMT2 BILL_AMT3 BILL_AMT4
## Min. :-165580 Min. :-69777 Min. :-157264 Min. :-170000
## 1st Qu.: 3559 1st Qu.: 2985 1st Qu.: 2666 1st Qu.: 2327
## Median : 22382 Median : 21200 Median : 20089 Median : 19052
## Mean : 51223 Mean : 49179 Mean : 47013 Mean : 43263
## 3rd Qu.: 67091 3rd Qu.: 64006 3rd Qu.: 60165 3rd Qu.: 54506
## Max. : 964511 Max. :983931 Max. :1664089 Max. : 891586
## BILL_AMT5 BILL_AMT6 PAY_AMT1 PAY_AMT2
## Min. :-81334 Min. :-339603 Min. : 0 Min. : 0
## 1st Qu.: 1763 1st Qu.: 1256 1st Qu.: 1000 1st Qu.: 833
## Median : 18105 Median : 17071 Median : 2100 Median : 2009
## Mean : 40311 Mean : 38872 Mean : 5664 Mean : 5921
## 3rd Qu.: 50191 3rd Qu.: 49198 3rd Qu.: 5006 3rd Qu.: 5000
## Max. :927171 Max. : 961664 Max. :873552 Max. :1684259
## PAY_AMT3 PAY_AMT4 PAY_AMT5 PAY_AMT6
## Min. : 0 Min. : 0 Min. : 0.0 Min. : 0.0
## 1st Qu.: 390 1st Qu.: 296 1st Qu.: 252.5 1st Qu.: 117.8
## Median : 1800 Median : 1500 Median : 1500.0 Median : 1500.0
## Mean : 5226 Mean : 4826 Mean : 4799.4 Mean : 5215.5
## 3rd Qu.: 4505 3rd Qu.: 4013 3rd Qu.: 4031.5 3rd Qu.: 4000.0
## Max. :896040 Max. :621000 Max. :426529.0 Max. :528666.0
## default.payment.next.month
## Min. :0.0000
## 1st Qu.:0.0000
## Median :0.0000
## Mean :0.2212
## 3rd Qu.:0.0000
## Max. :1.0000
Desc(dat1$SEX)
## ------------------------------------------------------------------------------
## dat1$SEX (integer - dichotomous)
##
## length n NAs unique
## 30'000 30'000 0 2
## 100.0% 0.0%
##
## freq perc lci.95 uci.95'
## 1 11'888 39.6% 39.1% 40.2%
## 2 18'112 60.4% 59.8% 60.9%
##
## ' 95%-CI (Wilson)
Desc(dat1$EDUCATION)
## ------------------------------------------------------------------------------
## dat1$EDUCATION (integer)
##
## length n NAs unique 0s mean meanCI'
## 30'000 30'000 0 7 14 1.85 1.84
## 100.0% 0.0% 0.0% 1.86
##
## .05 .10 .25 median .75 .90 .95
## 1.00 1.00 1.00 2.00 2.00 3.00 3.00
##
## range sd vcoef mad IQR skew kurt
## 6.00 0.79 0.43 1.48 1.00 0.97 2.08
##
##
## level freq perc cumfreq cumperc
## 1 0 14 0.0% 14 0.0%
## 2 1 10'585 35.3% 10'599 35.3%
## 3 2 14'030 46.8% 24'629 82.1%
## 4 3 4'917 16.4% 29'546 98.5%
## 5 4 123 0.4% 29'669 98.9%
## 6 5 280 0.9% 29'949 99.8%
## 7 6 51 0.2% 30'000 100.0%
##
## ' 95%-CI (classic)
Desc(dat1$LIMIT_BAL)
## ------------------------------------------------------------------------------
## dat1$LIMIT_BAL (numeric)
##
## length n NAs unique 0s mean'
## 30'000 30'000 0 81 0 167'484.32
## 100.0% 0.0% 0.0%
##
## .05 .10 .25 median .75 .90
## 20'000.00 30'000.00 50'000.00 140'000.00 240'000.00 360'000.00
##
## range sd vcoef mad IQR skew
## 990'000.00 129'747.66 0.77 133'434.00 190'000.00 0.99
##
## meanCI
## 166'016.06
## 168'952.59
##
## .95
## 430'000.00
##
## kurt
## 0.54
##
## lowest : 10'000.0 (493), 16'000.0 (2), 20'000.0 (1'976), 30'000.0 (1'610), 40'000.0 (230)
## highest: 750'000.0 (4), 760'000.0, 780'000.0 (2), 800'000.0 (2), 1'000'000.0
##
## heap(?): remarkable frequency (11.2%) for the mode(s) (= 50000)
##
## ' 95%-CI (classic)
Desc(dat1$MARRIAGE)
## ------------------------------------------------------------------------------
## dat1$MARRIAGE (integer)
##
## length n NAs unique 0s mean meanCI'
## 30'000 30'000 0 4 54 1.55 1.55
## 100.0% 0.0% 0.2% 1.56
##
## .05 .10 .25 median .75 .90 .95
## 1.00 1.00 1.00 2.00 2.00 2.00 2.00
##
## range sd vcoef mad IQR skew kurt
## 3.00 0.52 0.34 0.00 1.00 -0.02 -1.36
##
##
## level freq perc cumfreq cumperc
## 1 0 54 0.2% 54 0.2%
## 2 1 13'659 45.5% 13'713 45.7%
## 3 2 15'964 53.2% 29'677 98.9%
## 4 3 323 1.1% 30'000 100.0%
##
## ' 95%-CI (classic)
Desc(dat1$AGE)
## ------------------------------------------------------------------------------
## dat1$AGE (integer)
##
## length n NAs unique 0s mean meanCI'
## 30'000 30'000 0 56 0 35.49 35.38
## 100.0% 0.0% 0.0% 35.59
##
## .05 .10 .25 median .75 .90 .95
## 23.00 25.00 28.00 34.00 41.00 49.00 53.00
##
## range sd vcoef mad IQR skew kurt
## 58.00 9.22 0.26 8.90 13.00 0.73 0.04
##
## lowest : 21 (67), 22 (560), 23 (931), 24 (1'127), 25 (1'186)
## highest: 72 (3), 73 (4), 74, 75 (3), 79
##
## heap(?): remarkable frequency (5.3%) for the mode(s) (= 29)
##
## ' 95%-CI (classic)
str(dat1$EDUCATION)
## int [1:30000] 2 2 2 2 2 1 1 2 3 3 ...
str(dat1$MARRIAGE)
## int [1:30000] 1 2 2 1 1 2 2 2 1 2 ...
str(dat1$SEX)
## int [1:30000] 2 2 2 2 1 1 1 2 2 1 ...
str(dat1$LIMIT_BAL)
## num [1:30000] 20000 120000 90000 50000 50000 50000 500000 100000 140000 20000 ...
str(dat1$AGE)
## int [1:30000] 24 26 34 37 57 37 29 23 28 35 ...
Recodificacion de las variables Para una mejor visualizacion de los resultados y debido a que EDUCATION, SEX, MARRIAGE tiene valores igual “0” que corresponde a missing values se realiza un recodificacion de las variables. Ademas SEX, EDUCATION, MARRIGE, AGE no estan codificados correctamente.
dat1$SEX <- as.factor(dat1$SEX)
dat1$EDUCATION <- as.factor(dat1$EDUCATION)
dat1$MARRIAGE <- as.factor(dat1$MARRIAGE)
dat1$AGE <- as.numeric(dat1$AGE)
dat1 <- dat1[complete.cases(dat1),]
dat1 <- dat1 %>% mutate(EDUCATION = case_when(EDUCATION == '1' ~ 'graduate school', EDUCATION == '2' ~ 'university', EDUCATION == '3' ~ 'high school', EDUCATION == '4' ~ 'others'))
dat1 <- dat1 %>% mutate(SEX = case_when(SEX == '1' ~ 'male', SEX == '2' ~ 'female'))
dat1 <- dat1 %>% mutate(MARRIAGE = case_when(MARRIAGE == '1' ~ 'married', MARRIAGE == '2' ~ 'single',
MARRIAGE == '3' ~ 'others'))
Pregunta 1 Haciendo uso del conjunto de datos UCI_Credit_Card para estimar un modelo de regresión logit eligiendo. Donde la variable dependiente es: default.payment.next.month
Seleccionamos las 5 primeras variables para la elaboracion del presente trabajo
fit1 <- glm(default.payment.next.month ~ LIMIT_BAL + SEX + EDUCATION + MARRIAGE + AGE, data = dat1, family = binomial)
Pregunta 2 Estime e interprete el efecto de sus parámetros sobre el ratio de odds.
summary(fit1)
##
## Call:
## glm(formula = default.payment.next.month ~ LIMIT_BAL + SEX +
## EDUCATION + MARRIAGE + AGE, family = binomial, data = dat1)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -0.9978 -0.7693 -0.6494 -0.4251 2.5552
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -8.408e-01 7.740e-02 -10.864 < 2e-16 ***
## LIMIT_BAL -3.366e-06 1.342e-07 -25.088 < 2e-16 ***
## SEXmale 1.720e-01 2.898e-02 5.933 2.98e-09 ***
## EDUCATIONhigh school 8.396e-03 4.471e-02 0.188 0.851025
## EDUCATIONothers -1.350e+00 3.912e-01 -3.452 0.000557 ***
## EDUCATIONuniversity 3.339e-02 3.341e-02 0.999 0.317607
## MARRIAGEothers -1.271e-01 1.305e-01 -0.974 0.330019
## MARRIAGEsingle -2.075e-01 3.275e-02 -6.335 2.37e-10 ***
## AGE 3.722e-03 1.763e-03 2.112 0.034716 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 31427 on 29600 degrees of freedom
## Residual deviance: 30524 on 29592 degrees of freedom
## (399 observations deleted due to missingness)
## AIC: 30542
##
## Number of Fisher Scoring iterations: 5
exp(coef(fit1))
## (Intercept) LIMIT_BAL SEXmale
## 0.4313440 0.9999966 1.1876248
## EDUCATIONhigh school EDUCATIONothers EDUCATIONuniversity
## 1.0084315 0.2591690 1.0339581
## MARRIAGEothers MARRIAGEsingle AGE
## 0.8806090 0.8126289 1.0037292
Interpretacion de los OR:
Intercept:
El intercepto no tiene una interpretacion muy util. Matematicamente es definido como el valor de la variable dependiente (Y) cuando las variables independientes son igual a 0.
En nuestro model se podria interpretar de la siguiente manera:
Para una persona con LIMIT_BAL = 0, SEX = female, EDUCATION = graduate school), MARRIAGE = married, AGE = 0, el odds ratio de default.payment.next.month es 0.4313440. Este coeficiente es estaditicamente significativo.
LIMIT_BAL: Por un incremento de 1 unidad de LIMIT_BAL de una persona, el odd ratio es 0.9999966, manteniendo las demas variables a un valor fijo. Este coeficiente es estaditicamente significativo.
En simples palabras: si solo tomo en cuenta el monto de credito y mantengo las otras variables iguales, a major monto de credito, la probabilidad de pago el siguiente mess es menor.
SEXmale: Para un hombre comparado en comparacion una mujer, el odd ratio es 1.1876248, manteniendo las demas variables a un valor fijo. Este coeficiente es estaditicamente significativo.
En simples palabras: si solo tomo en cuenta el sexo de la persona que saca el credito y mantengo las otras variables iguales, los hombres comparados con las mujeres tienen mas probabilidad de realizar un pago el siguiente mes.
EDUCATIONhigh school: Para un persona con educacion tipo high school comparado una educacion tipo graduate school el odd ratio es 1.0084315, manteniendo las demas variables a un valor fijo.
En simples palabras: si solo tomo en cuenta el tipo de educacion de la persona que saca el credito, y mantengo las otras variables iguales, las personas con educacion secudaria comparadas con las peronas con education tecnica tienen mayor probabilidad de realizar un pago el siguiente mes.
EDUCATIONother: Para un persona con educacion tipo other comparado una educacion tipo graduate school el odd ratio es 0.2591690, manteniendo las demas variables a un valor fijo. Este coeficiente es estaditicamente significativo.
En simples palabras: si solo tomo en cuenta el tipo de educacion de la persona que saca el credito, y mantengo las otras variables iguales, las personas con educacion tipo “others” comparadas con las peronas con educacion tecnica tienen menor probabilidad de realizar un pago el siguiente mes.
EDUCATIONuniversity: Para un persona con educacion tipo university comparado una educacion tipo graduate school el odd ratio es 1.0339581, manteniendo las demas variables a un valor fijo.
En simples palabras: si solo tomo en cuenta el tipo de educacion de la persona que saca el credito, y mantengo las otras variables iguales, las personas con educacion universitaria comparadas con las peronas con educacion tecnica tienen mayor probabilidad de realizar un pago el siguiente mes.
MARRIAGEothers: Para una persona con estatus marital tipo others comparado con estatus marital married el odd ratio es 0.8806090.
En simples palabras: si solo tomo en cuenta el status marital de la persona que saca el credito y mantento las otras variables iguales, las personas con estatus marital tipo “others” comparadas con las casadas, tienen menor probabilidad de realizar un pago el siguiente mes.
MARRIAGEsingle: Para una persona con estatus marital tipo single comparado con estatus marital married el odd ratio es 0.8126289, manteniendo las demas variables a un valor fijo. Este coeficiente es estaditicamente significativo.
En simples palabras: si solo tomo en cuenta el status marital de la persona que saca el credito y mantento las otras variables iguales, las personas solteras comparadas con las casadas, tienen menor probabilidad de realizar un pago el siguiente mes.
AGE: Por un incremento de 1 unidad en la edad de una persona, el odd ratio es 1.0037292, manteniendo las demas variables a un valor fijo. Este coeficiente es estaditicamente significativo.
En simples palabras: si solo tomo en cuenta la edad de la persona que saca el credito y mantengo las otras variables iguales, las personas con mayor edad tienen mayor probabilidad de pago el siguiente mes.
Pregunta 3
Estime la probabilidad de caer en default asumiendo valores específicos para sus covariables.
Basado en el modelo creado: Un persona con LIMIT_BAL=20000, EDUCATION = ‘university’, MARRIAGE = ‘married’, AGE=21 tiene La probabilidad de realizar un pago el sigiente mes de 0.348724.
ND <- with(dat1, expand.grid(LIMIT_BAL=20000,
SEX = 'male', EDUCATION = 'university',
MARRIAGE = 'married', AGE=21))
ND$rankP <- predict(fit1, newdata = ND, type = "response")
ND
Effect plot Para una mejor visualizacion del modelo se realiza un diagrama de efecto con el 95% intervalo de confianza segun sexo y edad. Podemos observar que la probabilidad de pago en el siguiente mes para hombres y para mujeres aumenta con la edad. Los hombres con respecto a mujeres tienen en promedio major probabilidad de pago.
plot(Effect(focal.predictors = c("SEX", 'AGE'), mod = fit1,
xlevels=list(AGE=21:79)))
Pregunta 4
Evalúe la capacidad predictiva del modelo.
Asi 0.046 de la variacion de la occurrencia de un pago en el siguiente mes (default.payment.next.month) es explicado por la variables del modelo.
mod1b <- lrm(default.payment.next.month ~ LIMIT_BAL + SEX + EDUCATION + MARRIAGE + AGE, data = dat1)
print(mod1b)
## Frequencies of Missing Values Due to Each Variable
## default.payment.next.month LIMIT_BAL
## 0 0
## SEX EDUCATION
## 0 345
## MARRIAGE AGE
## 54 0
##
## Logistic Regression Model
##
## lrm(formula = default.payment.next.month ~ LIMIT_BAL + SEX +
## EDUCATION + MARRIAGE + AGE, data = dat1)
##
##
## Model Likelihood Discrimination Rank Discrim.
## Ratio Test Indexes Indexes
## Obs 29601 LR chi2 902.74 R2 0.046 C 0.623
## 0 22996 d.f. 8 g 0.511 Dxy 0.246
## 1 6605 Pr(> chi2) <0.0001 gr 1.667 gamma 0.246
## max |deriv| 0.0007 gp 0.081 tau-a 0.085
## Brier 0.168
##
## Coef S.E. Wald Z Pr(>|Z|)
## Intercept -0.8408 0.0774 -10.86 <0.0001
## LIMIT_BAL 0.0000 0.0000 -25.09 <0.0001
## SEX=male 0.1720 0.0290 5.93 <0.0001
## EDUCATION=high school 0.0084 0.0447 0.19 0.8510
## EDUCATION=others -1.3503 0.3912 -3.45 0.0006
## EDUCATION=university 0.0334 0.0334 1.00 0.3176
## MARRIAGE=others -0.1271 0.1305 -0.97 0.3300
## MARRIAGE=single -0.2075 0.0328 -6.34 <0.0001
## AGE 0.0037 0.0018 2.11 0.0347
##
Pregunta 5 Punto de corte para los valores ajustados.
El punto optimo del modelo que maximiza la sensibilidad y especificidad es 0.238.
prob <- predict(fit1, type='response')
ROC(test = prob, stat = fit1$y)
Pregunta 6 Especificidad / Sensibilidad
basado en el ROC del modelo la sensibilidad es 60.1% y la especificidad es 58.9%.
Pregunta 7 Análisis de la curva ROC
El area bajo la curva del modelo es 0.623