Objetivo

Realizar e interpretar regresión logísitica con datos de personas e ingresos de USA

Descripción

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.

Proceso

1. Cargar librerías

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

3. Las variables

4. Crear datos de entrenamiento y validación

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

5. Crear modelo de regresión logística

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')

6. Analizar y/o describir el modelo

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

7. Evaluar el modelo con matriz de confusión

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")
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

8. Realizar predicciones con el conjunto de datos de validación

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

Evaluar el modelo de predicción

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 de confusión de las predicciones

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

Evaluando la matriz de confusión de las predicciones

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

Realizar predicciones con al menos tres nuevos registros

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

Segundo: Escalar la edad y las horas trabajadas

  • Se escalan la edad y las horas trabajadas
  • Los primeros 3 valores escalados son los que interesan porque son 3 nuevos registros
  • Se escalan de igual forma con los valores originales de datos tanto en age como en hours.per.week
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

Tercero. Realizar la predicción con los nuevos registros

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

Interpretacion

De acuerdo a los datos de ingresos de adultos en el cual tenemos datos como edad de la persona, nivel educativo, estado civil, raza, horas por semana trabajadas, entre otras, se realizo una regresion logistica en el cual 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.

De acuerdo ala matriz de confucion nos dice que 24135 personas ganan por debajo de $50 mil de los 34,190 registros totales, se atino a 70.59 a los VP verdaderos positivos ademas de los registros totales anteriores nos dice el modelo que 4107 ganan por encima de $50 mil, se atino a 12.01 a los FP verdaderos negativos.

Por ultimo se realizo la prediccion de 3 nuevas observaciones aleatorias en las cuales se analizaron clase de empleo, nivel de edcuacion, estado civil , raza, genero y horas trabajadas. El resultado obtenido de las observaciones es que dos personas ganan menos de $50 mil y una gana mas de $50 mil de acuerdo a las variables analizadas.