Este modulo de datos de muestra rastrea el abandono de clientes de una empresa de telecomunicaciones ficticia en funcion de varios factores.
La columna de abandono indica si el cliente se fue en el ultimo mes. Otras columnas incluyen genero, dependientes, cargos mensuales y muchas con informacion sobre los tipos de servicios que tiene cada cliente.
Se predice el abandono o no del cliente de la empresa.
A posteriori se realizara una evaluacion del negocio en base a la realizacion de una campaña de fidelizacion de clientes(clientes retenidos), en la cual se tomaran datos ficticios sobre la campaña y se mostraran varias opciones para determinar su mejor rendimiento posible y finalmente el ROI de la misma.
customerID: Recoge el ID de cada cliente
VARIABLES PERSONALES:
gender: Indica si el cliente es hombre o mujer
SeniorCitizen: Indica si el cliente tiene 65 años o mas: Si, No
Partner: Indica si el cliente tiene pareja o no
Dependents: Indica si el cliente vive con algunn dependiente: Si, No. Los dependientes pueden ser hijos, padres, abuelos, etc.
VARIABLE QUE INDICA ANTIGUEDAD:
tenure: Indica el nº de meses que el cliente ha permanecido en la cia, o que permanece
VARIABLES ASOCIADAS AL SERVICIO SUSCRITO:
PhoneService: Indica si el cliente tiene servicio telefonico o no
MultipleLines: Indica si el cliente tiene multiples lineas o no
InternetService: Indica el proveedor de servicios de internet del cliente
OnlineSecurity: Indica si el cliente tiene seguridad en linea o no
OnlineBackup: Indica si el cliente tiene respaldo en linea o no
DeviceProtection: Indica si el cliente tiene proteccion del dispositivo o no
TechSupport: Indica si el cliente tiene soporte tecnico o no
StreamingTV: Indica si el cliente tiene servicio de TV en streaming o no
StreamingMovies: Indica si el cliente dispone del servicio de peliculas en streaming o no.
VARIABLES ASOCIADAS AL TIPO DE CONTRATO:
Contract: Indica el plazo del contrato
PaperlessBilling: Indica si el cliente dispone de facturacion electronica o no
PaymentMethod: Indica el metodo de pago del cliente
VARIABLES ASOCIADAS AL IMPORTE Y TOTAL CARGADO:
MonthlyCharges:Indica el importe cobrado mensualmente al cliente
TotalCharges:Indica la cantidad total cargada al cliente
VARIABLE OBJETIVO:
Churn: Indica si el cliente abandona o no, sera la variable a predecir
#Carga de librerias:
library(doParallel)
library(data.table)
library(tidyverse)
library(ROCR)
library(rpart)
library(rpart.plot)
library(smbinning)
library(randomForest)
library(corrplot)
library(ggpubr)
library(cowplot)
library(tidymodels)
library(rattle)
library(RColorBrewer)
library(DataExplorer)
library(inspectdf)
library(kableExtra)
options(scipen=999)#Desactiva la notacion cientifica
#Carga de datos:
dt <- fread('Telco-Customer-Churn.csv')
as.data.frame(sort(names(dt))) # Visualizacion nombre de variables datos
## sort(names(dt))
## 1 Churn
## 2 Contract
## 3 customerID
## 4 Dependents
## 5 DeviceProtection
## 6 gender
## 7 InternetService
## 8 MonthlyCharges
## 9 MultipleLines
## 10 OnlineBackup
## 11 OnlineSecurity
## 12 PaperlessBilling
## 13 Partner
## 14 PaymentMethod
## 15 PhoneService
## 16 SeniorCitizen
## 17 StreamingMovies
## 18 StreamingTV
## 19 TechSupport
## 20 tenure
## 21 TotalCharges
str(dt) # Visualizacion estructura datos 1
## Classes 'data.table' and 'data.frame': 7043 obs. of 21 variables:
## $ customerID : chr "7590-VHVEG" "5575-GNVDE" "3668-QPYBK" "7795-CFOCW" ...
## $ gender : chr "Female" "Male" "Male" "Male" ...
## $ SeniorCitizen : int 0 0 0 0 0 0 0 0 0 0 ...
## $ Partner : chr "Yes" "No" "No" "No" ...
## $ Dependents : chr "No" "No" "No" "No" ...
## $ tenure : int 1 34 2 45 2 8 22 10 28 62 ...
## $ PhoneService : chr "No" "Yes" "Yes" "No" ...
## $ MultipleLines : chr "No phone service" "No" "No" "No phone service" ...
## $ InternetService : chr "DSL" "DSL" "DSL" "DSL" ...
## $ OnlineSecurity : chr "No" "Yes" "Yes" "Yes" ...
## $ OnlineBackup : chr "Yes" "No" "Yes" "No" ...
## $ DeviceProtection: chr "No" "Yes" "No" "Yes" ...
## $ TechSupport : chr "No" "No" "No" "Yes" ...
## $ StreamingTV : chr "No" "No" "No" "No" ...
## $ StreamingMovies : chr "No" "No" "No" "No" ...
## $ Contract : chr "Month-to-month" "One year" "Month-to-month" "One year" ...
## $ PaperlessBilling: chr "Yes" "No" "Yes" "No" ...
## $ PaymentMethod : chr "Electronic check" "Mailed check" "Mailed check" "Bank transfer (automatic)" ...
## $ MonthlyCharges : num 29.9 57 53.9 42.3 70.7 ...
## $ TotalCharges : num 29.9 1889.5 108.2 1840.8 151.7 ...
## $ Churn : chr "No" "No" "Yes" "No" ...
## - attr(*, ".internal.selfref")=<externalptr>
glimpse(dt) # Visualizacion estructura datos 2
## Rows: 7,043
## Columns: 21
## $ customerID <chr> "7590-VHVEG", "5575-GNVDE", "3668-QPYBK", "7795-CFOCW~
## $ gender <chr> "Female", "Male", "Male", "Male", "Female", "Female",~
## $ SeniorCitizen <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,~
## $ Partner <chr> "Yes", "No", "No", "No", "No", "No", "No", "No", "Yes~
## $ Dependents <chr> "No", "No", "No", "No", "No", "No", "Yes", "No", "No"~
## $ tenure <int> 1, 34, 2, 45, 2, 8, 22, 10, 28, 62, 13, 16, 58, 49, 2~
## $ PhoneService <chr> "No", "Yes", "Yes", "No", "Yes", "Yes", "Yes", "No", ~
## $ MultipleLines <chr> "No phone service", "No", "No", "No phone service", "~
## $ InternetService <chr> "DSL", "DSL", "DSL", "DSL", "Fiber optic", "Fiber opt~
## $ OnlineSecurity <chr> "No", "Yes", "Yes", "Yes", "No", "No", "No", "Yes", "~
## $ OnlineBackup <chr> "Yes", "No", "Yes", "No", "No", "No", "Yes", "No", "N~
## $ DeviceProtection <chr> "No", "Yes", "No", "Yes", "No", "Yes", "No", "No", "Y~
## $ TechSupport <chr> "No", "No", "No", "Yes", "No", "No", "No", "No", "Yes~
## $ StreamingTV <chr> "No", "No", "No", "No", "No", "Yes", "Yes", "No", "Ye~
## $ StreamingMovies <chr> "No", "No", "No", "No", "No", "Yes", "No", "No", "Yes~
## $ Contract <chr> "Month-to-month", "One year", "Month-to-month", "One ~
## $ PaperlessBilling <chr> "Yes", "No", "Yes", "No", "Yes", "Yes", "Yes", "No", ~
## $ PaymentMethod <chr> "Electronic check", "Mailed check", "Mailed check", "~
## $ MonthlyCharges <dbl> 29.85, 56.95, 53.85, 42.30, 70.70, 99.65, 89.10, 29.7~
## $ TotalCharges <dbl> 29.85, 1889.50, 108.15, 1840.75, 151.65, 820.50, 1949~
## $ Churn <chr> "No", "No", "Yes", "No", "Yes", "Yes", "No", "No", "Y~
# Visualicemos la estructura de datos y sus dimensiones de forma mas grafica.
plot_intro(dt)
# Inspecionamos visualmente variables categoricas(niveles):
dt %>%
inspect_cat() %>%
show_plot()
Reservamos variables independientes y variable objetivo susceptibles a pasar a factores.
dt_factor <- c('gender','SeniorCitizen', 'Partner', 'Dependents', 'PhoneService', 'MultipleLines', 'InternetService', 'OnlineSecurity','OnlineBackup', 'DeviceProtection','TechSupport', 'StreamingTV', 'StreamingMovies','Contract',
'PaperlessBilling','PaymentMethod','Churn')
lapply(dt,summary)
## $customerID
## Length Class Mode
## 7043 character character
##
## $gender
## Length Class Mode
## 7043 character character
##
## $SeniorCitizen
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 0.0000 0.0000 0.0000 0.1621 0.0000 1.0000
##
## $Partner
## Length Class Mode
## 7043 character character
##
## $Dependents
## Length Class Mode
## 7043 character character
##
## $tenure
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 0.00 9.00 29.00 32.37 55.00 72.00
##
## $PhoneService
## Length Class Mode
## 7043 character character
##
## $MultipleLines
## Length Class Mode
## 7043 character character
##
## $InternetService
## Length Class Mode
## 7043 character character
##
## $OnlineSecurity
## Length Class Mode
## 7043 character character
##
## $OnlineBackup
## Length Class Mode
## 7043 character character
##
## $DeviceProtection
## Length Class Mode
## 7043 character character
##
## $TechSupport
## Length Class Mode
## 7043 character character
##
## $StreamingTV
## Length Class Mode
## 7043 character character
##
## $StreamingMovies
## Length Class Mode
## 7043 character character
##
## $Contract
## Length Class Mode
## 7043 character character
##
## $PaperlessBilling
## Length Class Mode
## 7043 character character
##
## $PaymentMethod
## Length Class Mode
## 7043 character character
##
## $MonthlyCharges
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 18.25 35.50 70.35 64.76 89.85 118.75
##
## $TotalCharges
## Min. 1st Qu. Median Mean 3rd Qu. Max. NA's
## 18.8 401.4 1397.5 2283.3 3794.7 8684.8 11
##
## $Churn
## Length Class Mode
## 7043 character character
CONCLUSION: A priori podemos comprobar que la variable TotalCharges contiene 11 NA’s.
-Visualizacion grafica de posibles NA’s.
plot_missing(dt)
data.frame(colSums(is.na(dt)))
## colSums.is.na.dt..
## customerID 0
## gender 0
## SeniorCitizen 0
## Partner 0
## Dependents 0
## tenure 0
## PhoneService 0
## MultipleLines 0
## InternetService 0
## OnlineSecurity 0
## OnlineBackup 0
## DeviceProtection 0
## TechSupport 0
## StreamingTV 0
## StreamingMovies 0
## Contract 0
## PaperlessBilling 0
## PaymentMethod 0
## MonthlyCharges 0
## TotalCharges 11
## Churn 0
CONCLUSION: Corroboramos que hay 11 nulos en nuestros datos en la variable TotalCharges que posteriormente hay que corregir en nuestros datos.
Examen preliminar de nuestros datos con la funcion considerada, para comprobar si existe un problema de calidad en nuestros datos respecto al numero de ceros considerados en las variables.
contar_ceros <- function(variable) {
temp <- transmute(dt,if_else(variable==0,1,0))
sum(temp)
}
num_ceros <- sapply(dt,contar_ceros)
num_ceros <- data.frame(VARIABLE=names(num_ceros),CEROS = as.numeric(num_ceros),stringsAsFactors = F)
num_ceros <- num_ceros %>%
arrange(desc(CEROS)) %>%
mutate(PORCENTAJE = CEROS / nrow(dt) * 100)
num_ceros
## VARIABLE CEROS PORCENTAJE
## 1 SeniorCitizen 5901 83.7853188
## 2 tenure 11 0.1561834
## 3 customerID 0 0.0000000
## 4 gender 0 0.0000000
## 5 Partner 0 0.0000000
## 6 Dependents 0 0.0000000
## 7 PhoneService 0 0.0000000
## 8 MultipleLines 0 0.0000000
## 9 InternetService 0 0.0000000
## 10 OnlineSecurity 0 0.0000000
## 11 OnlineBackup 0 0.0000000
## 12 DeviceProtection 0 0.0000000
## 13 TechSupport 0 0.0000000
## 14 StreamingTV 0 0.0000000
## 15 StreamingMovies 0 0.0000000
## 16 Contract 0 0.0000000
## 17 PaperlessBilling 0 0.0000000
## 18 PaymentMethod 0 0.0000000
## 19 MonthlyCharges 0 0.0000000
## 20 Churn 0 0.0000000
## 21 TotalCharges NA NA
CONCLUSION: No apreciamos problemas de calidad de datos salvando los NA’s citados de TotalCharges que corregiremos con posterioridad.
out <- function(variable){
t(t(head(sort(variable,decreasing = T),20)))
}
lapply(dt,function(x){
if(is.double(x)) out(x)
})
## $customerID
## NULL
##
## $gender
## NULL
##
## $SeniorCitizen
## NULL
##
## $Partner
## NULL
##
## $Dependents
## NULL
##
## $tenure
## NULL
##
## $PhoneService
## NULL
##
## $MultipleLines
## NULL
##
## $InternetService
## NULL
##
## $OnlineSecurity
## NULL
##
## $OnlineBackup
## NULL
##
## $DeviceProtection
## NULL
##
## $TechSupport
## NULL
##
## $StreamingTV
## NULL
##
## $StreamingMovies
## NULL
##
## $Contract
## NULL
##
## $PaperlessBilling
## NULL
##
## $PaymentMethod
## NULL
##
## $MonthlyCharges
## [,1]
## [1,] 118.75
## [2,] 118.65
## [3,] 118.60
## [4,] 118.60
## [5,] 118.35
## [6,] 118.20
## [7,] 117.80
## [8,] 117.60
## [9,] 117.50
## [10,] 117.45
## [11,] 117.35
## [12,] 117.20
## [13,] 117.15
## [14,] 116.95
## [15,] 116.85
## [16,] 116.80
## [17,] 116.75
## [18,] 116.60
## [19,] 116.60
## [20,] 116.55
##
## $TotalCharges
## [,1]
## [1,] 8684.80
## [2,] 8672.45
## [3,] 8670.10
## [4,] 8594.40
## [5,] 8564.75
## [6,] 8547.15
## [7,] 8543.25
## [8,] 8529.50
## [9,] 8496.70
## [10,] 8477.70
## [11,] 8477.60
## [12,] 8476.50
## [13,] 8468.20
## [14,] 8456.75
## [15,] 8443.70
## [16,] 8436.25
## [17,] 8425.30
## [18,] 8425.15
## [19,] 8424.90
## [20,] 8405.00
##
## $Churn
## NULL
CONCLUSION: No apreciamos nada relevante en lo referente a los atipicos numericos.
out <- function(variable){
t(t(table(variable)))
}
lapply(dt,function(x){
if(is.integer(x)) out(x)
})
## $customerID
## NULL
##
## $gender
## NULL
##
## $SeniorCitizen
##
## variable [,1]
## 0 5901
## 1 1142
##
## $Partner
## NULL
##
## $Dependents
## NULL
##
## $tenure
##
## variable [,1]
## 0 11
## 1 613
## 2 238
## 3 200
## 4 176
## 5 133
## 6 110
## 7 131
## 8 123
## 9 119
## 10 116
## 11 99
## 12 117
## 13 109
## 14 76
## 15 99
## 16 80
## 17 87
## 18 97
## 19 73
## 20 71
## 21 63
## 22 90
## 23 85
## 24 94
## 25 79
## 26 79
## 27 72
## 28 57
## 29 72
## 30 72
## 31 65
## 32 69
## 33 64
## 34 65
## 35 88
## 36 50
## 37 65
## 38 59
## 39 56
## 40 64
## 41 70
## 42 65
## 43 65
## 44 51
## 45 61
## 46 74
## 47 68
## 48 64
## 49 66
## 50 68
## 51 68
## 52 80
## 53 70
## 54 68
## 55 64
## 56 80
## 57 65
## 58 67
## 59 60
## 60 76
## 61 76
## 62 70
## 63 72
## 64 80
## 65 76
## 66 89
## 67 98
## 68 100
## 69 95
## 70 119
## 71 170
## 72 362
##
## $PhoneService
## NULL
##
## $MultipleLines
## NULL
##
## $InternetService
## NULL
##
## $OnlineSecurity
## NULL
##
## $OnlineBackup
## NULL
##
## $DeviceProtection
## NULL
##
## $TechSupport
## NULL
##
## $StreamingTV
## NULL
##
## $StreamingMovies
## NULL
##
## $Contract
## NULL
##
## $PaperlessBilling
## NULL
##
## $PaymentMethod
## NULL
##
## $MonthlyCharges
## NULL
##
## $TotalCharges
## NULL
##
## $Churn
## NULL
a1 <- ggplot(dt, aes(y = tenure)) +
geom_boxplot(fill="blue")+ggtitle("TENURE")
a2 <- ggplot(dt, aes(y = MonthlyCharges)) +
geom_boxplot(fill="orange")+ggtitle("MONTTHLYCHARGES")
a3 <- ggplot(dt,aes(y = TotalCharges)) +
geom_boxplot(fill="green")+ggtitle("TOTALCHARGES")
plot_grid(a1, a2, a3,labels = "AUTO")
No se muestran visualmente en los graficos datos atipicos.
CONCLUSION: Nada a destacar fuera de lo normal.
CONCLUSION: Nada a destacar puesto que no hay mismas variables con datos a lo largo del tiempo para analizarlas de forma longitudinal.
Hemos detectado en el analisis exploratorio preliminar variables categoricas y niveles de las mismas poco coherentes con su contenido a pares o inclusive entre mas variables:
Por ejemplo se detecta en la variable ‘Multiplelines’ que existe un nivel que ya esta definido como variable categorica anteriormente como es PhoneService con el indicativo de ‘No Phone Service’ cuando no deberia esta asociado a la variable categorica ‘Multiplelines’ sino a ‘PhoneService’.
Por otro lado las variables categoricas’OnlineSecurity’,‘OnlineBackup’, ‘DeviceProtection’, ‘TechSupport, ’StreamingMovies’,StreamingTV, hacen referencia a productos que el cliente solo puede contratar, no tiene mucho sentido tener un nivel que diga ‘No internet Service’ en todas ellas.
CONCLUSION:Lo comentado produce multicolinealidad entre las variables categoricas respectivas, para evitarla habria que corregir el nivel de ‘Multiplelines’ y cambiar No Phone Service a un nivel ‘No’ de forma aditiva al nivel de ‘No’ ya existente,y exactamente lo mismo con ‘No internet Service’ para las variables categoricas ‘OnlineSecurity’, ‘OnlineBackup’,‘DeviceProtection’,‘TechSupport’,‘StreamingMovies’, ‘StreamingTV’.
Se considera pasar la variable factor SeniorCitizen de caracter(numerico) a string, para igualarla con las demas variables independientes.
Comprobaremos con posterioridad la correlacion entre variables continuas para ver si nos sobra alguna variable que ya este bien informada y realizaremos en el caso de haberla las correcciones oportunas.
En primer lugar vamos a ocuparnos de los NA’s suscitados en la variable ‘TotalCharges’, concretamente 11 NA’s.
En este caso podemos hacer 2 cosas: eliminar las observaciones que son 11 respecto al total de los datos que son 7043 observaciones y seria bastante licito.
o bien imputarlos como voy a hacer yo sustituyendo estos valores perdidos NA’s en este caso con la mediana, de tal forma que:
dt$TotalCharges[is.na(dt$TotalCharges)] <- median(dt$TotalCharges,na.rm = TRUE)
Si comprobamos nuevamente, vemos que ya no tenemos datos NA´s
plot_missing(dt)
Vamos a convertir a factores las variables independientes y la variable objetivo que hemos # apartado previamente.
Ademas de cambiar la variable factor ‘SeniorCitizen’ de caracter(numerico) a string. De tal manera que:
dt <- dt %>%
mutate_at(dt_factor,.funs = factor)
levels(dt$SeniorCitizen) <- c("No","Yes")
levels(dt$Churn) <- c("0","1")
# Comprobamos los cambios a factor:
glimpse(dt)
## Rows: 7,043
## Columns: 21
## $ customerID <chr> "7590-VHVEG", "5575-GNVDE", "3668-QPYBK", "7795-CFOCW~
## $ gender <fct> Female, Male, Male, Male, Female, Female, Male, Femal~
## $ SeniorCitizen <fct> No, No, No, No, No, No, No, No, No, No, No, No, No, N~
## $ Partner <fct> Yes, No, No, No, No, No, No, No, Yes, No, Yes, No, Ye~
## $ Dependents <fct> No, No, No, No, No, No, Yes, No, No, Yes, Yes, No, No~
## $ tenure <int> 1, 34, 2, 45, 2, 8, 22, 10, 28, 62, 13, 16, 58, 49, 2~
## $ PhoneService <fct> No, Yes, Yes, No, Yes, Yes, Yes, No, Yes, Yes, Yes, Y~
## $ MultipleLines <fct> No phone service, No, No, No phone service, No, Yes, ~
## $ InternetService <fct> DSL, DSL, DSL, DSL, Fiber optic, Fiber optic, Fiber o~
## $ OnlineSecurity <fct> No, Yes, Yes, Yes, No, No, No, Yes, No, Yes, Yes, No ~
## $ OnlineBackup <fct> Yes, No, Yes, No, No, No, Yes, No, No, Yes, No, No in~
## $ DeviceProtection <fct> No, Yes, No, Yes, No, Yes, No, No, Yes, No, No, No in~
## $ TechSupport <fct> No, No, No, Yes, No, No, No, No, Yes, No, No, No inte~
## $ StreamingTV <fct> No, No, No, No, No, Yes, Yes, No, Yes, No, No, No int~
## $ StreamingMovies <fct> No, No, No, No, No, Yes, No, No, Yes, No, No, No inte~
## $ Contract <fct> Month-to-month, One year, Month-to-month, One year, M~
## $ PaperlessBilling <fct> Yes, No, Yes, No, Yes, Yes, Yes, No, Yes, No, Yes, No~
## $ PaymentMethod <fct> Electronic check, Mailed check, Mailed check, Bank tr~
## $ MonthlyCharges <dbl> 29.85, 56.95, 53.85, 42.30, 70.70, 99.65, 89.10, 29.7~
## $ TotalCharges <dbl> 29.85, 1889.50, 108.15, 1840.75, 151.65, 820.50, 1949~
## $ Churn <fct> 0, 0, 1, 0, 1, 1, 0, 0, 1, 0, 0, 0, 0, 1, 0, 0, 0, 0,~
dt$MultipleLines <- recode(dt$MultipleLines, "No phone service"="No")
vect <- c('OnlineSecurity', 'OnlineBackup','DeviceProtection','TechSupport', 'StreamingMovies','StreamingTV')
for (i in vect){
dt[[i]] <- recode(dt[[i]], "No internet service"="No")}
# Comprobamos que ya esta corregidos los niveles en nuestros datos, y nos quedarian como aparecen:
dt %>%
inspect_cat() %>%
show_plot()
glimpse(dt)
## Rows: 7,043
## Columns: 21
## $ customerID <chr> "7590-VHVEG", "5575-GNVDE", "3668-QPYBK", "7795-CFOCW~
## $ gender <fct> Female, Male, Male, Male, Female, Female, Male, Femal~
## $ SeniorCitizen <fct> No, No, No, No, No, No, No, No, No, No, No, No, No, N~
## $ Partner <fct> Yes, No, No, No, No, No, No, No, Yes, No, Yes, No, Ye~
## $ Dependents <fct> No, No, No, No, No, No, Yes, No, No, Yes, Yes, No, No~
## $ tenure <int> 1, 34, 2, 45, 2, 8, 22, 10, 28, 62, 13, 16, 58, 49, 2~
## $ PhoneService <fct> No, Yes, Yes, No, Yes, Yes, Yes, No, Yes, Yes, Yes, Y~
## $ MultipleLines <fct> No, No, No, No, No, Yes, Yes, No, Yes, No, No, No, Ye~
## $ InternetService <fct> DSL, DSL, DSL, DSL, Fiber optic, Fiber optic, Fiber o~
## $ OnlineSecurity <fct> No, Yes, Yes, Yes, No, No, No, Yes, No, Yes, Yes, No,~
## $ OnlineBackup <fct> Yes, No, Yes, No, No, No, Yes, No, No, Yes, No, No, N~
## $ DeviceProtection <fct> No, Yes, No, Yes, No, Yes, No, No, Yes, No, No, No, Y~
## $ TechSupport <fct> No, No, No, Yes, No, No, No, No, Yes, No, No, No, No,~
## $ StreamingTV <fct> No, No, No, No, No, Yes, Yes, No, Yes, No, No, No, Ye~
## $ StreamingMovies <fct> No, No, No, No, No, Yes, No, No, Yes, No, No, No, Yes~
## $ Contract <fct> Month-to-month, One year, Month-to-month, One year, M~
## $ PaperlessBilling <fct> Yes, No, Yes, No, Yes, Yes, Yes, No, Yes, No, Yes, No~
## $ PaymentMethod <fct> Electronic check, Mailed check, Mailed check, Bank tr~
## $ MonthlyCharges <dbl> 29.85, 56.95, 53.85, 42.30, 70.70, 99.65, 89.10, 29.7~
## $ TotalCharges <dbl> 29.85, 1889.50, 108.15, 1840.75, 151.65, 820.50, 1949~
## $ Churn <fct> 0, 0, 1, 0, 1, 1, 0, 0, 1, 0, 0, 0, 0, 1, 0, 0, 0, 0,~
corr_cont <- dt%>%select(c("tenure", "MonthlyCharges", "TotalCharges"))
dat.cor <- cor(corr_cont, method = "pearson")
corrplot(dat.cor, method = "shade",
shade.col = NA, tl.col = "black",
tl.srt = 90,
addCoef.col = "black", addcolorlabel = "no",
order = "AOE")
CONCLUSION: Vemos sobre todo que hay una correlacion alta entre ’ ‘TotalCharges’ y las demas variables continuas ‘MonthlyCharges’ y ‘tenure’. El motivo es que la variable ‘TotalCharges’ parece ser combinacion lineal de las otras dos ‘Tenure’ y ‘MonthlyCharges’, ya que su producto da como resultado la variable ‘TotalCharges’. Vamos a tomar la accion de eliminar la variable ‘TotalCharges’ de nuestros datos porque consideramos que no la necesitamos por estar ya suficientemente informada.
dte <- dt # Copiamos el dataframe antes de hacer cambios, para posterior
# visualizacion grafica de variables al completo.
# Eliminamos variable
dt$TotalCharges=NULL
# Comprobamos los cambios:
glimpse(dt)
## Rows: 7,043
## Columns: 20
## $ customerID <chr> "7590-VHVEG", "5575-GNVDE", "3668-QPYBK", "7795-CFOCW~
## $ gender <fct> Female, Male, Male, Male, Female, Female, Male, Femal~
## $ SeniorCitizen <fct> No, No, No, No, No, No, No, No, No, No, No, No, No, N~
## $ Partner <fct> Yes, No, No, No, No, No, No, No, Yes, No, Yes, No, Ye~
## $ Dependents <fct> No, No, No, No, No, No, Yes, No, No, Yes, Yes, No, No~
## $ tenure <int> 1, 34, 2, 45, 2, 8, 22, 10, 28, 62, 13, 16, 58, 49, 2~
## $ PhoneService <fct> No, Yes, Yes, No, Yes, Yes, Yes, No, Yes, Yes, Yes, Y~
## $ MultipleLines <fct> No, No, No, No, No, Yes, Yes, No, Yes, No, No, No, Ye~
## $ InternetService <fct> DSL, DSL, DSL, DSL, Fiber optic, Fiber optic, Fiber o~
## $ OnlineSecurity <fct> No, Yes, Yes, Yes, No, No, No, Yes, No, Yes, Yes, No,~
## $ OnlineBackup <fct> Yes, No, Yes, No, No, No, Yes, No, No, Yes, No, No, N~
## $ DeviceProtection <fct> No, Yes, No, Yes, No, Yes, No, No, Yes, No, No, No, Y~
## $ TechSupport <fct> No, No, No, Yes, No, No, No, No, Yes, No, No, No, No,~
## $ StreamingTV <fct> No, No, No, No, No, Yes, Yes, No, Yes, No, No, No, Ye~
## $ StreamingMovies <fct> No, No, No, No, No, Yes, No, No, Yes, No, No, No, Yes~
## $ Contract <fct> Month-to-month, One year, Month-to-month, One year, M~
## $ PaperlessBilling <fct> Yes, No, Yes, No, Yes, Yes, Yes, No, Yes, No, Yes, No~
## $ PaymentMethod <fct> Electronic check, Mailed check, Mailed check, Bank tr~
## $ MonthlyCharges <dbl> 29.85, 56.95, 53.85, 42.30, 70.70, 99.65, 89.10, 29.7~
## $ Churn <fct> 0, 0, 1, 0, 1, 1, 0, 0, 1, 0, 0, 0, 0, 1, 0, 0, 0, 0,~
Guardamos sesion para poder recuperar el trabajo:
save.image(file='sesion1.RData')
#load(file='sesion1.RData')--para cargar
dtc <- dte%>%select(c(tenure,MonthlyCharges,TotalCharges))
s1 <-ggplot(data = dtc, aes(x=tenure))+
geom_histogram(color="white", fill="blue")+
theme_bw() +
theme(legend.position = "none")
s2 <- ggplot(data = dtc, aes(x=MonthlyCharges))+
geom_histogram(color="white", fill="orange") +
theme_bw() +
theme(legend.position = "none")
s3 <- ggplot(data = dtc, aes(x=TotalCharges))+
geom_histogram(color="white", fill="green") +
theme_bw() +
theme(legend.position = "none")
plot_grid(s1, s2, s3,labels = "AUTO")
plot_density(
dtc ,
ncol = 2,
title = "Distribucion variables continuas",
ggtheme = theme_bw(),
theme_config = list(
plot.title = element_text(size = 16, face = "bold"),
strip.text = element_text(colour = "black", size = 12, face = 2)
)
)
Podemos destacar en nuestros datos que en la variable ‘tenure’ que representa la permanencia en meses del cliente, que la mayor frecuencia de clientes la encontramos en los primeros meses de permanencia y en lo ultimos meses que son los que llevan mas tiempo (apareciendo como picos maximos en las graficas).
Tambien podemos destacar la larga cola de asimetria a la derecha de la variable ‘TotalCharges’ que indica la cantidad total cargada al cliente y que nos indica que hay pocos clientes con cantidades cargadas elevadas o muy elevadas en su totalidad, manteniendose la mayoria en un primer nivel de rango de cantidad cargada.
dts <- dt%>%select(-customerID)
plot_bar(
dts,
nrow = 5, ncol = 3,
title = "Numero de observaciones por grupo",
ggtheme = theme_bw(),
theme_config = list(plot.background = element_rect(fill = "yellow"),
plot.title = element_text(size = 16, face = "bold"),
strip.text = element_text(colour = "black", size = 10, face = 2),
legend.position = "none"
)
)
plot_bar(dts, by="Churn", nrow=5, ncol=3)
Un dato a destacar que podria ser interesante viendo las graficas de penetracion de la target, es que hay un abandono notable de los clientes que utilizan internet de fibra optica, con respecto a los que deciden no abandonar y respecto a otros servicios de internet como DSL. Seria un componente importante a retener tratando de fidelizar al cliente.
dt_ind <- dt %>%select(-Churn,-customerID)%>%names()
pre_rf <- randomForest(formula = reformulate(dt_ind,'Churn'), data= dt,mtry=2,ntree=50, importance = T)
imp_rf <- randomForest::importance(pre_rf)[,4]
imp_rf <- data.frame(VARIABLE = names(imp_rf), IMP_RF = imp_rf)
imp_rf <- imp_rf %>% arrange(desc(IMP_RF)) %>% mutate(RANKING_RF = 1:nrow(imp_rf))
temp <- mutate(dt,Churn = as.numeric(as.character(Churn))) %>% as.data.frame()
imp_iv <- smbinning.sumiv(temp[c(dt_ind,'Churn')],y="Churn")
##
##
|
| | 0%
|
|--- | 5%
|
|----- | 11%
|
|-------- | 16%
|
|----------- | 21%
|
|------------- | 26%
|
|---------------- | 32%
|
|------------------ | 37%
|
|--------------------- | 42%
|
|------------------------ | 47%
|
|-------------------------- | 53%
|
|----------------------------- | 58%
|
|-------------------------------- | 63%
|
|---------------------------------- | 68%
|
|------------------------------------- | 74%
|
|--------------------------------------- | 79%
|
|------------------------------------------ | 84%
|
|--------------------------------------------- | 89%
|
|----------------------------------------------- | 95%
|
|--------------------------------------------------| 100%
##
imp_iv <- imp_iv %>% mutate(Ranking = 1:nrow(imp_iv), IV = ifelse(is.na(.$IV),0,IV)) %>% select(-Process)
names(imp_iv) <- c('VARIABLE','IMP_IV','RANKING_IV')
imp_final <- inner_join(imp_rf,imp_iv,by='VARIABLE') %>%
select(VARIABLE,IMP_RF,IMP_IV,RANKING_RF,RANKING_IV) %>%
mutate(RANKING_TOT = RANKING_RF + RANKING_IV) %>%
arrange(RANKING_TOT)
imp_final
## VARIABLE IMP_RF IMP_IV RANKING_RF RANKING_IV RANKING_TOT
## 1 tenure 264.66089 0.8659 1 2 3
## 2 Contract 221.47715 1.2386 2 1 3
## 3 MonthlyCharges 206.25964 0.4842 3 4 7
## 4 InternetService 116.67450 0.6179 5 3 8
## 5 PaymentMethod 119.62190 0.4571 4 5 9
## 6 OnlineSecurity 48.17474 0.1719 6 7 13
## 7 PaperlessBilling 40.64527 0.2030 7 6 13
## 8 TechSupport 40.47674 0.1575 8 8 16
## 9 Dependents 29.04432 0.1555 10 9 19
## 10 SeniorCitizen 30.67409 0.1057 9 11 20
## 11 Partner 28.48917 0.1188 11 10 21
## 12 OnlineBackup 27.87895 0.0359 12 12 24
## 13 DeviceProtection 26.54301 0.0230 13 13 26
## 14 StreamingTV 25.29591 0.0202 15 14 29
## 15 MultipleLines 25.64500 0.0082 14 16 30
## 16 StreamingMovies 23.53669 0.0191 17 15 32
## 17 gender 24.78697 0.0004 16 18 34
## 18 PhoneService 13.63273 0.0008 18 17 35
Vamos a ver si los metodos son fiables para ello hacemos una correlacion entre ellos para ver si dan cosas similares.
cor(imp_final$IMP_RF,imp_final$IMP_IV)
## [1] 0.9145929
Vemos que da una correlacion alta por tanto lo consideramos confiable.
Vamos ahora a descartar las variables que no hayan salido entre las 11 primeras y la extraeremos añadiendo a mayores la variable objetivo y customerID:
dt1 <- dt %>% select(c('tenure','Contract', 'MonthlyCharges', 'InternetService', 'PaymentMethod' , 'PaperlessBilling', 'TechSupport', 'OnlineSecurity', 'Dependents', 'Partner', 'OnlineBackup', 'Churn', 'customerID'))
y guardamos nuestro dataframe.
saveRDS(dt1,'cache2.rds')
# dt1 <- readRDS('cache2.rds')--carga
discretizar <- function(vi,target){
temp_df <- data.frame(vi = vi, target = target)
temp_df$target <- as.numeric(as.character(temp_df$target))
disc <- smbinning(temp_df, y = 'target', x = 'vi')
return(disc)
}
# MONTHLYCHARGES:
disc_temp_MonthlyCharges <- discretizar(dt1$MonthlyCharges,dt1$Churn)
dt1_temp <- select(dt1,MonthlyCharges,Churn)
dt1_temp <- smbinning.gen(dt1_temp,disc_temp_MonthlyCharges,chrname = 'MONTHLYCHARGES_DISC')
dt1 <- cbind(dt1,dt1_temp[,3]) %>% select(-MonthlyCharges)
Escogemos para discretizacion manual la variable ‘tenure’ por su alto nivel de importancia y significaco de negocio.
Vemos la distribucion que tiene la variable:
OBJETIVOS:
-Conseguir una discretizacion de la variable lo mas similar a la distribucion original.
-Que la penetracion de la target con respecto a la variable sea monotonica(crece o decrece en todo su dominio)
# Densidad:
ggplot(dt1,aes(tenure))+geom_density()+scale_x_continuous(limits = c(0, 72))
#Histograma:
ggplot(dt1,aes(tenure))+geom_histogram(bins=25)+scale_x_continuous(limits = c(0, 72))
La variable decrece hasta un cierto punto, donde se mantiene de forma continua y al final vuelve a crecer.
Comportamiento de la variable con respecto a la target:
ggplot(dt1,aes(tenure,fill=Churn))+geom_histogram(bins=50, position='fill')+scale_x_continuous(limits = c(0, 72))
Aqui vemos al analizar en el grafico el comportamiento de la variable con respecto a la target tiene un comportamiento monotonico decreciente,vemos que cuando mas aumenta la permanencia(tenure), mas disminuye el abandono(churn).
Una vez hemos verificado ambas cosas vamos a apoyarnos en los deciles para selecionar los mejores cortes.
as.data.frame(quantile(dt1$tenure,prob = seq(0, 1, length = 11)))
## quantile(dt1$tenure, prob = seq(0, 1, length = 11))
## 0% 0
## 10% 2
## 20% 6
## 30% 12
## 40% 20
## 50% 29
## 60% 40
## 70% 50
## 80% 60
## 90% 69
## 100% 72
Hacemos los cortes en funcion de lo anteriormente visto:
dt1 <- dt1 %>% mutate(TENURE_DISC = as.factor(case_when(
between(tenure,-Inf,10) ~ '01_MENOR_10',
between(tenure,10,30) ~ '02_DE_10_A_30',
between(tenure,30,45) ~ '03_DE_30_A_45',
between(tenure,45,60) ~ '04_DE_45_A_60',
between(tenure,60,Inf) ~ '05_MAYOR_60',
TRUE ~ '00_ERROR'))
)
Eliminamos la variable ‘tenure’ original sobrante.
dt1$tenure <- NULL
Despues de probar diferentes cortes halle distribucion semejante en los cortes a la original.
ggplot(dt1,aes(TENURE_DISC)) + geom_bar()
Comprobamos si la penetracion de la target es monotonica.
ggplot(dt1,aes(TENURE_DISC, fill=Churn))+geom_bar(position='fill')
y efectivamente es monotona decreciente.
Ahora inspeccionamos visualmente todas las variables que son factores:
dt1 %>%
select_if(is.factor) %>%
gather() %>%
ggplot(aes(value)) +
geom_bar() +
facet_wrap(~ key, scales = "free") +
theme(axis.text=element_text(size=4))
Vemos ahora la penetracion de la target en cada categoria para ver las variables que han salido monotonicas.
a <- function(var1,var2) {
dt1_temp <- data.frame(var1 = dt1[[var1]],var2 = dt1[[var2]])
dt1_temp %>%
group_by(var1) %>%
summarise(Conteo = n(), Porc = mean(as.numeric(as.character(var2)))) %>%
ggplot(aes(var1,Porc)) + geom_bar(stat='identity') + xlab(var1)
}
dt1_nombres <- dt1 %>% select_if(is.factor) %>% names()
lapply(dt1_nombres,function(x){a(x,'Churn')})
## [[1]]
##
## [[2]]
##
## [[3]]
##
## [[4]]
##
## [[5]]
##
## [[6]]
##
## [[7]]
##
## [[8]]
##
## [[9]]
##
## [[10]]
##
## [[11]]
##
## [[12]]
Como podemos ver las unicas que no tienen una penetracion de la target monotonica son ‘PaymentMethod’, ‘MonthlyCharges’ y ‘InternetService’.
Lo dejamos aqui aunque se podria hacer el mismo proceso que hicimos con tenure para las otras variables en intentar ajustarlas para que sean monotonicas.
Ahora guardamos las discretizaciones en un objeto R,por que las necesitaremos despues si queremos poner el modelo en produccion.
discretizaciones <- list(
disc_temp_MonthlyCharges = disc_temp_MonthlyCharges
)
saveRDS(discretizaciones,'02_CortesDiscretizaciones.rds')
Reordenamos variables por simplemente por estetica y guardamos el dataframe con el que vamos a trabajar en el modelado:
centrales <- setdiff(names(dt1),c('customerID', 'Churn'))
dt1 <- dt1%>%select(
customerID,
one_of(centrales),
Churn)
Vemos como ha quedado el fichero antes de la modelizacion y guardamos.
glimpse(dt1)
## Rows: 7,043
## Columns: 13
## $ customerID <chr> "7590-VHVEG", "5575-GNVDE", "3668-QPYBK", "7795-CF~
## $ Contract <fct> Month-to-month, One year, Month-to-month, One year~
## $ InternetService <fct> DSL, DSL, DSL, DSL, Fiber optic, Fiber optic, Fibe~
## $ PaymentMethod <fct> Electronic check, Mailed check, Mailed check, Bank~
## $ PaperlessBilling <fct> Yes, No, Yes, No, Yes, Yes, Yes, No, Yes, No, Yes,~
## $ TechSupport <fct> No, No, No, Yes, No, No, No, No, Yes, No, No, No, ~
## $ OnlineSecurity <fct> No, Yes, Yes, Yes, No, No, No, Yes, No, Yes, Yes, ~
## $ Dependents <fct> No, No, No, No, No, No, Yes, No, No, Yes, Yes, No,~
## $ Partner <fct> Yes, No, No, No, No, No, No, No, Yes, No, Yes, No,~
## $ OnlineBackup <fct> Yes, No, Yes, No, No, No, Yes, No, No, Yes, No, No~
## $ MONTHLYCHARGES_DISC <fct> 02 <= 55.95, 03 <= 68.8, 02 <= 55.95, 02 <= 55.95,~
## $ TENURE_DISC <fct> 01_MENOR_10, 03_DE_30_A_45, 01_MENOR_10, 03_DE_30_~
## $ Churn <fct> 0, 0, 1, 0, 1, 1, 0, 0, 1, 0, 0, 0, 0, 1, 0, 0, 0,~
saveRDS(dt1,'cache3.rds')
#dt1 <- readRDS('cache3.rds')
Funcion para crear una matriz de confusion
confusion<-function(real,scoring,umbral){
conf<-table(real,scoring>=umbral)
if(ncol(conf)==2) return(conf) else return(NULL)
}
Funcion para calcular las metricas de los modelos: acierto, precision, cobertura y F1
metricas<-function(matriz_conf){
acierto <- (matriz_conf[1,1] + matriz_conf[2,2]) / sum(matriz_conf) *100
precision <- matriz_conf[2,2] / (matriz_conf[2,2] + matriz_conf[1,2]) *100
cobertura <- matriz_conf[2,2] / (matriz_conf[2,2] + matriz_conf[2,1]) *100
F1 <- 2*precision*cobertura/(precision+cobertura)
salida<-c(acierto,precision,cobertura,F1)
return(salida)
}
Funcion para probar distintos umbrales y ver el efecto sobre precision y cobertura.
umbrales<-function(real,scoring){
umbrales<-data.frame(umbral=rep(0,times=19),acierto=rep(0,times=19),precision=rep(0,times=19),cobertura=rep(0,times=19),F1=rep(0,times=19))
cont <- 1
for (cada in seq(0.05,0.95,by = 0.05)){
datos<-metricas(confusion(real,scoring,cada))
registro<-c(cada,datos)
umbrales[cont,]<-registro
cont <- cont + 1
}
return(umbrales)
}
Funciones que calculan la curva ROC y el AUC.
roc<-function(prediction){
r<-performance(prediction,'tpr','fpr')
plot(r)
}
auc<-function(prediction){
a<-performance(prediction,'auc')
return(a@y.values[[1]])
}
En primer lugar vamos a ver los datos que nos arroja en proporcion nuestra variable objetivo, por si baja de un nivel minimo aceptable alguno de sus niveles o factores y hubiera que balancear los datos (entendiendo cuestionarselo por debajo del 10% de algun nivel).
barplot(prop.table(table(dt1$Churn)),col=c("orange","blue"))
y comprobamos que no se da el caso, lo consideramos suficiente y seguimos con nuestro proceso.
Vamos a utilizar para ello el paquete tidymodels para ayudarnos en esta tarea.
Vamos a eliminar de nuestros datos ‘CustomerID’, puesto que es un identificador que no tiene ningun poder predictivo ni descriptivo.
dt3 <- dt1%>%select(-customerID)
y vamos a dividir nuestros datos en un 70 % para entrenamiento y un 30% para test. ademas para todos los modelos consideraremos un cross-validation con vfolds=5.
set.seed(1971)
split_inicial <- initial_split(
data = dt3,
prop = 0.7,
strata = Churn
)
train <- training(split_inicial)
test <- testing(split_inicial)
mod_rl <- logistic_reg(mode = 'classification') %>%
set_engine(engine = "glm",family=binomial(link='logit'))
mod_rcp <- recipe(Churn~.,
data = train)
set.seed(1981)
cv_folds <- vfold_cv(
data = train,
v = 5,
strata = Churn
)
wf_m <- workflow()%>%
add_recipe(mod_rcp)%>%
add_model(mod_rl)
my_metrics <- metric_set(accuracy,roc_auc,kap)
set.seed(1981)
rlog <- wf_m%>%
fit_resamples(
resamples = cv_folds,
metrics = my_metrics,
control = control_resamples(save_pred = TRUE)
)
collect_metrics(rlog)
## # A tibble: 3 x 6
## .metric .estimator mean n std_err .config
## <chr> <chr> <dbl> <int> <dbl> <chr>
## 1 accuracy binary 0.799 5 0.00270 Preprocessor1_Model1
## 2 kap binary 0.445 5 0.00972 Preprocessor1_Model1
## 3 roc_auc binary 0.843 5 0.00417 Preprocessor1_Model1
show_best(rlog, metric = "roc_auc")
## # A tibble: 1 x 6
## .metric .estimator mean n std_err .config
## <chr> <chr> <dbl> <int> <dbl> <chr>
## 1 roc_auc binary 0.843 5 0.00417 Preprocessor1_Model1
best_rlog1 <- rlog%>%
select_best("roc_auc")
wf_final <-
wf_m%>%
finalize_workflow(best_rlog1)
set.seed(2100)
final_model_rl <- fit(object = wf_final, data = train)
model_extr<- final_model_rl%>%extract_fit_engine()
summary(model_extr)
##
## Call:
## stats::glm(formula = ..y ~ ., family = ~binomial(link = "logit"),
## data = data)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -1.8189 -0.6882 -0.2902 0.6885 3.1644
##
## Coefficients:
## Estimate Std. Error z value
## (Intercept) -0.29338 0.30450 -0.964
## ContractOne year -0.66256 0.12655 -5.236
## ContractTwo year -1.80472 0.22776 -7.924
## InternetServiceFiber optic 0.75268 0.20826 3.614
## InternetServiceNo -1.11629 0.31148 -3.584
## PaymentMethodCredit card (automatic) -0.26035 0.13568 -1.919
## PaymentMethodElectronic check 0.36789 0.11147 3.300
## PaymentMethodMailed check -0.13329 0.13697 -0.973
## PaperlessBillingYes 0.42361 0.08852 4.785
## TechSupportYes -0.34383 0.10520 -3.268
## OnlineSecurityYes -0.47564 0.10211 -4.658
## DependentsYes -0.16178 0.10590 -1.528
## PartnerYes 0.06647 0.09184 0.724
## OnlineBackupYes -0.10219 0.09253 -1.104
## MONTHLYCHARGES_DISC02 <= 55.95 -0.02616 0.29964 -0.087
## MONTHLYCHARGES_DISC03 <= 68.8 -0.47830 0.34362 -1.392
## MONTHLYCHARGES_DISC04 <= 106.75 0.06797 0.35656 0.191
## MONTHLYCHARGES_DISC05 > 106.75 0.63719 0.43667 1.459
## TENURE_DISC02_DE_10_A_30 -1.02598 0.10322 -9.939
## TENURE_DISC03_DE_30_A_45 -1.05110 0.13295 -7.906
## TENURE_DISC04_DE_45_A_60 -1.27203 0.15178 -8.381
## TENURE_DISC05_MAYOR_60 -1.69245 0.19899 -8.505
## Pr(>|z|)
## (Intercept) 0.335296
## ContractOne year 0.00000016453097794 ***
## ContractTwo year 0.00000000000000231 ***
## InternetServiceFiber optic 0.000301 ***
## InternetServiceNo 0.000339 ***
## PaymentMethodCredit card (automatic) 0.055009 .
## PaymentMethodElectronic check 0.000965 ***
## PaymentMethodMailed check 0.330488
## PaperlessBillingYes 0.00000170711003889 ***
## TechSupportYes 0.001081 **
## OnlineSecurityYes 0.00000318916428306 ***
## DependentsYes 0.126595
## PartnerYes 0.469234
## OnlineBackupYes 0.269420
## MONTHLYCHARGES_DISC02 <= 55.95 0.930439
## MONTHLYCHARGES_DISC03 <= 68.8 0.163947
## MONTHLYCHARGES_DISC04 <= 106.75 0.848821
## MONTHLYCHARGES_DISC05 > 106.75 0.144508
## TENURE_DISC02_DE_10_A_30 < 0.0000000000000002 ***
## TENURE_DISC03_DE_30_A_45 0.00000000000000266 ***
## TENURE_DISC04_DE_45_A_60 < 0.0000000000000002 ***
## TENURE_DISC05_MAYOR_60 < 0.0000000000000002 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 5703.8 on 4928 degrees of freedom
## Residual deviance: 4101.7 on 4907 degrees of freedom
## AIC: 4145.7
##
## Number of Fisher Scoring iterations: 6
Podemos observar a la hora de interpretar los coeficientes de estimacion y signos de los mismos que en general todos se encuentran en un mismo rango de valores sin coeficientes muy desproporcionados que se salgan de ese mismo rango.
La clase predictora son los que abandonan como clientes, por tanto podemos decir que segun los resultados los que presentan menor tasa de abandono son los clientes que tienen contrato de 2 años y los clientes con una permanencia mayor de 60 meses. Entendemos que esto tiene sentido debido a la fidelidad del cliente con mas tiempo de permanencia y mayor contrato en la empresa.
Asi mismo se interpreta que los que mayor tasa de abandono presentan son los que tienen mayores cargos mensuales, esto, tambien tiene sentido, que los que soportan mayores gastos consideren abandonar por considerar caro el servicio u considerar ofertas mejores de la competencia.
y tambien los clientes con servicio de internet de fibra optica, por algun motivo de negocio que desconocemos.
Despues de analizar los p-valores vamos a escoger las siguientes variables basandonos en el criterio de significacion del 99,9 % de intervalo de confianza(***) considerando las mas predictivas de todas y eliminar las restantes, de tal forma nos quedarian:
‘Contract’,‘InternetService’,‘PaymentMethod’,‘PaperlessBilling’,‘OnlineSecurity’,‘TENURE_DISC’
Ahora aplicamos las variables escogidas en un nuevo modelo de regresion logistica.
dt4 <- dt3%>%select(-TechSupport,-Dependents,-Partner,-OnlineBackup,-MONTHLYCHARGES_DISC)
set.seed(1981)
split_inicial <- initial_split(
data = dt4,
prop = 0.7,
strata = Churn
)
train <- training(split_inicial)
test <- testing(split_inicial)
mod_rl <- logistic_reg(mode = 'classification') %>%
set_engine(engine = "glm",family=binomial(link='logit'))
mod_rcp <- recipe(Churn~.,
data = train)
set.seed(1981)
cv_folds <- vfold_cv(
data = train,
v = 5,
strata = Churn
)
wf_m <- workflow()%>%
add_recipe(mod_rcp)%>%
add_model(mod_rl)
my_metrics <- metric_set(accuracy,roc_auc,kap)
set.seed(1981)
rlog <- wf_m%>%
fit_resamples(
resamples = cv_folds,
metrics = my_metrics,
control = control_resamples(save_pred = TRUE)
)
collect_metrics(rlog)
## # A tibble: 3 x 6
## .metric .estimator mean n std_err .config
## <chr> <chr> <dbl> <int> <dbl> <chr>
## 1 accuracy binary 0.799 5 0.00463 Preprocessor1_Model1
## 2 kap binary 0.433 5 0.0158 Preprocessor1_Model1
## 3 roc_auc binary 0.839 5 0.00532 Preprocessor1_Model1
show_best(rlog, metric = "roc_auc")
## # A tibble: 1 x 6
## .metric .estimator mean n std_err .config
## <chr> <chr> <dbl> <int> <dbl> <chr>
## 1 roc_auc binary 0.839 5 0.00532 Preprocessor1_Model1
best_rlog2 <- rlog%>%
select_best("roc_auc")
wf_final <-
wf_m%>%
finalize_workflow(best_rlog2)
set.seed(2100)
final_model_rl <- fit(object = wf_final, data = train)
model_extr1 <- final_model_rl%>%extract_fit_engine()
summary(model_extr1)
##
## Call:
## stats::glm(formula = ..y ~ ., family = ~binomial(link = "logit"),
## data = data)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -1.7274 -0.6874 -0.2986 0.7139 3.0315
##
## Coefficients:
## Estimate Std. Error z value
## (Intercept) -0.51965 0.14146 -3.673
## ContractOne year -0.80412 0.12543 -6.411
## ContractTwo year -1.67381 0.20787 -8.052
## InternetServiceFiber optic 0.94559 0.09213 10.263
## InternetServiceNo -0.98636 0.14874 -6.632
## PaymentMethodCredit card (automatic) -0.06093 0.13529 -0.450
## PaymentMethodElectronic check 0.43820 0.11364 3.856
## PaymentMethodMailed check -0.08777 0.13830 -0.635
## PaperlessBillingYes 0.37299 0.08783 4.247
## OnlineSecurityYes -0.40244 0.09967 -4.038
## TENURE_DISC02_DE_10_A_30 -0.98493 0.09959 -9.890
## TENURE_DISC03_DE_30_A_45 -1.19957 0.12972 -9.247
## TENURE_DISC04_DE_45_A_60 -1.31725 0.14502 -9.083
## TENURE_DISC05_MAYOR_60 -1.76452 0.18698 -9.437
## Pr(>|z|)
## (Intercept) 0.000239 ***
## ContractOne year 0.000000000144884445 ***
## ContractTwo year 0.000000000000000813 ***
## InternetServiceFiber optic < 0.0000000000000002 ***
## InternetServiceNo 0.000000000033227057 ***
## PaymentMethodCredit card (automatic) 0.652437
## PaymentMethodElectronic check 0.000115 ***
## PaymentMethodMailed check 0.525678
## PaperlessBillingYes 0.000021669630671992 ***
## OnlineSecurityYes 0.000054004096979827 ***
## TENURE_DISC02_DE_10_A_30 < 0.0000000000000002 ***
## TENURE_DISC03_DE_30_A_45 < 0.0000000000000002 ***
## TENURE_DISC04_DE_45_A_60 < 0.0000000000000002 ***
## TENURE_DISC05_MAYOR_60 < 0.0000000000000002 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 5703.8 on 4928 degrees of freedom
## Residual deviance: 4157.8 on 4915 degrees of freedom
## AIC: 4185.8
##
## Number of Fisher Scoring iterations: 6
y vemos que ahora las variables de la relacion tienen al menos una categoria con 3 estrellas(***).
La Prediccion del scoring del modelo (conjunto test) sera:
rl_predict <- predict(model_extr1, test, type = "response")
head(rl_predict,20)
## 2 7 10 12 15 17 18
## 0.04681705 0.43857893 0.02957676 0.01440497 0.46250901 0.02377296 0.03183484
## 19 24 32 33 34 35 36
## 0.35879783 0.03921674 0.67655702 0.05738301 0.18153105 0.37293345 0.03183484
## 41 42 49 55 57 58
## 0.19599097 0.01715155 0.02656879 0.06114540 0.09681750 0.20890767
El resultado es el scoring de la probabilidad de un cliente de que abandone la empresa.
Veamos graficamente la pinta que tiene:
plot(rl_predict~test$Churn)
Calculamos ahora el pseudo R cuadrado, visualizando con anterioridad los coeficientes y significancia de acuerdo con lo que hemos decidido escoger para nuestro modelo:
psdo_R <- 1-(model_extr1$deviance/model_extr1$null.deviance)
psdo_R
## [1] 0.2710531
El valor que nos da no es muy elevado, este valor nos indica que con las variables que tenemos somos capaces de explicar el 27% de los motivos por los que un cliente abandona, hay mas de un 70% que no somos capaces de explicar, pero con la información disponible es lo mejor hemos podido obtener.
Ahora tenemos que transformar la probabilidad en una decision de si el cliente va a abandonar o se mantiene.
Con la funcion umbrales probamos diferentes cortes y para este caso vamos a determinar el punto de corte que maximiza a la F1(Punto intermedio entre precision y cobertura)
umb_rl<-umbrales(test$Churn,rl_predict)
umb_rl
## umbral acierto precision cobertura F1
## 1 0.05 52.69631 35.79288 98.57398 52.51662
## 2 0.10 62.25166 40.91954 95.18717 57.23473
## 3 0.15 66.55629 43.77133 91.44385 59.20369
## 4 0.20 71.14475 47.58621 86.09626 61.29442
## 5 0.25 73.93567 50.56306 80.03565 61.97378
## 6 0.30 76.96310 54.75578 75.93583 63.62957
## 7 0.35 77.38884 55.72414 72.01426 62.83048
## 8 0.40 78.85525 59.28339 64.88414 61.95745
## 9 0.45 78.76064 60.64639 56.86275 58.69365
## 10 0.50 79.09177 64.76427 46.52406 54.14938
## 11 0.55 79.23368 66.66667 43.49376 52.64293
## 12 0.60 77.90918 70.79646 28.52050 40.66074
## 13 0.65 77.57805 70.61611 26.55971 38.60104
## 14 0.70 77.01041 74.50980 20.32086 31.93277
## 15 0.75 76.44276 75.20000 16.75579 27.40525
## 16 0.80 0.80000 0.80000 0.80000 0.80000
## 17 0.85 0.85000 0.85000 0.85000 0.85000
## 18 0.90 0.90000 0.90000 0.90000 0.90000
## 19 0.95 0.95000 0.95000 0.95000 0.95000
Seleccionamos la formula para hallarlo, aunque tambien podemos deducirlo de la tabla que seria el numero mas alto de la tabla F1 del cual el siguiente baja, el que de el umbral corresponderia a esa fila.
umbral_final_rl<-umb_rl[which.max(umb_rl$F1),1]
umbral_final_rl
## [1] 0.3
Hallamos la matriz de confusion del umbral optimizado.
confusion(test$Churn,rl_predict,0.3)
##
## real FALSE TRUE
## 0 1201 352
## 1 135 426
rl_metricas<-filter(umb_rl,umbral==umbral_final_rl)
rl_metricas
## umbral acierto precision cobertura F1
## 1 0.3 76.9631 54.75578 75.93583 63.62957
Evaluamos la ROC
#creamos el objeto prediction
rl_prediction<-prediction(rl_predict,test$Churn)
#visualizamos la ROC
roc(rl_prediction)
#Sacamos las metricas definitivas incluyendo el AUC
rl_metricas<-cbind(rl_metricas,AUC=round(auc(rl_prediction),2)*100)
print(t(rl_metricas))
## [,1]
## umbral 0.30000
## acierto 76.96310
## precision 54.75578
## cobertura 75.93583
## F1 63.62957
## AUC 84.00000
set.seed(3001)
split_inicial <- initial_split(
data = dt3,
prop = 0.7,
strata = Churn
)
train <- training(split_inicial)
test <- testing(split_inicial)
mod_rp <- decision_tree(
mode = 'classification',
cost_complexity=0.00001,
tree_depth = tune(),
min_n = tune()
)%>%
set_engine(engine = "rpart",parms = list(split="information"))
mod_rcp <- recipe(Churn~.,
data = train)
set.seed(3001)
cv_folds <- vfold_cv(
data = train,
v = 5,
strata = Churn
)
wf_rp <- workflow()%>%
add_recipe(mod_rcp)%>%
add_model(mod_rp)
registerDoParallel(cores = parallel::detectCores() - 1)
set.seed(2341)
rp_fit <- wf_rp%>%tune_bayes(
resamples = cv_folds,
initial = 15,
iter = 20,
metrics = metric_set(roc_auc),
control = control_bayes(save_pred=TRUE,no_improve = 10,
verbose = FALSE))
stopImplicitCluster()
show_best(rp_fit, metric = "roc_auc")
## # A tibble: 5 x 9
## tree_depth min_n .metric .estimator mean n std_err .config .iter
## <int> <int> <chr> <chr> <dbl> <int> <dbl> <chr> <int>
## 1 8 28 roc_auc binary 0.809 5 0.00708 Iter7 7
## 2 7 11 roc_auc binary 0.808 5 0.00448 Preprocessor1_M~ 0
## 3 8 38 roc_auc binary 0.808 5 0.00732 Iter6 6
## 4 10 24 roc_auc binary 0.808 5 0.00403 Preprocessor1_M~ 0
## 5 13 23 roc_auc binary 0.808 5 0.00362 Iter17 17
mejor_hiper1 <- select_best(rp_fit, metric="roc_auc")
mejor_hiper1
## # A tibble: 1 x 3
## tree_depth min_n .config
## <int> <int> <chr>
## 1 8 28 Iter7
wf_finalrp <-
wf_rp%>%
finalize_workflow(mejor_hiper1)
set.seed(3100)
final_model_rp <- fit(object = wf_finalrp, data = train)
model_ext <- final_model_rp%>%
extract_fit_engine()
Visualizamos el arbol del modelo.
fancyRpartPlot(model_ext)
Ahora revisamos la tabla donde podemos ver en primer termino el parametro de complejidad, error relativo,error absoluto, y desviacion estandar.
printcp(model_ext)
##
## Classification tree:
## rpart::rpart(formula = ..y ~ ., data = data, parms = ~list(split = "information"),
## cp = ~0.00001, maxdepth = ~8, minsplit = min_rows(28, data))
##
## Variables actually used in tree construction:
## [1] Contract Dependents InternetService
## [4] MONTHLYCHARGES_DISC OnlineBackup OnlineSecurity
## [7] PaperlessBilling Partner PaymentMethod
## [10] TechSupport TENURE_DISC
##
## Root node error: 1308/4929 = 0.26537
##
## n= 4929
##
## CP nsplit rel error xerror xstd
## 1 0.06090724 0 1.00000 1.00000 0.023699
## 2 0.01490826 3 0.81728 0.81728 0.022121
## 3 0.00324924 5 0.78746 0.79281 0.021877
## 4 0.00267584 10 0.76988 0.82798 0.022224
## 5 0.00229358 12 0.76453 0.81498 0.022098
## 6 0.00152905 15 0.75765 0.81040 0.022053
## 7 0.00114679 17 0.75459 0.81728 0.022121
## 8 0.00076453 19 0.75229 0.81957 0.022143
## 9 0.00050968 22 0.75000 0.82187 0.022165
## 10 0.00045872 25 0.74847 0.82569 0.022202
## 11 0.00025484 36 0.74312 0.83028 0.022247
## 12 0.00001000 39 0.74235 0.83410 0.022283
En ella nos fijamos en el parametro de complejidad desde el que tomamos de partida 0.00001 abajo del todo hasta arriba de esa relacion, tenemos que mirar en la columna xerror(error absoluto), que empieza a bajar desde arriba, En el momento que encontramos un numero donde empieza a cambiar la tendencia (a aumentar) ese seria, el parametro de complejidad que tenemos que tomar correspondiente a ese valor, y por ahi vamos a podar el arbol y ejecutar otra vez el modelo.
En este caso el parametro de complejidad que corresponderia a lo dicho anteriormente parece que minimiza aproximadamente en 0.001.
Tambien lo podemos interpretar ayudandonos viendo la representacion grafica como vemos a continuacion.
plotcp(model_ext)
Cambiamos el parametro de complejidad de partida a 0.001 e incluimos en el arbol que no tenga mas de 5 niveles.
set.seed(3001)
split_inicial <- initial_split(
data = dt3,
prop = 0.7,
strata = Churn)
mod_rp <- decision_tree(
mode = 'classification',
cost_complexity=0.001,
tree_depth =5,
min_n = tune()
)%>%
set_engine(engine = "rpart",parms = list(split="information"))
mod_rcp <- recipe(Churn~.,
data = train)
set.seed(3001)
cv_folds <- vfold_cv(
data = train,
v = 5,
strata = Churn)
wf_rp <- workflow()%>%
add_recipe(mod_rcp)%>%
add_model(mod_rp)
registerDoParallel(cores = parallel::detectCores() - 1)
set.seed(2201)
rp_fit <- wf_rp%>%tune_bayes(
resamples = cv_folds,
initial = 15,
iter = 30,
metrics = metric_set(roc_auc),
control = control_bayes(save_pred=TRUE,no_improve = 10, verbose = FALSE))
stopImplicitCluster()
show_best(rp_fit, metric = "roc_auc")
## # A tibble: 5 x 8
## min_n .metric .estimator mean n std_err .config .iter
## <int> <chr> <chr> <dbl> <int> <dbl> <chr> <int>
## 1 31 roc_auc binary 0.793 5 0.00922 Preprocessor1_Model10 0
## 2 29 roc_auc binary 0.793 5 0.00922 Preprocessor1_Model13 0
## 3 30 roc_auc binary 0.793 5 0.00922 Iter1 1
## 4 27 roc_auc binary 0.793 5 0.00922 Iter2 2
## 5 28 roc_auc binary 0.793 5 0.00922 Iter3 3
mejor_hiper2 <- select_best(rp_fit, metric="roc_auc")
mejor_hiper2
## # A tibble: 1 x 2
## min_n .config
## <int> <chr>
## 1 31 Preprocessor1_Model10
wf_finalrp <-
wf_rp%>%
finalize_workflow(mejor_hiper2)
set.seed(3100)
final_model_rp <- fit(object = wf_finalrp, data = train)
model_extr2 <- final_model_rp%>%
extract_fit_engine()
Visualizamos tabla complejidad-errores.
printcp(model_extr2)
##
## Classification tree:
## rpart::rpart(formula = ..y ~ ., data = data, parms = ~list(split = "information"),
## cp = ~0.001, maxdepth = ~5, minsplit = min_rows(31, data))
##
## Variables actually used in tree construction:
## [1] Contract InternetService OnlineSecurity PaymentMethod
## [5] TENURE_DISC
##
## Root node error: 1308/4929 = 0.26537
##
## n= 4929
##
## CP nsplit rel error xerror xstd
## 1 0.0609072 0 1.00000 1.00000 0.023699
## 2 0.0149083 3 0.81728 0.81728 0.022121
## 3 0.0022936 5 0.78746 0.80352 0.021985
## 4 0.0010000 7 0.78287 0.80581 0.022007
Podemos ver una clara tendencia descendente de xerror(error absoluto) hasta llegar al parametro de complejidad 0.0007 de partida y tanto su error relativo como el absoluto son bastane proximos en dicho parametro y en los demas.
Por tanto consideramos como buena dicha aproximacion.
Visualizamos grafica error validacion cruzada y complejidad.
plotcp(model_extr2)
Visualizamos el arbol del modelo.
fancyRpartPlot(model_extr2)
Reglas que podrian ser utiizadas para hacer implantacion del arbol.(y es Churn–la variable objetivo(dependiente) en las reglas del arbol)
rpart.rules(model_extr2, style = 'tall', cover = T,roundint=FALSE)
## ..y is 0.07 with cover 45% when
## Contract is One year or Two year
##
## ..y is 0.25 with cover 0% when
## Contract is Month-to-month
## InternetService is Fiber optic
## TENURE_DISC is 01_MENOR_10
## PaymentMethod is Mailed check
## OnlineSecurity is Yes
##
## ..y is 0.29 with cover 25% when
## Contract is Month-to-month
## InternetService is DSL or No
##
## ..y is 0.34 with cover 8% when
## Contract is Month-to-month
## InternetService is Fiber optic
## TENURE_DISC is 02_DE_10_A_30 or 03_DE_30_A_45 or 04_DE_45_A_60 or 05_MAYOR_60
## PaymentMethod is Bank transfer (automatic) or Credit card (automatic) or Mailed check
##
## ..y is 0.42 with cover 5% when
## Contract is Month-to-month
## InternetService is Fiber optic
## TENURE_DISC is 03_DE_30_A_45 or 04_DE_45_A_60 or 05_MAYOR_60
## PaymentMethod is Electronic check
##
## ..y is 0.54 with cover 1% when
## Contract is Month-to-month
## InternetService is Fiber optic
## TENURE_DISC is 01_MENOR_10
## PaymentMethod is Bank transfer (automatic) or Credit card (automatic) or Electronic check
## OnlineSecurity is Yes
##
## ..y is 0.58 with cover 5% when
## Contract is Month-to-month
## InternetService is Fiber optic
## TENURE_DISC is 02_DE_10_A_30
## PaymentMethod is Electronic check
##
## ..y is 0.73 with cover 11% when
## Contract is Month-to-month
## InternetService is Fiber optic
## TENURE_DISC is 01_MENOR_10
## OnlineSecurity is No
Vemos el numero de nodos final por cliente.
rp_numnod <- rpart.predict(model_extr2,test,nn=T)
head(rp_numnod)
## 0 1 nn
## 4 0.9279238 0.07207616 2
## 9 0.4235294 0.57647059 59
## 12 0.9279238 0.07207616 2
## 13 0.9279238 0.07207616 2
## 15 0.4235294 0.57647059 59
## 18 0.9279238 0.07207616 2
Calculo de prediccion del modelo (scorings).
rp_predict <- predict(model_extr2, test, type = 'prob')[,2]
head(rp_predict,20)
## 4 9 12 13 15 18 21
## 0.07207616 0.57647059 0.07207616 0.07207616 0.57647059 0.07207616 0.28662930
## 22 23 24 28 33 35 38
## 0.07207616 0.28662930 0.07207616 0.28662930 0.07207616 0.28662930 0.33587786
## 41 42 53 54 57 60
## 0.07207616 0.07207616 0.28662930 0.72916667 0.07207616 0.07207616
Veamos graficamente como se ve:
plot(rp_predict~test$Churn)
Probamos cortes con la funcion umbrales, visualizamos la tabla respectiva.
umb_rp<-umbrales(test$Churn,rp_predict)
umb_rp
## umbral acierto precision cobertura F1
## 1 0.05 0.05000 0.05000 0.05000 0.05000
## 2 0.10 66.84011 43.92361 90.19608 59.07764
## 3 0.15 66.84011 43.92361 90.19608 59.07764
## 4 0.20 66.84011 43.92361 90.19608 59.07764
## 5 0.25 66.84011 43.92361 90.19608 59.07764
## 6 0.30 77.43614 56.46154 65.41889 60.61107
## 7 0.35 79.37559 63.10273 53.65419 57.99615
## 8 0.40 79.37559 63.10273 53.65419 57.99615
## 9 0.45 79.89593 67.43590 46.88057 55.31020
## 10 0.50 79.89593 67.43590 46.88057 55.31020
## 11 0.55 79.80132 67.91444 45.27629 54.33155
## 12 0.60 78.80795 74.04255 31.01604 43.71859
## 13 0.65 78.80795 74.04255 31.01604 43.71859
## 14 0.70 78.80795 74.04255 31.01604 43.71859
## 15 0.75 0.75000 0.75000 0.75000 0.75000
## 16 0.80 0.80000 0.80000 0.80000 0.80000
## 17 0.85 0.85000 0.85000 0.85000 0.85000
## 18 0.90 0.90000 0.90000 0.90000 0.90000
## 19 0.95 0.95000 0.95000 0.95000 0.95000
Lo hallamos directamente con la siguiente sintaxis.
umbral_final_rp<-umb_rp[which.max(umb_rp$F1),1]
umbral_final_rp
## [1] 0.3
Hallamos la matriz de confusion del umbral optimizado.
confusion(test$Churn,rp_predict,0.3)
##
## real FALSE TRUE
## 0 1270 283
## 1 194 367
rp_metricas<-filter(umb_rp,umbral==umbral_final_rp)
rp_metricas
## umbral acierto precision cobertura F1
## 1 0.3 77.43614 56.46154 65.41889 60.61107
Evaluamos la ROC
#creamos el objeto prediction
rp_prediction<-prediction(rp_predict,test$Churn)
#visualizamos la ROC
roc(rp_prediction)
Sacamos las metricas definitivas incluyendo el AUC
rp_metricas<-cbind(rp_metricas,AUC=round(auc(rp_prediction),2)*100)
print(t(rp_metricas))
## [,1]
## umbral 0.30000
## acierto 77.43614
## precision 56.46154
## cobertura 65.41889
## F1 60.61107
## AUC 82.00000
mod_rf <- rand_forest(
mode = "classification",
mtry = tune(),
trees = 500,
min_n = tune()
) %>%
set_engine(engine = "randomForest",importance=T)
set.seed(3001)
cv_folds <- vfold_cv(
data = train,
v = 5,
strata = Churn)
mod_rcp <- recipe(Churn~.,
data = train)
wf_rf <- workflow()%>%
add_recipe(mod_rcp)%>%
add_model(mod_rf)
hiperpar_grid=grid_max_entropy(
# Rango de busqueda para cada hiperparametro
mtry(range = c(1L, 10L), trans = NULL),
min_n(range = c(2L, 100L), trans = NULL),
# Numero de combinaciones totales
size = 100
)
registerDoParallel(cores = parallel::detectCores() - 1)
set.seed(2051)
rf_fit <- wf_rf%>%tune_grid(
resamples = cv_folds,
metrics = metric_set(roc_auc),
control = control_resamples(save_pred = TRUE),
grid = hiperpar_grid)
stopImplicitCluster()
show_best(rf_fit, metric = "roc_auc")
## # A tibble: 5 x 8
## mtry min_n .metric .estimator mean n std_err .config
## <int> <int> <chr> <chr> <dbl> <int> <dbl> <chr>
## 1 1 42 roc_auc binary 0.827 5 0.00325 Preprocessor1_Model70
## 2 1 19 roc_auc binary 0.826 5 0.00346 Preprocessor1_Model80
## 3 1 92 roc_auc binary 0.826 5 0.00347 Preprocessor1_Model89
## 4 1 81 roc_auc binary 0.825 5 0.00318 Preprocessor1_Model77
## 5 2 46 roc_auc binary 0.824 5 0.00354 Preprocessor1_Model96
mejor_hiper3 <- select_best(rf_fit, metric="roc_auc")
mejor_hiper3
## # A tibble: 1 x 3
## mtry min_n .config
## <int> <int> <chr>
## 1 1 42 Preprocessor1_Model70
wf_finalrf <-
wf_rf%>%
finalize_workflow(mejor_hiper3)
set.seed(3100)
final_model_rf <- fit(object = wf_finalrf, data = train)
model_xtr <- final_model_rf%>%
extract_fit_engine()
Visualizamos las variables de nuestros datos mas importantes en un grafico.
varImpPlot(model_xtr)
Como hay 2 criterios, creamos una unica variable agregada para tener una mejor idea de la importancia de cada variable.
importancia <- randomForest::importance(model_xtr)[,3:4]
#normalizamos para poner las dos variables en la misma escala. los valores negativos son las que menos predicen y los positivos las que mas
importancia_norm <- as.data.frame(scale(importancia))
#creamos una única variable como suma de las otras
importancia_norm <- importancia_norm %>% mutate(
Variable = rownames(importancia_norm),
Imp_tot = MeanDecreaseAccuracy + MeanDecreaseGini) %>%
mutate(Imp_tot = Imp_tot + abs(min(Imp_tot))) %>%
arrange(desc(Imp_tot)) %>%
select(Variable,Imp_tot,MeanDecreaseAccuracy,MeanDecreaseGini)
#hacemos un gráfico para ver la curva de caÃÂda de importancia
ggplot(importancia_norm, aes(reorder(Variable,-Imp_tot),Imp_tot)) + geom_bar(stat = "identity") + theme(axis.text.x = element_text(angle = 90,size = 7))
importancia_norm
## Variable Imp_tot MeanDecreaseAccuracy
## TENURE_DISC TENURE_DISC 5.4236211 1.49575077
## Contract Contract 5.2713893 0.67705055
## PaymentMethod PaymentMethod 4.5010238 0.79539834
## InternetService InternetService 4.2146883 0.71651601
## MONTHLYCHARGES_DISC MONTHLYCHARGES_DISC 3.5369080 0.26410281
## OnlineSecurity OnlineSecurity 2.3577771 0.28371930
## PaperlessBilling PaperlessBilling 2.0165402 -0.19271414
## TechSupport TechSupport 2.0141544 -0.09635295
## OnlineBackup OnlineBackup 1.0460069 -0.66157242
## Dependents Dependents 0.6261723 -1.37779260
## Partner Partner 0.0000000 -1.90410567
## MeanDecreaseGini
## TENURE_DISC 1.1089356
## Contract 1.7754041
## PaymentMethod 0.8866908
## InternetService 0.6792376
## MONTHLYCHARGES_DISC 0.4538705
## OnlineSecurity -0.7448769
## PaperlessBilling -0.6096803
## TechSupport -0.7084273
## OnlineBackup -1.1113553
## Dependents -0.8149698
## Partner -0.9148290
Vemos el grafico y nuestros datos en la tabla de forma decreciente segun nivel de importancia, y vemos donde puede estar un corte claro diferencia de una con respecto a la anterior, y vemos de la variable ‘TechSupport’ a ‘OnlineBackup’ un corte interesante.
Vamos a deseleccionar de nuestros datos esas 3 variables sobrantes por debajo del corte de ‘TechSupport’ y nos quedaremos con esta y las de la parte superior.
dt5 <- dt3%>%select(-OnlineBackup,-Dependents,-Partner)
y ahora vamos a realizar de nuevo el modelo con las nuevas variables.
set.seed(5001)
split_inicial <- initial_split(
data = dt5,
prop = 0.7,
strata = Churn
)
train <- training(split_inicial)
test <- testing(split_inicial)
mod_rf <- rand_forest(
mode = "classification",
mtry = tune(),
trees = 500,
min_n = tune()
) %>%
set_engine(engine = "randomForest",importance=T)
mod_rcp <- recipe(Churn~.,
data = train)
set.seed(3001)
cv_folds <- vfold_cv(
data = train,
v = 5,
strata = Churn)
wf_rf <- workflow()%>%
add_recipe(mod_rcp)%>%
add_model(mod_rf)
hiperpar_grid=grid_max_entropy(
# Rango de busqueda para cada hiperparametro
mtry(range = c(1L, 10L), trans = NULL),
min_n(range = c(2L, 100L), trans = NULL),
# Numero de combinaciones totales
size = 100
)
registerDoParallel(cores = parallel::detectCores() - 1)
set.seed(2051)
rf_fit <- wf_rf%>%tune_grid(
resamples = cv_folds,
metrics = metric_set(roc_auc),
control = control_resamples(save_pred = TRUE),
grid = hiperpar_grid)
stopImplicitCluster()
show_best(rf_fit, metric = "roc_auc")
## # A tibble: 5 x 8
## mtry min_n .metric .estimator mean n std_err .config
## <int> <int> <chr> <chr> <dbl> <int> <dbl> <chr>
## 1 2 21 roc_auc binary 0.833 5 0.00664 Preprocessor1_Model60
## 2 2 37 roc_auc binary 0.831 5 0.00631 Preprocessor1_Model66
## 3 1 81 roc_auc binary 0.831 5 0.00611 Preprocessor1_Model77
## 4 2 29 roc_auc binary 0.831 5 0.00591 Preprocessor1_Model50
## 5 2 46 roc_auc binary 0.831 5 0.00552 Preprocessor1_Model96
mejor_hiper4 <- select_best(rf_fit, metric="roc_auc")
mejor_hiper4
## # A tibble: 1 x 3
## mtry min_n .config
## <int> <int> <chr>
## 1 2 21 Preprocessor1_Model60
wf_finalrf <-
wf_rf%>%
finalize_workflow(mejor_hiper4)
set.seed(3100)
final_model_rf <- fit(object = wf_finalrf, data = train)
model_extr3 <- final_model_rf%>%
extract_fit_engine()
Calculo de prediccion del modelo (scorings).
rf_predict <- predict(model_extr3, test, type = 'prob')[,2]
head(rf_predict,20)
## 1 3 5 7 8 9 13 15 16 17 25 26 27
## 0.798 0.150 0.988 0.382 0.054 0.514 0.012 0.286 0.000 0.000 0.012 0.054 0.508
## 28 31 34 37 40 42 45
## 0.584 0.000 0.024 0.988 0.428 0.000 0.068
Veamos graficamente como se ve:
plot(rf_predict~test$Churn)
Probamos cortes con la funcion umbrales, visualizamos la tabla respectiva.
umb_rf<-umbrales(test$Churn,rf_predict)
umb_rf
## umbral acierto precision cobertura F1
## 1 0.05 72.18543 48.53420 79.67914 60.32389
## 2 0.10 75.82781 53.22997 73.44029 61.72285
## 3 0.15 77.05771 55.41311 69.34046 61.59937
## 4 0.20 78.28761 57.79817 67.37968 62.22222
## 5 0.25 78.99716 59.73378 63.99287 61.79002
## 6 0.30 79.37559 60.98418 61.85383 61.41593
## 7 0.35 79.75402 62.29205 60.07130 61.16152
## 8 0.40 79.65941 63.50515 54.90196 58.89101
## 9 0.45 79.94324 64.73118 53.65419 58.67446
## 10 0.50 79.99054 65.06550 53.11943 58.48871
## 11 0.55 79.37559 66.06684 45.81105 54.10526
## 12 0.60 79.56481 68.69565 42.24599 52.31788
## 13 0.65 79.28098 68.35821 40.81996 51.11607
## 14 0.70 79.32829 68.67470 40.64171 51.06383
## 15 0.75 79.13907 69.10828 38.68093 49.60000
## 16 0.80 78.85525 69.93007 35.65062 47.22550
## 17 0.85 78.85525 70.21277 35.29412 46.97509
## 18 0.90 77.38884 75.77640 21.74688 33.79501
## 19 0.95 76.77389 75.73529 18.36007 29.55524
Lo hallamos directamente con la siguiente sintaxis.
umbral_final_rf<-umb_rf[which.max(umb_rf$F1),1]
umbral_final_rf
## [1] 0.2
Hallamos la matriz de confusion del umbral optimizado.
confusion(test$Churn,rf_predict,0.2)
##
## real FALSE TRUE
## 0 1277 276
## 1 183 378
rf_metricas<-filter(umb_rf,umbral==umbral_final_rf)
rf_metricas
## umbral acierto precision cobertura F1
## 1 0.2 78.28761 57.79817 67.37968 62.22222
Evaluamos la ROC
#creamos el objeto prediction
rf_prediction<-prediction(rf_predict,test$Churn)
#visualizamos la ROC
roc(rf_prediction)
Sacamos las metricas definitivas incluyendo el AUC
rf_metricas<-cbind(rf_metricas,AUC=round(auc(rf_prediction),2)*100)
print(t(rf_metricas))
## [,1]
## umbral 0.20000
## acierto 78.28761
## precision 57.79817
## cobertura 67.37968
## F1 62.22222
## AUC 83.00000
set.seed(6000)
split_inicial <- initial_split(
data = dt3,
prop = 0.7,
strata = Churn
)
train <- training(split_inicial)
test <- testing(split_inicial)
mod_xg <- boost_tree(mtry = tune(),
min_n = tune(),
tree_depth = tune(),
trees = 500,
learn_rate = 0.1,
sample_size = 0.8) %>%
set_engine("xgboost") %>%
set_mode("classification")
set.seed(3001)
cv_folds <- vfold_cv(
data = train,
v = 5,
strata = Churn)
wf_xg<- workflow()%>%
add_formula(Churn~.)%>%
add_model(mod_xg)
# Parametros para tuning
params_xgb <- parameters(
finalize(mtry(), x = train[, -1]),
min_n(range = c(2L, 50L)),
tree_depth(range = c(3L, 8L))
)
# Grid
set.seed(2000)
grid_xgb <- params_xgb %>%
grid_max_entropy(size = 10)
registerDoParallel(cores = parallel::detectCores() - 1)
set.seed(2631)
xg_fit <- wf_xg%>%tune_grid(
resamples = cv_folds,
metrics = metric_set(roc_auc),
grid = grid_xgb,
control = control_grid(save_pred = TRUE))
stopImplicitCluster()
show_best(xg_fit, metric = "roc_auc")
## # A tibble: 5 x 9
## mtry min_n tree_depth .metric .estimator mean n std_err .config
## <int> <int> <int> <chr> <chr> <dbl> <int> <dbl> <chr>
## 1 1 37 7 roc_auc binary 0.840 5 0.00210 Preprocessor1_M~
## 2 2 36 4 roc_auc binary 0.838 5 0.00241 Preprocessor1_M~
## 3 1 8 3 roc_auc binary 0.838 5 0.00240 Preprocessor1_M~
## 4 7 48 8 roc_auc binary 0.837 5 0.00186 Preprocessor1_M~
## 5 9 29 3 roc_auc binary 0.835 5 0.00188 Preprocessor1_M~
mejor_hiper5 <- select_best(xg_fit, metric="roc_auc")
mejor_hiper5
## # A tibble: 1 x 4
## mtry min_n tree_depth .config
## <int> <int> <int> <chr>
## 1 1 37 7 Preprocessor1_Model10
wf_finalxg <-
wf_rf%>%
finalize_workflow(mejor_hiper5)
set.seed(3300)
final_model_xg <- fit(object = wf_finalxg, data = train)
model_extr4 <- final_model_xg%>%
extract_fit_engine()
Calculo de prediccion del modelo (scorings).
xg_predict <- predict(model_extr4, test, type = 'prob')[,2]
head(xg_predict,20)
## 1 2 3 4 6 8 9 17 19 24 25 28 31
## 0.486 0.000 0.118 0.000 0.766 0.058 0.316 0.000 0.128 0.002 0.002 0.344 0.000
## 32 37 39 47 50 51 58
## 0.534 0.766 0.472 0.174 0.000 0.472 0.036
Veamos graficamente como se ve:
plot(xg_predict~test$Churn)
Probamos cortes con la funcion umbrales, visualizamos la tabla respectiva.
umb_xg<-umbrales(test$Churn,xg_predict)
umb_xg
## umbral acierto precision cobertura F1
## 1 0.05 69.77294 46.14625 83.24421 59.37699
## 2 0.10 75.07096 52.16837 72.90553 60.81784
## 3 0.15 77.57805 56.50224 67.37968 61.46341
## 4 0.20 78.94986 60.35714 60.24955 60.30330
## 5 0.25 78.61873 60.96579 54.01070 57.27788
## 6 0.30 79.13907 63.27434 50.98039 56.46594
## 7 0.35 79.47020 65.91479 46.88057 54.79167
## 8 0.40 79.51750 66.84211 45.27629 53.98512
## 9 0.45 79.61211 68.36158 43.13725 52.89617
## 10 0.50 79.32829 72.79412 35.29412 47.53902
## 11 0.55 78.66604 72.17742 31.90731 44.25216
## 12 0.60 77.15232 75.65789 20.49911 32.25806
## 13 0.65 76.72658 76.74419 17.64706 28.69565
## 14 0.70 76.72658 76.74419 17.64706 28.69565
## 15 0.75 76.72658 76.74419 17.64706 28.69565
## 16 0.80 0.80000 0.80000 0.80000 0.80000
## 17 0.85 0.85000 0.85000 0.85000 0.85000
## 18 0.90 0.90000 0.90000 0.90000 0.90000
## 19 0.95 0.95000 0.95000 0.95000 0.95000
Lo hallamos directamente con la siguiente sintaxis.
umbral_final_xg<-umb_xg[which.max(umb_xg$F1),1]
umbral_final_xg
## [1] 0.15
Hallamos la matriz de confusion del umbral optimizado.
confusion(test$Churn,xg_predict,0.15)
##
## real FALSE TRUE
## 0 1261 292
## 1 183 378
xg_metricas<-filter(umb_xg,umbral==umbral_final_xg)
xg_metricas
## umbral acierto precision cobertura F1
## 1 0.15 77.57805 56.50224 67.37968 61.46341
Evaluamos la ROC
#creamos el objeto prediction
xg_prediction<-prediction(xg_predict,test$Churn)
#visualizamos la ROC
roc(xg_prediction)
Sacamos las metricas definitivas incluyendo el AUC
xg_metricas<-cbind(xg_metricas,AUC=round(auc(xg_prediction),2)*100)
print(t(xg_metricas))
## [,1]
## umbral 0.15000
## acierto 77.57805
## precision 56.50224
## cobertura 67.37968
## F1 61.46341
## AUC 82.00000
comparativa <- rbind(rl_metricas,rp_metricas,rf_metricas,xg_metricas)
rownames(comparativa) <- c('Regresion Logistica', 'Arbol Decision', 'Random Forest', 'Extreme Gradient Boosting')
t(comparativa)
## Regresion Logistica Arbol Decision Random Forest
## umbral 0.30000 0.30000 0.20000
## acierto 76.96310 77.43614 78.28761
## precision 54.75578 56.46154 57.79817
## cobertura 75.93583 65.41889 67.37968
## F1 63.62957 60.61107 62.22222
## AUC 84.00000 82.00000 83.00000
## Extreme Gradient Boosting
## umbral 0.15000
## acierto 77.57805
## precision 56.50224
## cobertura 67.37968
## F1 61.46341
## AUC 82.00000
rl_auc <- rlog%>%collect_predictions(summarize=TRUE,parameters=best_rlog2)%>%
roc_curve(Churn,.pred_0)%>%
mutate(model = "Regresion Logistica")
rp_auc <- rp_fit%>%collect_predictions(summarize=TRUE,parameters=mejor_hiper2)%>%
roc_curve(Churn,.pred_0)%>%
mutate(model = "Arbol de Decision")
rf_auc <- rf_fit%>%collect_predictions(summarize=TRUE,parameters=mejor_hiper4)%>%
roc_curve(Churn,.pred_0)%>%
mutate(model = "Random Forest")
xg_auc <- xg_fit%>%collect_predictions(summarize=TRUE,parameters=mejor_hiper5)%>%
roc_curve(Churn,.pred_0)%>%
mutate(model = "Xgboost")
bind_rows(rl_auc,rp_auc,rf_auc,xg_auc)%>% # Dibuja las 4 curvas AUC juntas
ggplot(aes(x=1-specificity, y = sensitivity, col = model))+ # Especifica el eje x e y dibuja la columna para usar el nombre de la metrica
geom_path(lwd = 1.5, alpha = 0.8) + # Conecta a las 4 AUC, lwd = anchura de linea, alpha = Color transparencia del valor.
geom_abline(lty = 3) + # abline del dibujo, lty= linea tipo
coord_equal() + # Asegura que los rangos de los ejes sean iguales
scale_color_viridis_d(option = "plasma", end = .6)
En base a la comparativa el modelo de Regresion Logistica ha obtenido levemente mejores resultados, asi que seleccionamos este modelo ganador en la comparativa.
Añadimos la variable scoring a nuestro dataframe en base a la prediccion del modelo seleccionado y guardamos el modelo final.
dt1$SCORING_CHURN <- predict(model_extr1,dt1,type='response')
# Visualizamos tabla con el scoring reflejado:
kable(head(dt1,6), booktabs = T) %>%
kable_styling(font_size=7)
| customerID | Contract | InternetService | PaymentMethod | PaperlessBilling | TechSupport | OnlineSecurity | Dependents | Partner | OnlineBackup | MONTHLYCHARGES_DISC | TENURE_DISC | Churn | SCORING_CHURN |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| 7590-VHVEG | Month-to-month | DSL | Electronic check | Yes | No | No | No | Yes | Yes | 02 <= 55.95 | 01_MENOR_10 | 0 | 0.5723726 |
| 5575-GNVDE | One year | DSL | Mailed check | No | No | Yes | No | No | No | 03 <= 68.8 | 03_DE_30_A_45 | 0 | 0.0468170 |
| 3668-QPYBK | Month-to-month | DSL | Mailed check | Yes | No | Yes | No | No | Yes | 02 <= 55.95 | 01_MENOR_10 | 1 | 0.3459541 |
| 7795-CFOCW | One year | DSL | Bank transfer (automatic) | No | Yes | Yes | No | No | No | 02 <= 55.95 | 03_DE_30_A_45 | 0 | 0.0508932 |
| 9237-HQITU | Month-to-month | Fiber optic | Electronic check | Yes | No | No | No | No | No | 04 <= 106.75 | 01_MENOR_10 | 1 | 0.7750630 |
| 9305-CDSKC | Month-to-month | Fiber optic | Electronic check | Yes | No | No | No | No | No | 04 <= 106.75 | 01_MENOR_10 | 1 | 0.7750630 |
#Hacemos los guardados que incluyen al modelo final.
saveRDS(model_extr1,'modelo_final.rds')
saveRDS(dt1,'cache4.rds')
Vamos a visualizar el abandono real por tramos de scoring. Este grafico es muy potente para ver que el modelo es consistente, ya que debe presentar una linea descendente en la tasa de abandono conforme se desciende en el scoring.
Creamos una funcion para visualizar la contratacion real por percentiles de scoring.
vis <- function(scoring,real) {
#Preparar el dataframe de visualizacion
vis_df <- data.frame(Scoring = scoring, Perc_Scoring = cut_number(scoring, 10), Real = real)
levels(vis_df$Perc_Scoring) <- seq(from = 100,to = 10,by = -10)
vis_gr <- vis_df %>% group_by(Perc_Scoring) %>% summarise(Tasa_Contr = mean(as.numeric(as.character(Real)))) %>% arrange(Perc_Scoring)
#Ordenar el factor para el grafico
vis_gr$Perc_Scoring <- factor(vis_gr$Perc_Scoring, levels = vis_gr$Perc_Scoring[order(vis_gr$Perc_Scoring, decreasing = T)])
#hacemos el grafico
ggplot(vis_gr,aes(Perc_Scoring, Tasa_Contr)) +
geom_col(fill='grey') +
geom_hline(aes(yintercept = mean(as.numeric(as.character(vis_df$Real)))), col = 'black') +
labs(title = 'Abandono real por tramo de scoring', x = 'Tramo de Scoring', y = 'Abandono Real')
}
vis(dt1$SCORING_CHURN,dt1$Churn)
Como vemos el modelo esta bien calibrado pues se aprecia esa tendencia descendente de la tasa de abandono conforme disminuye el scoring.
Vemos tambien en la grafica que a partir del 40% el abandono es inferior a la media, a partir de este punto no nos interesaria lanzar las campañas, aunque para asegurarnos un mejor ROI, se tratara de coger tramos anteriores a ese 40%, para mayor confiabilidad de los resultados.
Nuestro caso sera lanzar una campaña de fidelizacion para los clientes haciendole diversas promociones y descuentos incentivando el engagement con el cliente.
La campaña tiene un presupuesto de 23000 euros, un coste unitario de 25 euros por cliente y sera realizada por call center.
Haciendo las cuentas el numero de clientes a contactar sera 23000/25=920
Filtramos pues los clientes a los que vamos a lanzar la campaña para intentar fidelizarlos.
y escojemos estos 920 clientes primeros ordenados por scoring que nos han salido al ajustar el coste del cliente en funcion del tamaño de la campaña.
y visualizamos los 50 primeros.
tamaño_campaña <- 920
bote_campaña <- dt1%>%
filter(Churn==0)%>%
arrange(desc(SCORING_CHURN))%>%
slice(1:tamaño_campaña)%>%
select(customerID,SCORING_CHURN)
#Previsualizamos la salida
head(bote_campaña,50)
## customerID SCORING_CHURN
## 1: 4482-EWFMI 0.775063
## 2: 4847-TAJYI 0.775063
## 3: 8266-VBFQL 0.775063
## 4: 2799-ARNLO 0.775063
## 5: 0021-IKXGC 0.775063
## 6: 4115-NZRKS 0.775063
## 7: 5168-MSWXT 0.775063
## 8: 1452-VOQCH 0.775063
## 9: 4234-XTNEA 0.775063
## 10: 8270-RKSAP 0.775063
## 11: 6786-OBWQR 0.775063
## 12: 2123-AGEEN 0.775063
## 13: 6630-UJZMY 0.775063
## 14: 8242-SOQUO 0.775063
## 15: 0455-XFASS 0.775063
## 16: 7439-DKZTW 0.775063
## 17: 7363-QTBIW 0.775063
## 18: 6077-BDPXA 0.775063
## 19: 4291-SHSBH 0.775063
## 20: 9957-YODKZ 0.775063
## 21: 7321-VGNKU 0.775063
## 22: 1219-NNDDO 0.775063
## 23: 6435-SRWBJ 0.775063
## 24: 2081-VEYEH 0.775063
## 25: 3030-YDNRM 0.775063
## 26: 3318-NMQXL 0.775063
## 27: 2018-QKYGT 0.775063
## 28: 1960-UYCNN 0.775063
## 29: 8087-LGYHQ 0.775063
## 30: 1393-IMKZG 0.775063
## 31: 4813-HQMGZ 0.775063
## 32: 7245-NIIWQ 0.775063
## 33: 8739-XNIKG 0.775063
## 34: 2789-HQBOU 0.775063
## 35: 9603-OAIHC 0.775063
## 36: 1475-VWVDO 0.775063
## 37: 1197-BVMVG 0.775063
## 38: 8622-ZLFKO 0.775063
## 39: 4927-WWOOZ 0.775063
## 40: 0722-SVSFK 0.775063
## 41: 8325-QRPZR 0.775063
## 42: 5542-TBBWB 0.775063
## 43: 4770-UEZOX 0.775063
## 44: 6860-YRJZP 0.775063
## 45: 5150-ITWWB 0.775063
## 46: 6818-DJXAA 0.775063
## 47: 9643-AVVWI 0.775063
## 48: 1640-PLFMP 0.775063
## 49: 0840-DFEZH 0.775063
## 50: 2545-EBUPK 0.775063
## customerID SCORING_CHURN
Vamos a ver graficamente si de esta forma estamos aprovechando el potencial de nuestro modelo
penetracion_target <- mean(as.numeric(as.character(dt1$Churn)))
dt1 %>%
arrange(desc(SCORING_CHURN)) %>%
ggplot(aes(y = SCORING_CHURN, x = seq_along(SCORING_CHURN))) +
geom_line() +
geom_vline(xintercept = tamaño_campaña, col = 'orange') +
geom_hline(yintercept = penetracion_target,col='blue') +
labs(x = 'CLIENTES ORDENADOS POR SCORING', y = 'SCORING')
Vemos que tenemos clientes cercanos a 3000 con un scorring mayor de la media(punto de corte linea azul hacia la izquierda), a partir de ahi hacia la derecha tendrian un scoring inferior a la media, si tomamos el triangulo que dibuja donde intersectan las rectas naranja del tamaño de la campaña(en nuestro caso 920 clientes) y azul de penetracion de la media de la target con nuestra linea del grafico de descenso del scoring tenemos el recorrido de clientes posibles a los cuales lanzar la campaña y tambien vemos que podriamos tratar de ampliar a mayor tamaño de campaña porque tendriamos margen para ello.
Suele ser normal coger aquellos clientes cuyo scoring este x2 o x3 por encima de la penetracion de la target
scoring_x2 <- penetracion_target * 2
scoring_x3 <- penetracion_target * 3
bote_campaña_x2 <- dt1 %>%
filter(Churn==0 & SCORING_CHURN > scoring_x2) %>%
select(customerID,SCORING_CHURN)
#Tamaño del bote_campaña_x2
nrow(bote_campaña_x2)
## [1] 401
bote_campaña_x3 <- dt1 %>%
filter(Churn==0 & SCORING_CHURN > scoring_x3) %>%
select(customerID,SCORING_CHURN)
#Tamaño del bote_campaña_x3
nrow(bote_campaña_x3)
## [1] 0
Vamos a ver graficamente si de esta forma estamos aprovechando el potencial de nuestro modelo
dt1 %>%
arrange(desc(SCORING_CHURN)) %>%
ggplot(aes(y = SCORING_CHURN, x = seq_along(SCORING_CHURN))) +
geom_line() +
geom_hline(yintercept = scoring_x2,col='blue') +
geom_hline(yintercept = scoring_x3,col='orange') +
geom_hline(yintercept = penetracion_target,col='black') +
labs(x = 'CLIENTES ORDENADOS POR SCORING', y = 'SCORING')
Al hacer el estudio suponiendo que no tenemos restricciones presupuestarias para esta opcion concreta, vemos la recomendacion del numero de clientes(tamaño de campaña) que dariamos para lanzar la campaña analizando el scoring 2 veces por encima de la media de penetracion de la target (scoringx2) y el scoring 3 veces por encima de esa misma media(scoringx3), y al hacerlo vemos que tenemos 401 clientes para lanzar la campaña con resultados optimos con un nivel de scoringx2, en cambio no somos capaces de seleccionar clientes por encima de la media 3 veces(scoringx3), nos da cero(numero de clientes) el resultado, y se puede observar esto visualmente tambien en el grafico en el que por encima de la linea naranja (ni tan siquiera llega) no supera la linea de descenso del scoring esa franja.
Por lo tanto los maximos resultados optimos que obtenemos son los que ofrece el nivel de scoringx2.
Para ello tenemos en cuenta el ingreso medio previsto y el coste medio por accion comercial.
Asi pues calculamos el ingreso medio del cliente retenido.
clientes_retenidos <- dt%>%filter(Churn==0)
ingreso_medio_cliente <- mean(clientes_retenidos$MonthlyCharges)
ingreso_medio_cliente
## [1] 61.26512
suponiendo un margen neto de un 15%
margen_medio_anual=61.26*0.15*12
margen_medio_anual
## [1] 110.268
mar_medio <- margen_medio_anual
Supuesto de margen medio anual por cliente = 110.26€
Supuesto de coste medio por accion comercial (call center cliente contactado) = 25€
Definicion de margen esperado = probabilidad de evento * margen evento
Definicion de margen neto = margen esperado - coste medio
mar_medio <- 110.26
coste_medio <- 25
#Calculamos el margen esperado de cada cliente
df_campaña <- dt1 %>%
filter(Churn==0) %>%
mutate(
ME = SCORING_CHURN * mar_medio,
MN = ME - coste_medio) %>%
arrange(desc(MN)) %>%
mutate(INDICE = 1:nrow(.)) %>%
select(customerID,INDICE,ME,MN)
head(df_campaña,50)
## customerID INDICE ME MN
## 1: 4482-EWFMI 1 85.45845 60.45845
## 2: 4847-TAJYI 2 85.45845 60.45845
## 3: 8266-VBFQL 3 85.45845 60.45845
## 4: 2799-ARNLO 4 85.45845 60.45845
## 5: 0021-IKXGC 5 85.45845 60.45845
## 6: 4115-NZRKS 6 85.45845 60.45845
## 7: 5168-MSWXT 7 85.45845 60.45845
## 8: 1452-VOQCH 8 85.45845 60.45845
## 9: 4234-XTNEA 9 85.45845 60.45845
## 10: 8270-RKSAP 10 85.45845 60.45845
## 11: 6786-OBWQR 11 85.45845 60.45845
## 12: 2123-AGEEN 12 85.45845 60.45845
## 13: 6630-UJZMY 13 85.45845 60.45845
## 14: 8242-SOQUO 14 85.45845 60.45845
## 15: 0455-XFASS 15 85.45845 60.45845
## 16: 7439-DKZTW 16 85.45845 60.45845
## 17: 7363-QTBIW 17 85.45845 60.45845
## 18: 6077-BDPXA 18 85.45845 60.45845
## 19: 4291-SHSBH 19 85.45845 60.45845
## 20: 9957-YODKZ 20 85.45845 60.45845
## 21: 7321-VGNKU 21 85.45845 60.45845
## 22: 1219-NNDDO 22 85.45845 60.45845
## 23: 6435-SRWBJ 23 85.45845 60.45845
## 24: 2081-VEYEH 24 85.45845 60.45845
## 25: 3030-YDNRM 25 85.45845 60.45845
## 26: 3318-NMQXL 26 85.45845 60.45845
## 27: 2018-QKYGT 27 85.45845 60.45845
## 28: 1960-UYCNN 28 85.45845 60.45845
## 29: 8087-LGYHQ 29 85.45845 60.45845
## 30: 1393-IMKZG 30 85.45845 60.45845
## 31: 4813-HQMGZ 31 85.45845 60.45845
## 32: 7245-NIIWQ 32 85.45845 60.45845
## 33: 8739-XNIKG 33 85.45845 60.45845
## 34: 2789-HQBOU 34 85.45845 60.45845
## 35: 9603-OAIHC 35 85.45845 60.45845
## 36: 1475-VWVDO 36 85.45845 60.45845
## 37: 1197-BVMVG 37 85.45845 60.45845
## 38: 8622-ZLFKO 38 85.45845 60.45845
## 39: 4927-WWOOZ 39 85.45845 60.45845
## 40: 0722-SVSFK 40 85.45845 60.45845
## 41: 8325-QRPZR 41 85.45845 60.45845
## 42: 5542-TBBWB 42 85.45845 60.45845
## 43: 4770-UEZOX 43 85.45845 60.45845
## 44: 6860-YRJZP 44 85.45845 60.45845
## 45: 5150-ITWWB 45 85.45845 60.45845
## 46: 6818-DJXAA 46 85.45845 60.45845
## 47: 9643-AVVWI 47 85.45845 60.45845
## 48: 1640-PLFMP 48 85.45845 60.45845
## 49: 0840-DFEZH 49 85.45845 60.45845
## 50: 2545-EBUPK 50 85.45845 60.45845
## customerID INDICE ME MN
Visualizamos las curvas
Localizamos el punto en el que el margen neto pasa a ser cero o menos
MN_cero <- df_campaña %>% filter(MN <= 0 ) %>% slice(1) %>% select(INDICE)
MN_cero <- MN_cero$INDICE
Hacemos el grafico
ggplot(df_campaña,aes(x = INDICE)) +
geom_line(aes(y = ME, col = "ME")) +
geom_line(aes(y = MN, col = "MN")) +
geom_hline(aes(yintercept = coste_medio, col = 'COSTE MEDIO')) +
geom_vline(aes(xintercept = MN_cero, col = 'MARGEN NETO CERO')) +
labs(x = 'TAMAÑO DE CAMPAÑA', y = 'MARGEN', colour = 'KPI')
print(paste('Tamaño maximo de campaña rentable: ',MN_cero))
## [1] "Tamaño maximo de campaña rentable: 1688"
Vemos que segun los parametros que nosotros le hemos definido donde la linea vertical verde se cruza hacia abajo con la azul(ME) y la rosa(MN) que coincide con el valor ‘0’ de la linea del MN(rosa) al cruzarse, a partir de ese punto (toda la linea verde que marca el equilibrio) hacia la derecha perdemos dinero, hacia la izquierda de la linea verde ganamos dinero, asi pues determinamos ese punto donde intersecta la linea verde como tamaño maximo de campaña rentable, que luego imprimimos con la sentencia print para verlo, ya que no se muestra en numero en la grafica y seria ese el tamaño maximo de la campaña que resultara rentable para nuestro caso que serian 1688 clientes. Lanzar la campaña a mas de 1688 clientes seria empezar a perder dinero.Mientras que cortar en ese punto (1688 clientes), seria lo que nos optimizaria el tamaño de campaña(punto de equilibrio).
Vamos a calcular 2 nuevas variables que sean un agregado de los ingresos agregados y de los gastos agregados en cada potencial tamaño de campaña, y el ROI como diferencia de las anteriores, y vamos a localizar el tamaño de la campaña que va a maximizar ese ROI y tambien cuanto vamos a ganar previsiblemente
Vamos a usar la funcion cumsum(), que hace una suma acumulada secuencial de la variable que le pases como parametro
df_campaña <- df_campaña %>%
mutate(
INGRESOS_AGRE = cumsum(ME),
COSTES_AGRE = INDICE * coste_medio,
ROI = INGRESOS_AGRE - COSTES_AGRE)
head(df_campaña,50)
## customerID INDICE ME MN INGRESOS_AGRE COSTES_AGRE ROI
## 1: 4482-EWFMI 1 85.45845 60.45845 85.45845 25 60.45845
## 2: 4847-TAJYI 2 85.45845 60.45845 170.91689 50 120.91689
## 3: 8266-VBFQL 3 85.45845 60.45845 256.37534 75 181.37534
## 4: 2799-ARNLO 4 85.45845 60.45845 341.83378 100 241.83378
## 5: 0021-IKXGC 5 85.45845 60.45845 427.29223 125 302.29223
## 6: 4115-NZRKS 6 85.45845 60.45845 512.75067 150 362.75067
## 7: 5168-MSWXT 7 85.45845 60.45845 598.20912 175 423.20912
## 8: 1452-VOQCH 8 85.45845 60.45845 683.66756 200 483.66756
## 9: 4234-XTNEA 9 85.45845 60.45845 769.12601 225 544.12601
## 10: 8270-RKSAP 10 85.45845 60.45845 854.58445 250 604.58445
## 11: 6786-OBWQR 11 85.45845 60.45845 940.04290 275 665.04290
## 12: 2123-AGEEN 12 85.45845 60.45845 1025.50135 300 725.50135
## 13: 6630-UJZMY 13 85.45845 60.45845 1110.95979 325 785.95979
## 14: 8242-SOQUO 14 85.45845 60.45845 1196.41824 350 846.41824
## 15: 0455-XFASS 15 85.45845 60.45845 1281.87668 375 906.87668
## 16: 7439-DKZTW 16 85.45845 60.45845 1367.33513 400 967.33513
## 17: 7363-QTBIW 17 85.45845 60.45845 1452.79357 425 1027.79357
## 18: 6077-BDPXA 18 85.45845 60.45845 1538.25202 450 1088.25202
## 19: 4291-SHSBH 19 85.45845 60.45845 1623.71046 475 1148.71046
## 20: 9957-YODKZ 20 85.45845 60.45845 1709.16891 500 1209.16891
## 21: 7321-VGNKU 21 85.45845 60.45845 1794.62736 525 1269.62736
## 22: 1219-NNDDO 22 85.45845 60.45845 1880.08580 550 1330.08580
## 23: 6435-SRWBJ 23 85.45845 60.45845 1965.54425 575 1390.54425
## 24: 2081-VEYEH 24 85.45845 60.45845 2051.00269 600 1451.00269
## 25: 3030-YDNRM 25 85.45845 60.45845 2136.46114 625 1511.46114
## 26: 3318-NMQXL 26 85.45845 60.45845 2221.91958 650 1571.91958
## 27: 2018-QKYGT 27 85.45845 60.45845 2307.37803 675 1632.37803
## 28: 1960-UYCNN 28 85.45845 60.45845 2392.83647 700 1692.83647
## 29: 8087-LGYHQ 29 85.45845 60.45845 2478.29492 725 1753.29492
## 30: 1393-IMKZG 30 85.45845 60.45845 2563.75336 750 1813.75336
## 31: 4813-HQMGZ 31 85.45845 60.45845 2649.21181 775 1874.21181
## 32: 7245-NIIWQ 32 85.45845 60.45845 2734.67026 800 1934.67026
## 33: 8739-XNIKG 33 85.45845 60.45845 2820.12870 825 1995.12870
## 34: 2789-HQBOU 34 85.45845 60.45845 2905.58715 850 2055.58715
## 35: 9603-OAIHC 35 85.45845 60.45845 2991.04559 875 2116.04559
## 36: 1475-VWVDO 36 85.45845 60.45845 3076.50404 900 2176.50404
## 37: 1197-BVMVG 37 85.45845 60.45845 3161.96248 925 2236.96248
## 38: 8622-ZLFKO 38 85.45845 60.45845 3247.42093 950 2297.42093
## 39: 4927-WWOOZ 39 85.45845 60.45845 3332.87937 975 2357.87937
## 40: 0722-SVSFK 40 85.45845 60.45845 3418.33782 1000 2418.33782
## 41: 8325-QRPZR 41 85.45845 60.45845 3503.79626 1025 2478.79626
## 42: 5542-TBBWB 42 85.45845 60.45845 3589.25471 1050 2539.25471
## 43: 4770-UEZOX 43 85.45845 60.45845 3674.71316 1075 2599.71316
## 44: 6860-YRJZP 44 85.45845 60.45845 3760.17160 1100 2660.17160
## 45: 5150-ITWWB 45 85.45845 60.45845 3845.63005 1125 2720.63005
## 46: 6818-DJXAA 46 85.45845 60.45845 3931.08849 1150 2781.08849
## 47: 9643-AVVWI 47 85.45845 60.45845 4016.54694 1175 2841.54694
## 48: 1640-PLFMP 48 85.45845 60.45845 4102.00538 1200 2902.00538
## 49: 0840-DFEZH 49 85.45845 60.45845 4187.46383 1225 2962.46383
## 50: 2545-EBUPK 50 85.45845 60.45845 4272.92227 1250 3022.92227
## customerID INDICE ME MN INGRESOS_AGRE COSTES_AGRE ROI
Visualizamos las curvas
ggplot(df_campaña,aes(x = INDICE)) +
geom_line(aes(y = INGRESOS_AGRE, col='INGRESOS_AGRE')) +
geom_line(aes(y = COSTES_AGRE, col='COSTES_AGRE')) +
geom_line(aes(y = ROI, col='ROI')) +
labs(y='EUROS', x = 'TAMAÑO DE CAMPAÑA')
Lo que hago a posteriori es amplificar la visualizacion de la linea del ROI azul hasta el punto de la curva en el que el ROI=0 a partir de ahi hacia la derecha cada accion comercial que realicemos estaremos perdiendo dinero. Con ello pretendemos dilucidar el pto optimo de la curva (que va a coincidir con el hallado en la opcion 3), la diferencia es que de forma añadida podremos entre otros datos sacar el ROI de la campaña en ese punto optimo que es lo que estamos buscando y no tenemos.
Vamos por tanto a visualizar un zoom sobre el ROI solo en los tamaños de campaña que son positivos para localizar el punto optimo
df_campaña %>%
filter(ROI > 0) %>%
ggplot(aes(x = INDICE)) +
geom_line(aes(y = ROI, col='ROI')) +
geom_vline(aes(xintercept = MN_cero, col = 'PUNTO OPTIMO')) +
labs(x = 'TAMAÑO DE CAMPAÑA', y = 'ROI', colour = 'KPI')
Esta curva(linea azul del ROI)interceptara con la linea roja vertical que corresponde al MN_cero para darnos el punto optimo descrito.
Imprimimos un pequeño informe donde nos dara el resultado final con los datos totales requeridos de forma resumida.
cat(
paste(
'El tamaño optimo de campaña para el ROI es de:', MN_cero, 'clientes',
'\nCon unos ingresos esperados de margen neto acumulado de:', round(df_campaña[which(df_campaña$INDICE == MN_cero),'INGRESOS_AGRE']), '€',
'\nY unos costes agregados de:',
df_campaña[which(df_campaña$INDICE == MN_cero),'COSTES_AGRE'], '€',
'\nQue van a generar un Retorno Neto de la Inversion de:',
round(df_campaña[which(df_campaña$INDICE == MN_cero),'ROI']),'€'
)
)
## El tamaño optimo de campaña para el ROI es de: 1688 clientes
## Con unos ingresos esperados de margen neto acumulado de: 78985 €
## Y unos costes agregados de: 42200 €
## Que van a generar un Retorno Neto de la Inversion de: 36785 €