library(readxl)
library(tidyverse)
library(lubridate)
library(dplyr)
library(data.table)
library(caTools)
library(tidyr)
library(ggplot2)
library(ggthemes)
library(caret)
Base <- read_xlsx("C:/Users/matir/OneDrive/Escritorio/Mati/Universidad/2022/Segundo Trimestre/Analitica aplicada a Negocios-Finanzas/Examen/data_ex1.xlsx")
#1- Concéntrese primero en las variables DURATION y AMOUNT. Explore estas variables en el conjunto de datos: ¿Hay valores atípicos y/o faltantes? Si encuentra valores atípicos, elimínelos.
sum(is.na(Base$DURATION)) # NO hay valores faltantes
## [1] 0
ggplot(Base, aes(x=DURATION))+ geom_boxplot() + theme_economist() # Valores atipicos pasado los 60
Base <- filter(Base, DURATION<=60)
sum(is.na(Base$AMOUNT)) # NO hay valores faltantes
## [1] 0
ggplot(Base, aes(x=AMOUNT))+ geom_boxplot() + theme_economist() # Valores atipicos pasado los 15000 app
Base <- filter(Base, AMOUNT<=17000)
#2- ¿Cómo varía la variable RESPONSE con los niveles de DURATION y AMOUNT? Ilustreestas relaciones con gráficos, por ejemplo, gráficos bivariados. Interprete
Base$RESPONSE<- haven::as_factor(Base$RESPONSE)
DURATION_1 <- ggplot(data = Base,
mapping = aes(x = DURATION, fill = RESPONSE, color = RESPONSE))
DURATION_1 + geom_density(alpha = 0.3)+
labs(y = "%",
x = "DURATION",
title = "Relación RESPONSE-DURATION") + theme_classic()
Con la variable duration se puede ver que existen 2 picks donde el ratio
crediticio es categorizado como bueno, donde después del segundo
(cercano a los 23 meses)comienza un declive del buen ratio,teniendo
mayor predominancia los no pagadores. De lo anterior se puede ver que se
tiende a tener buen ratio crediticio con valores de Duration más
bajos.
AMOUNT_1 <- ggplot(data = Base,
mapping = aes(x = AMOUNT, fill = RESPONSE, color = RESPONSE))
AMOUNT_1 + geom_density(alpha = 0.3)+
labs(y = "%",
x = "AMOUNT",
title = "Relación RESPONSE-AMOUNT") + theme_classic()
La variable Amount logra reflejar una predominancia de un buen ratio
crediticio hasta un valor de 2500 del Amount, donde luego se estabilizan
los buenos y malos de manera marginal, con una pequeña preferencia por
los no pagadores de manera sostenida. Lo anterior se interpreta como que
entre mayor sea el monto del crédito, se tiende a tener peor ratio de
crédito.
#3- Produzca análisis gráficos (por ejemplo, histogramas) para las variables que considere más relevantes para el problema que estamos estudiando. Justifique la elección de dichas variables.
Considerando que tenemos dentro de nuestros objetivos saber la probabilidad que un cliente vaya a pagar o no, por lo que las variables a seleccionar deben tener pertinencia con elementos que afecten al hecho de cumplir con su compromiso crediticio o no.
Algunas de estas variables, y las que se consideran más influyentes son:
JOB: Esta variable categórica nos indica la naturaleza del empleo del cliente, por lo que nos puede dar una intuición del tipo de empleo del individuo y su capacidad adquisitiva para poder pagar el crédito.
NUM_CREDITS: El cliente al poseer una gran concentración de créditos en un mismo banco, nos puede dar la intuición de que este tiene concentrado sus préstamos en un mismo banco y no en otros bancos.
AGE: Dado el rango etario de las personas podemos analizar los distintos tipos de crédito solicitados. Por ejemplo, podemos intuir que la gente en sus 20 años tiende a solicitar créditos de consumo, mientras que la gente tiende a optar por su primer crédito hipotecario. Dado el análisis del histograma de esta variable, podemos decir que sus datos se distribuyen con un sesgo positivo a la derecha, una tendencia modal que va para los 30 años.
REAL_ESTATE: Esta variable binaria nos indica si es que el cliente posee o no bienes de raíces. El cliente al poseer una propiedad, nos indica de qué puede hacer pasar esta como activo de garantía frente a una concesión de un crédito.
DURATION: Está variable numérica nos indica en meses el tiempo del crédito otorgado. A través del análisis gráfico podemos ver que tiene un leve sesgo a la derecha, por lo que podemos concluir que la gran mayoría de los créditos son menores a 20 meses.
par(mfcol = c(1, 2))
hist(Base$JOB, main = "JOB")
hist(Base$NUM_CREDITS, main = "NUM_CREDITS")
hist(Base$AGE, main = "AGE")
hist(Base$REAL_ESTATE, main = "REAL_ESTATE")
hist(Base$DURATION, main = "DURATION")
#B- Construcción del modelo 1: LDA (15 puntos)
#4. Divida el conjunto de datos en dos subconjuntos, el conjunto de entrenamiento y el conjunto de prueba, los cuales representan el 60% y el 40% de las observaciones, respectivamente. Establezca una semilla igual a 100.
set.seed(100)
Entrenamiento <- createDataPartition(y = Base$RESPONSE, p = 0.6, list = FALSE, times = 1)
Datos_Entrenamiento <- Base[Entrenamiento,]
Datos_Prueba <- Base[-Entrenamiento,]
#5. Considerando el conjunto de entrenamiento, implemente el LDA en R, con solo dos predictores DURATION y AMOUNT. Siga la metodología aprendida en la magistral, es decir, no use la función LDA.
#+
#6. Con base en los parámetros estimados en el conjunto de entrenamiento, identifique aquellos solicitantes de crédito en el conjunto de prueba que no habrían obtenido el préstamo. Como punto de corte (cut-off), considere que es el punto medio entre los vectores xA y xB, es decir, el cut-off que no considera ni priors ni costos de errores.
LDA <- as.data.frame(cbind(Datos_Entrenamiento$RESPONSE,Datos_Entrenamiento$AMOUNT, Datos_Entrenamiento$DURATION))
names(LDA) = c("Tipo", "X1", "X2")
LDA$Tipo = LDA$Tipo - 1 # corrección de la trasformación
Base1 <- filter(LDA, Tipo == 0)
Base0 <- filter(LDA, Tipo == 1)
cov1Base1 <- cov(Base1[ ,2:3] )
cov2Base1 <- cov(Base0[ ,2:3] )
na <- nrow(Datos_Entrenamiento)
nb <- nrow(Datos_Prueba)
escalarA <- (na-1)/(na+nb-2)
escalarB <- (nb-1)/(na+nb-2)
covarianzaT = (cov1Base1*escalarA) + (cov2Base1*escalarB)
invcovarianzaT = solve(covarianzaT)
#Diferenciales / promedios
Xa1 = mean(Base1$X1)
Xa2 = mean(Base1$X2)
Xb1 = mean(Base0$X1)
Xb2 = mean(Base0$X2)
MXAB12 <- matrix(c(Xa1-Xb1,Xa2-Xb2),nrow = 2, ncol = 1,byrow = TRUE)
#matriz varianzaa covarianza totales
ALPHA <- invcovarianzaT%*%MXAB12
LDA$Zi <- LDA$X1*ALPHA[1,1]+ LDA$X2*ALPHA[2,1] # SCORE (Zi)
Summean_M <- matrix(c(Xa1+Xb1,Xa2+Xb2),nrow = 2, ncol = 1,byrow = TRUE)
#alpha no priors
tALPHA <- t(ALPHA)
alphao_M <- (tALPHA%*%Summean_M)/2
#probabilidad (entendiendo que se pretenden buscar probabiliades, ES NECESARIO USAR LOS PRIORS)
Priors <- matrix(c(na/(na+nb),nb/(na+nb)),nrow = 2, ncol = 1,byrow = TRUE)
LDA$probabilidad <- as.vector(as.vector(1/(1+(Priors[1, 1]/Priors[2, 1])*exp(LDA$Zi-alphao_M))))
#decision otorgamiento en base al SCORE
LDA$Otorgamiento <- as.numeric(LDA$Zi > alphao_M[1, 1])
#7. Con base en los parámetros estimados en el conjunto de entrenamiento, calcule la probabilidad de que un solicitante de crédito sea de alto riesgo en el conjunto de prueba.
Datos_Prueba <- as.data.frame(cbind(Datos_Prueba$RESPONSE,Datos_Prueba$AMOUNT, Datos_Prueba$DURATION))
names(Datos_Prueba) = c("RESPONSE", "AMOUNT", "DURATION")
Datos_Prueba$RESPONSE = Datos_Prueba$RESPONSE - 1
Datos_Prueba$Zi<-Datos_Prueba$AMOUNT*ALPHA[1,1]+Datos_Prueba$DURATION*ALPHA[2,1]
Datos_Prueba$OTORGAMIENTO <- as.numeric(Datos_Prueba$Zi > alphao_M[1, 1])
Datos_Prueba$probabilidad <- as.vector(as.vector(1/(1+(Priors[1, 1]/Priors[2, 1])*exp(Datos_Prueba$Zi-alphao_M))))
FINAL = Datos_Prueba # Creamos una Base igual a Datos de prueba para su manipulación.
Solicitante_de_Malos <- sum(FINAL$OTORGAMIENTO == 0) # no paga
Solicitante_de_Buenos <- sum(FINAL$OTORGAMIENTO == 1) # paga
Probabilidad_alto_riesgo <- Solicitante_de_Malos/(Solicitante_de_Buenos+Solicitante_de_Malos)
print(paste("Probabilidad alto riesgo: ", round(Solicitante_de_Malos/(Solicitante_de_Buenos+Solicitante_de_Malos)*100,2), "%"))
## [1] "Probabilidad alto riesgo: 58.9 %"
#8. Suponga ahora que el banco cambia su regla de decisión y decide que deniega el crédito si la probabilidad que obtiene en el punto 7 es superior a 0,5. Identifique aquellos solicitantes de crédito en el conjunto de prueba que no habrían obtenido el préstamo de acuerdo con esta nueva regla
FINAL$REGLA_0.5 <- ifelse(FINAL$probabilidad>0.5 ,0,1)
#9. Compare el número de solicitantes de crédito en el conjunto de prueba que no habrían obtenido el préstamo de acuerdo con las dos reglas definidas en el punto 6 y 7. Interprete.
FINAL$DOBLE_UMBRAL <- ifelse(((FINAL$probabilidad>alphao_M[1, 1]) & (FINAL$REGLA_0.5==1) | (FINAL$probabilidad>alphao_M[1, 1]) & (FINAL$REGLA_0.5==0)
| (FINAL$probabilidad<alphao_M[1, 1]) | (FINAL$REGLA_0.5==1)),1,0)
#10. Con base al conjunto de entrenamiento, estimar una regresión logística para la probabilidad de que un solicitante de crédito sea de alto riesgo. Considere los mismos dos predictores que para LDA, es decir, DURATION y AMOUNT.
#+
#11. Con base en los parámetros estimados en el conjunto de entrenamiento, calcule la probabilidad de ser un solicitante de alto riesgo en el conjunto de prueba.
LDA_logistic <- as.data.frame(cbind(Datos_Entrenamiento$RESPONSE,Datos_Entrenamiento$AMOUNT, Datos_Entrenamiento$DURATION))
names(LDA_logistic) = c("RESPONSE", "AMOUNT", "DURATION")
LDA_logistic$RESPONSE = LDA_logistic$RESPONSE - 1
#Estimamos
Modelo_logistico <- glm(LDA_logistic$RESPONSE ~ LDA_logistic$AMOUNT+LDA_logistic$DURATION,
data=LDA_logistic, family="binomial")
#Predecimos
LDA_logistic$PREDICCION <- predict(Modelo_logistico, newdata=LDA_logistic, type="response")
#12. Como antes, suponga que el banco toma una regla conservadora y deniega el crédito si la probabilidad que obtiene en el punto anterior es superior a 0,5. Identifique aquellos solicitantes de crédito en el conjunto de prueba que no habrían obtenido el préstamo de acuerdo a esta regla.
LDA_logistic$OBTIENE_CREDITO <- ifelse(LDA_logistic$PREDICCION > 0.50, 0, 1)
#13. Compare las predicciones del modelo de LDA calculadas en el punto 8 y las predicciones de la regresión logística mediante tablas de contingencia. En cada una de las tablas de contingencia (deberá construir dos tablas), considere la variable RESPONSE observada y la predicha por el modelo correspondiente. Analice los resultados. ¿Cuál es el mejor modelo?
tabla <- table(LDA_logistic$OBTIENE_CREDITO,LDA_logistic$RESPONSE)
par(mfrow = c(1, 2))
colores <- c("#80FFFF", "#FFFFFF")
tabla_2 <- prop.table(tabla, margin = 2)
barplot(tabla, col = colores)
barplot(tabla_2, col = colores)
# Volvemos a la ventana gráfica original
par(mfrow = c(1, 1))
##
tabla3 <- table(FINAL$OTORGAMIENTO,FINAL$RESPONSE)
par(mfrow = c(1, 2))
colores <- c("#80FFFF", "#FFFFFF")
tabla_4 <- prop.table(tabla3, margin = 2)
barplot(tabla3, col = colores)
barplot(tabla_4, col = colores)
# Volvemos a la ventana gráfica original
par(mfrow = c(1, 1))
#D- Reconsiderando el modelo 2: Regresión Logística (15 puntos)
#14. Permita que más regresores expliquen potencialmente la probabilidad de tener un alto riesgo crediticio en el modelo de regresión logística. ¿Qué variables considera que son los predictores más importantes de la probabilidad de tener un alto riesgo crediticio? Explique cómo identificó estas variables
Variables elegidas:
HISTORY : El historial crediticio del sujeto es importante, puesto que permite saber si ha solicitado anteriormente un instrumento de este tipo.
SAV_ACCT: El hecho de tener o no ahorros inlfuye en el riesgo, pues pueden servir en caso de emergencia para solventar una deuda.
CHK_ACCT: El nivel de la cuenta corriente de la persona es sin duda el mayor indicador de solvencia ante una deuda.
EMPLOYMENT : Tener o no tener empleo es un factor relevante a la hora de tener una fuente de ingresos mensual que pueda ayudar a pagar un credito.
REAL_ESTATE : En el caso de no tener, la probabilidad de riesgo crediticio aumenta.
La elección de estas 4 variables es principalmente por el principio de solvencia. Asimimso, podemos demostrar que las variables elegidas guardan relacion con “RESPONSE”.
library(corrplot)
library(ggcorrplot)
correlacion <- Base[, c("RESPONSE",
"HISTORY",
"SAV_ACCT",
"EMPLOYMENT",
"REAL_ESTATE",
"CHK_ACCT")]
model.matrix(~0+., data=correlacion) %>%
cor(use="pairwise.complete.obs") %>%
ggcorrplot(show.diag = F, type="lower", lab=TRUE, lab_size=3.5)
LDA_logistic2 <- as.data.frame(cbind(Datos_Entrenamiento$RESPONSE,Datos_Entrenamiento$HISTORY,Datos_Entrenamiento$SAV_ACCT, Datos_Entrenamiento$EMPLOYMENT,Datos_Entrenamiento$REAL_ESTATE, Datos_Entrenamiento$CHK_ACCT))
names(LDA_logistic2) = c("RESPONSE","HISTORY", "SAV_ACCT", "EMPLOYMENT","REAL_ESTATE", "CHK_ACCT")
LDA_logistic2$RESPONSE = LDA_logistic2$RESPONSE - 1
#Estimamos
Modelo_logistico2 <- glm(LDA_logistic2$RESPONSE ~ LDA_logistic2$HISTORY+LDA_logistic2$SAV_ACCT+
LDA_logistic2$EMPLOYMENT+LDA_logistic2$REAL_ESTATE+LDA_logistic2$CHK_ACCT,
data=LDA_logistic2, family="binomial")
#Predecimos
LDA_logistic2$PREDICCION <- predict(Modelo_logistico2, newdata=LDA_logistic2, type="response")
LDA_logistic2$OBTIENE_CREDITO <- ifelse(LDA_logistic2$PREDICCION > 0.50, 0, 1)
Dado que estamos en presencia de una
library(caret)
mejores <- glm(RESPONSE~ .,data=Base,family=binomial())
summary(mejores)
##
## Call:
## glm(formula = RESPONSE ~ ., family = binomial(), data = Base)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -2.6412 -0.7119 0.3913 0.7058 2.3346
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 1.033e+00 8.680e-01 1.191 0.233797
## CHK_ACCT 5.633e-01 7.246e-02 7.774 7.63e-15 ***
## DURATION -2.584e-02 9.152e-03 -2.823 0.004754 **
## HISTORY 3.971e-01 8.983e-02 4.420 9.85e-06 ***
## NEW_CAR -8.213e-01 3.888e-01 -2.113 0.034642 *
## USED_CAR 7.855e-01 4.865e-01 1.615 0.106411
## FURNITURE -6.522e-02 4.027e-01 -0.162 0.871331
## `RADIO/TV` 4.921e-02 3.920e-01 0.126 0.900099
## EDUCATION -8.894e-01 5.030e-01 -1.768 0.077016 .
## RETRAINING -1.189e-01 4.452e-01 -0.267 0.789464
## AMOUNT -1.161e-04 4.323e-05 -2.686 0.007239 **
## SAV_ACCT 2.479e-01 6.061e-02 4.091 4.30e-05 ***
## EMPLOYMENT 1.185e-01 7.468e-02 1.587 0.112550
## INSTALL_RATE -3.253e-01 8.656e-02 -3.758 0.000171 ***
## MALE_DIV -3.526e-01 3.814e-01 -0.925 0.355217
## MALE_SINGLE 5.296e-01 2.053e-01 2.580 0.009877 **
## MALE_MAR_or_WID 1.292e-01 3.073e-01 0.421 0.674067
## `CO-APPLICANT` -3.619e-01 3.992e-01 -0.907 0.364664
## GUARANTOR 9.302e-01 4.194e-01 2.218 0.026561 *
## PRESENT_RESIDENT -1.338e-02 8.401e-02 -0.159 0.873469
## REAL_ESTATE 2.070e-01 2.093e-01 0.989 0.322577
## PROP_UNKN_NONE -5.657e-01 3.726e-01 -1.518 0.128953
## AGE 1.154e-02 8.653e-03 1.333 0.182515
## OTHER_INSTALL -6.177e-01 2.043e-01 -3.024 0.002494 **
## RENT -6.611e-01 4.595e-01 -1.439 0.150171
## OWN_RES -2.400e-01 4.346e-01 -0.552 0.580761
## NUM_CREDITS -2.312e-01 1.660e-01 -1.393 0.163643
## JOB -2.594e-02 1.422e-01 -0.182 0.855260
## NUM_DEPENDENTS -2.554e-01 2.454e-01 -1.041 0.297978
## TELEPHONE 3.514e-01 1.950e-01 1.802 0.071519 .
## FOREIGN 1.553e+00 6.608e-01 2.351 0.018735 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 1216.90 on 997 degrees of freedom
## Residual deviance: 907.83 on 967 degrees of freedom
## AIC: 969.83
##
## Number of Fisher Scoring iterations: 5
mejores_pred <- varImp(mejores, scale = FALSE)
library(tibble)
mejores_pred2 <- tibble::rownames_to_column(mejores_pred, "VARIABLE")
head(mejores_pred2 %>% arrange(desc(Overall)),3)
## VARIABLE Overall
## 1 CHK_ACCT 7.773575
## 2 HISTORY 4.420454
## 3 SAV_ACCT 4.090776
Según los tests de significancia, y gracias a la funcion “varImp”, que en pocas palabras hace un calculo del peso de un preductor según vaya o no en el modelo, las tres variables mas importantes serán: “CHK_ACCT”, “HISTORY”, “SAV_ACCT”.
library(MASS)
table(LDA_logistic$OBTIENE_CREDITO, LDA_logistic$RESPONSE)
##
## 0 1
## 0 161 402
## 1 18 18
prop.table(table(LDA_logistic$OBTIENE_CREDITO, LDA_logistic$RESPONSE))
##
## 0 1
## 0 0.26878130 0.67111853
## 1 0.03005008 0.03005008
table(LDA_logistic2$OBTIENE_CREDITO, LDA_logistic2$RESPONSE)
##
## 0 1
## 0 113 375
## 1 66 45
prop.table(table(LDA_logistic2$OBTIENE_CREDITO, LDA_logistic2$RESPONSE))
##
## 0 1
## 0 0.18864775 0.62604341
## 1 0.11018364 0.07512521
table(LDA_logistic$OBTIENE_CREDITO, LDA_logistic2$OBTIENE_CREDITO)
##
## 0 1
## 0 461 102
## 1 27 9
prop.table(table(LDA_logistic$OBTIENE_CREDITO, LDA_logistic2$OBTIENE_CREDITO))
##
## 0 1
## 0 0.76961603 0.17028381
## 1 0.04507513 0.01502504
library(caret)
conf_matrix<-table(LDA_logistic$OBTIENE_CREDITO, LDA_logistic$RESPONSE)
show(confusionMatrix(conf_matrix))
## Confusion Matrix and Statistics
##
##
## 0 1
## 0 161 402
## 1 18 18
##
## Accuracy : 0.2988
## 95% CI : (0.2624, 0.3373)
## No Information Rate : 0.7012
## P-Value [Acc > NIR] : 1
##
## Kappa : -0.0357
##
## Mcnemar's Test P-Value : <2e-16
##
## Sensitivity : 0.89944
## Specificity : 0.04286
## Pos Pred Value : 0.28597
## Neg Pred Value : 0.50000
## Prevalence : 0.29883
## Detection Rate : 0.26878
## Detection Prevalence : 0.93990
## Balanced Accuracy : 0.47115
##
## 'Positive' Class : 0
##
conf_matrix2<-table(LDA_logistic2$OBTIENE_CREDITO, LDA_logistic2$RESPONSE)
show(confusionMatrix(conf_matrix2))
## Confusion Matrix and Statistics
##
##
## 0 1
## 0 113 375
## 1 66 45
##
## Accuracy : 0.2638
## 95% CI : (0.2289, 0.301)
## No Information Rate : 0.7012
## P-Value [Acc > NIR] : 1
##
## Kappa : -0.1749
##
## Mcnemar's Test P-Value : <2e-16
##
## Sensitivity : 0.6313
## Specificity : 0.1071
## Pos Pred Value : 0.2316
## Neg Pred Value : 0.4054
## Prevalence : 0.2988
## Detection Rate : 0.1886
## Detection Prevalence : 0.8147
## Balanced Accuracy : 0.3692
##
## 'Positive' Class : 0
##
confusion_matrix3<-table(FINAL$OTORGAMIENTO, FINAL$RESPONSE)
show(confusionMatrix(confusion_matrix3))
## Confusion Matrix and Statistics
##
##
## 0 1
## 0 60 175
## 1 59 105
##
## Accuracy : 0.4135
## 95% CI : (0.3648, 0.4636)
## No Information Rate : 0.7018
## P-Value [Acc > NIR] : 1
##
## Kappa : -0.0944
##
## Mcnemar's Test P-Value : 5.571e-14
##
## Sensitivity : 0.5042
## Specificity : 0.3750
## Pos Pred Value : 0.2553
## Neg Pred Value : 0.6402
## Prevalence : 0.2982
## Detection Rate : 0.1504
## Detection Prevalence : 0.5890
## Balanced Accuracy : 0.4396
##
## 'Positive' Class : 0
##
PARA MODELO1:
library(ROCR)
pred1 = predict(Modelo_logistico, newdata=LDA_logistic, type="response")
pred1 = prediction(pred1, LDA_logistic$RESPONSE)
perf1 = performance(pred1, "acc")
plot(perf1)
roc = performance(pred1,"tpr","fpr")
plot(roc, colorize = T, lwd = 2)
abline(a = 0, b = 1)
PARA MODELO 2:
pred = predict(Modelo_logistico2, newdata=LDA_logistic2, type="response")
pred = prediction(pred, LDA_logistic2$RESPONSE)
perf = performance(pred, "acc")
plot(perf)
roc = performance(pred,"tpr","fpr")
plot(roc, colorize = T, lwd = 2)
abline(a = 0, b = 1)