Realizar e interpretar regresión logísitica con datos de personas e ingresos de USA
Construir un modelo de regresión logísitca aplicado a datos de personas y sus ingresos en USA La variable dependiente es los ingresos identificado por 0 y 1, los ganan por debajo o igual a 50 Mil y los que ganan por encima de 50 Mil.
Posiblemente se utilicen algunas de ellas
library(ggplot2)
library(dplyr)
##
## Attaching package: 'dplyr'
## The following objects are masked from 'package:stats':
##
## filter, lag
## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, union
library(knitr)
library(caret)
## Warning: package 'caret' was built under R version 4.0.3
## Loading required package: lattice
library(readr)
library(knitr)
library(DT) # install.packages('DT')
## Warning: package 'DT' was built under R version 4.0.3
library(scales)
##
## Attaching package: 'scales'
## The following object is masked from 'package:readr':
##
## col_factor
datos <- read.csv('C:/Users/Blue/Documents/adultos_clean.csv', encoding = "UTF-8")
datatable(datos, caption = "Los datos", options = list(pageLength = 10))
## Warning in instance$preRenderHook(instance): It seems your data is too big
## for client-side DataTables. You may consider server-side processing: https://
## rstudio.github.io/DT/server.html
70 % datos de entrenamiento 30 % datos de validación
set.seed(2020)
entrena <- createDataPartition(y = datos$income10, p = 0.7, list = FALSE, times = 1)
# Datos entrenamiento
datos.entrenamiento <- datos[entrena, ] # [renglones, columna]
# Datos validación
datos.validacion <- datos[-entrena, ]
# kable(head(datos.entrenamiento, 10), caption = "Datos de entrenamiento (primeros diez)", row.names = 1:nrow(datos.entrenamiento))
# kable(head(datos.validacion, 10), caption = "Datos de validación (primeros diez)", row.names = 1:nrow(datos.entrenamiento))
datatable(datos.entrenamiento, caption = "Datos de entrenamiento", options = list(pageLength = 10))
## Warning in instance$preRenderHook(instance): It seems your data is too big
## for client-side DataTables. You may consider server-side processing: https://
## rstudio.github.io/DT/server.html
datatable(datos.validacion, caption = "Datos de entrenamiento",options = list(pageLength = 10))
## Warning in instance$preRenderHook(instance): It seems your data is too big
## for client-side DataTables. You may consider server-side processing: https://
## rstudio.github.io/DT/server.html
Con la regresión logística, dado un conjunto particular de valores de las variables independientes elegidas, se estima la probabilidad de los ingresos de una persona ‘<=50’ o ‘>50’ o lo que es lo mismo ingresos con valores de 0 y 1 en la variable income10.
Por medio de la función gml() se contruye un modelo de regresión logística
Variable dependiente o predictiva es ‘income10’, ya que depende de todas las demás variables
Variables independientes o predictoras, todas las demás: “age”, “workclass”, “education”, “educational.num”, “marital.status”, “race”, “gender”, “hours.per.week”
Se utiliza el conjunto de datos de entrenamiento
La finalidad de consruir el modelo de regresión logística es entre otros cosas, para conocer los coeficienes y el nivel de significación de cada variable independiente o predictora así como las pruebas t y F
La ecuación income10=age.scale+workclass+education+marital.status+race+gender+hours.per.week.scale
Se asigna a una variable formula y se utiliza para construir el modelo,
Significa que la variable ‘income10’ depende o es dependiente de todas las demás variables
La fórmula dada por la ecuacuación:
formula = income10 ~ age.scale + workclass + education + marital.status + race + gender + hours.per.week.scale
modelo <- glm(formula, data = datos.entrenamiento, family = 'binomial')
summary(modelo)
##
## Call:
## glm(formula = formula, family = "binomial", data = datos.entrenamiento)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -2.7548 -0.5778 -0.2603 -0.0684 3.3249
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -2.76730 0.23517 -11.767 < 2e-16 ***
## age.scale 2.04720 0.10474 19.546 < 2e-16 ***
## workclassFederal-gov 1.39346 0.12545 11.107 < 2e-16 ***
## workclassLocal-gov 0.76117 0.11093 6.862 6.79e-12 ***
## workclassNever-worked -7.98118 98.02680 -0.081 0.93511
## workclassPrivate 0.88906 0.09705 9.161 < 2e-16 ***
## workclassSelf-emp-inc 1.32433 0.11899 11.130 < 2e-16 ***
## workclassSelf-emp-not-inc 0.31448 0.10762 2.922 0.00348 **
## workclassState-gov 0.48978 0.12316 3.977 6.99e-05 ***
## workclassWithout-pay 0.15627 0.85256 0.183 0.85457
## educationCommunity -1.00153 0.04434 -22.589 < 2e-16 ***
## educationDropout -2.69532 0.07543 -35.731 < 2e-16 ***
## educationHighGrad -1.58616 0.04511 -35.164 < 2e-16 ***
## educationMaster 0.59823 0.06103 9.803 < 2e-16 ***
## educationPhD 1.22350 0.13962 8.763 < 2e-16 ***
## marital.statusNot_married -2.53858 0.05407 -46.948 < 2e-16 ***
## marital.statusSeparated -2.10873 0.05619 -37.531 < 2e-16 ***
## marital.statusWidow -2.10698 0.12348 -17.064 < 2e-16 ***
## raceAsian-Pac-Islander 0.32012 0.21854 1.465 0.14298
## raceBlack 0.34440 0.20941 1.645 0.10005
## raceOther 0.26768 0.28964 0.924 0.35540
## raceWhite 0.58693 0.20059 2.926 0.00343 **
## genderMale 0.10151 0.04450 2.281 0.02252 *
## hours.per.week.scale 3.05128 0.13834 22.056 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 37584 on 34189 degrees of freedom
## Residual deviance: 25090 on 34166 degrees of freedom
## AIC: 25138
##
## Number of Fisher Scoring iterations: 11
Para evaluar el rendimiento del modelo, se crea la matriz de confusión
Una matriz de confusión es una herramienta que permite la visualización del desempeño de un algoritmo que se emplea en aprendizaje supervisado.
Cada columna de la matriz representa el número de predicciones de cada clase, mientras que cada fila representa a las instancias en la clase real.
Uno de los beneficios de las matrices de confusión es que facilitan ver si el sistema está confundiendo las diferentes clases o resultados.
Matriz de Confusión
include_graphics(“../imagenes/matriz confusion.jpg”)
include_graphics('C:/Users/Blue/Documents/MatrizC.jpg')
Comparar los valores de income10 VS valores ajustados Utilizando los datos de entrenamiento
Tres variables income10 original con valores 0 y 1 s
valores ajustados valores ajustados codificados 0 y 1 s aquellos cuya probabilidad sea > 0.5 o al 50% Con las columnas1 y 3 se puede generar la matriz de confusión ¿Qué son los valores ajustados?, es la probabilidad, si es mayor al 50% entonces es 1, si es menor o igual al 50% entoces es 0 fitted.values
Se observan los primeros diez y últimos diez registros
comparar <- data.frame(datos.entrenamiento$income10, as.vector(modelo$fitted.values) )
comparar <- comparar %>%
mutate(income10ajustados = if_else (modelo$fitted.values > 0.5, 1, 0))
colnames(comparar) <- c("income10", "ajuste", 'income10ajustados')
# kable(head(comparar, 10), caption = "Comparar valores, primeros diez")
# kable(tail(comparar, 10), caption = "Comparar valores, útlimos diez")
datatable(comparar, caption = "Comparar valores", options = list(pageLength = 10))
## Warning in instance$preRenderHook(instance): It seems your data is too big
## for client-side DataTables. You may consider server-side processing: https://
## rstudio.github.io/DT/server.html
### Visualizando los valores ajustados
ggplot(data = comparar, aes(x =row.names(comparar), y = ajuste, )) +
geom_point(aes(colour = factor(income10))) +
geom_hline(yintercept = 0.50)
Matriz de confusión Con los datos generados en la variable comparar se genera la matriz de confusión
matriz_confusion <- table(comparar$income10, comparar$income10ajustados, dnn = c("income10", "income10ajustados para predicciones"))
kable(matriz_confusion, caption = "Matriz de confusión")
| 0 | 1 | |
|---|---|---|
| 0 | 24135 | 1892 |
| 1 | 4056 | 4107 |
¿Qué significa la matriz de confusión?
El modelo dice que 24135 ganan por debajo de 50 Mil de los 34190 registros totales. De tal forma que le atinó a 70.59 a los VP verdaderos positivos
El modelo dice que 4107 ganan por encima de 50 Mil de los 34190 registros totales. De tal forma que le atinó a 12.01 a los FP verdaderos negativos
De acuerdo a estas condiciones:
include_graphics('C:/Users/Blue/Documents/M.jpg')
Entonces en cuanto a exactitud (accuracy) como una forma de evaluar el modelo en los datos de entrenamiento se tiene un: 82.6, “%” que significa:
El modelo es capaz de predecir y clasificar con exactidud al 82.6, “%”, o sea que se puede equivocar en 17.4% de los casos
predicciones <- predict(object = modelo, newdata = datos.validacion, se.fit = TRUE)
predicciones.ajustadas <- predicciones$fit
kable(head(predicciones.ajustadas, 10))
| x | |
|---|---|
| 3 | -0.7964486 |
| 4 | -0.4623941 |
| 6 | -5.0440275 |
| 12 | 1.0617147 |
| 19 | -3.8320026 |
| 25 | 0.2488330 |
| 26 | 0.2057424 |
| 29 | -0.6797358 |
| 41 | 2.2801743 |
| 42 | 0.3710815 |
include_graphics('C:/Users/Blue/Documents/F.jpg')
predicciones_prob <- exp(predicciones$fit) / (1 + exp(predicciones$fit))
kable(head(predicciones_prob))
| x | |
|---|---|
| 3 | 0.3107857 |
| 4 | 0.3864180 |
| 6 | 0.0064064 |
| 12 | 0.7430181 |
| 19 | 0.0212067 |
| 25 | 0.5618892 |
las.predicciones <- cbind(datos.validacion, predicciones_prob)
las.predicciones <- las.predicciones %>%
mutate(income10.prediccion = if_else(predicciones_prob > 0.5, 1, 0))
kable(head(las.predicciones))
| X | age | workclass | education | educational.num | marital.status | race | gender | hours.per.week | income | age.scale | educational.num.scale | hours.per.week.scale | income10 | predicciones_prob | income10.prediccion |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| 3 | 28 | Local-gov | Community | 12 | Married | White | Male | 40 | >50K | 0.1506849 | 0.7333333 | 0.3979592 | 1 | 0.3107857 | 0 |
| 4 | 44 | Private | Community | 10 | Married | Black | Male | 40 | >50K | 0.3698630 | 0.6000000 | 0.3979592 | 1 | 0.3864180 | 0 |
| 6 | 34 | Private | Dropout | 6 | Not_married | White | Male | 30 | <=50K | 0.2328767 | 0.3333333 | 0.2959184 | 0 | 0.0064064 | 0 |
| 12 | 36 | Federal-gov | Bachelors | 13 | Married | White | Male | 40 | <=50K | 0.2602740 | 0.8000000 | 0.3979592 | 0 | 0.7430181 | 1 |
| 19 | 37 | Private | HighGrad | 9 | Widow | White | Female | 20 | <=50K | 0.2739726 | 0.5333333 | 0.1938776 | 0 | 0.0212067 | 0 |
| 25 | 25 | Private | Bachelors | 13 | Married | White | Male | 40 | <=50K | 0.1095890 | 0.8000000 | 0.3979592 | 0 | 0.5618892 | 1 |
matriz_confusion <- table(las.predicciones$income10, las.predicciones$income10.prediccion, dnn = c("income10", "predicciones"))
matriz_confusion
## predicciones
## income10 0 1
## 0 10361 767
## 1 1752 1772
matriz_confusion <- table(las.predicciones$income10, las.predicciones$income10.prediccion, dnn = c("income10", "predicciones"))
matriz_confusion
## predicciones
## income10 0 1
## 0 10361 767
## 1 1752 1772
edad <- c(43,54,66)
clase.empleo <- c('Federal-gov', 'State-gov','Never-worked')
nivel.educacion <- c('HighGrad', 'Community', 'Bachelors')
edo.civil <- c('Married', 'Married', 'Not_married')
raza <- c('White', 'White', 'Black')
genero <- c('Female', 'Male', 'Female')
horas <- c(54,49,57)
a.predecir <- data.frame(rbind(cbind(edad, clase.empleo, nivel.educacion, edo.civil,raza, genero, horas)))
colnames(a.predecir) <- c('age.scale', 'workclass', 'education', 'marital.status', 'race', 'gender', 'hours.per.week.scale')
a.predecir
## age.scale workclass education marital.status race gender
## 1 43 Federal-gov HighGrad Married White Female
## 2 54 State-gov Community Married White Male
## 3 66 Never-worked Bachelors Not_married Black Female
## hours.per.week.scale
## 1 54
## 2 49
## 3 57
edad; horas
## [1] 43 54 66
## [1] 54 49 57
edad.escalada <- rescale(c(edad, min(datos$age), max(datos$age)))
edad.escalada <- edad.escalada[1:3]
# Escalando las horas por semana
horas.escalada<- rescale(c(horas, min(datos$hours.per.week), max(datos$hours.per.week)))
horas.escalada <- horas.escalada[1:3]
a.predecir <- a.predecir %>%
mutate(age.scale = edad.escalada,
hours.per.week.scale = horas.escalada)
a.predecir
## age.scale workclass education marital.status race gender
## 1 0.3561644 Federal-gov HighGrad Married White Female
## 2 0.5068493 State-gov Community Married White Male
## 3 0.6712329 Never-worked Bachelors Not_married Black Female
## hours.per.week.scale
## 1 0.5408163
## 2 0.4897959
## 3 0.5714286
Realizar la predicción Establecer la probabilidad de predicción Determinar si es 0 a 1 la predicción
prediccion <- predict(modelo, a.predecir, se.fit = TRUE)
prediccion
## $fit
## 1 2 3
## 0.006243832 -0.058486973 -9.824937653
##
## $se.fit
## 1 2 3
## 0.09942861 0.08658633 98.02677753
##
## $residual.scale
## [1] 1
prediccion_prob <- exp(prediccion$fit) / (1 + exp(prediccion$fit))
prediccion_prob
## 1 2 3
## 5.015610e-01 4.853824e-01 5.408294e-05
# Predecir si será 0 o 1 conforme a la probabilidad de predicción
las.predicciones <- if_else(prediccion_prob > 0.5, 1, 0)
cat("Son las predicciones para las personas con esas características ")
## Son las predicciones para las personas con esas características
las.predicciones
## [1] 1 0 0
print("Las predicciones en la columna final")
## [1] "Las predicciones en la columna final"
cbind(a.predecir ,las.predicciones)
## age.scale workclass education marital.status race gender
## 1 0.3561644 Federal-gov HighGrad Married White Female
## 2 0.5068493 State-gov Community Married White Male
## 3 0.6712329 Never-worked Bachelors Not_married Black Female
## hours.per.week.scale las.predicciones
## 1 0.5408163 1
## 2 0.4897959 0
## 3 0.5714286 0