Health Insurence

Preguntas Guía

i) ¿Qué es Supervised Machine Learning y cuáles son algunas de sus aplicaciones en Inteligencia de Negocios?

Supervised learning, also known as supervised machine learning, is a subcategory of machine learning and artificial intelligence. It is defined by its use of labeled datasets to train algorithms that to classify data or predict outcomes accurately. As input data is fed into the model, it adjusts its weights until the model has been fitted appropriately, which occurs as part of the cross validation process. Supervised learning helps organizations solve for a variety of real-world problems at scale, such as classifying spam in a separate folder from your inbox.

ii) ¿Cuáles son los principales algoritmos de Supervised Machine Learning?

En los procesos de machine learning supervisado se utilizan varios algoritmos y técnicas de cálculo. A continuación, se incluyen breves explicaciones de algunos de los métodos de aprendizaje más utilizados, generalmente calculados mediante el uso de programas como R o Python:

1. Redes Neuronales:

Las redes neuronales son modelos utilizados en el ámbito del aprendizaje profundo para procesar datos imitando la interconexión de neuronas en el cerebro humano. Estas redes consisten en capas de nodos, donde cada nodo recibe entradas, aplica ponderaciones, un sesgo y produce una salida. A través del aprendizaje supervisado, estas redes ajustan sus conexiones basadas en una función de pérdida, buscando minimizar errores y mejorar la precisión del modelo.

2. Naive Bayes:

Naive Bayes es un método de clasificación que se basa en el principio de independencia condicional de clases, según el Teorema de Bayes. Este enfoque asume que las características predictoras son independientes entre sí, lo que lo hace eficaz en la clasificación de textos, sistemas de recomendación y detección de correo no deseado. Se divide en tres tipos: Multinomial, Bernoulli y Gaussiano.

3. Regresión Lineal:

La regresión lineal busca entender la relación entre una variable dependiente y una o más variables independientes para realizar predicciones. Se ajusta a través del método de mínimos cuadrados, intentando trazar una línea recta que mejor se ajuste a los datos. Si hay una variable independiente, se llama regresión lineal simple; si hay múltiples variables independientes, es una regresión lineal múltiple.

4. Regresión Logística:

La regresión logística se emplea cuando la variable dependiente es categórica, es decir, tiene resultados binarios como “verdadero” o “falso”. A diferencia de la regresión lineal, que se utiliza con variables continuas, la regresión logística se aplica principalmente en problemas de clasificación binaria, como la detección de correo no deseado.

5. Máquinas de Vectores de Soporte (SVM):

Las SVM son modelos de aprendizaje supervisado utilizados para clasificación y regresión. Buscan un hiperplano que maximice la distancia entre dos clases de puntos de datos, actuando como límite de decisión. Esto permite separar eficientemente las clases de puntos de datos en el espacio.

6. K Vecino Más Cercano (KNN):

El algoritmo KNN clasifica los puntos de datos según la proximidad a otros puntos de datos conocidos. Se basa en la suposición de que puntos similares están cerca unos de otros. Utiliza la distancia euclidiana para calcular la similitud y asignar una categoría basada en la mayoría de las categorías de los vecinos más cercanos.

7. Bosque Aleatorio:

El bosque aleatorio es un algoritmo flexible de aprendizaje supervisado que combina múltiples árboles de decisión para mejorar la precisión de las predicciones. Cada árbol de decisión opera de manera independiente y su resultado se combina para reducir la varianza y mejorar la capacidad predictiva del modelo.

iii) ¿Qué es la R2 Ajustada? ¿Qué es la métrica RMSE? ¿Cuál es la diferencia entre la R2 Ajustada y la métrica RMSE?

** El R cuadrado ajustado (o coeficiente de determinación ajustado) se utiliza en la regresión múltiple para ver el grado de intensidad o efectividad que tienen las variables independientes en explicar la variable dependiente.

** El error cuadrático medio (RMSE) es una regla de puntuación cuadrática que también mide la magnitud media del error. Es la raíz cuadrada del promedio de diferencias cuadradas entre la predicción y la observación real.

** En resumen, la R2 ajustada evalúa la capacidad explicativa del modelo considerando el número de predictores, mientras que el RMSE evalúa la precisión de las predicciones del modelo sin considerar la complejidad del mismo. Ambas métricas son importantes en la evaluación de modelos de regresión, ya que proporcionan información sobre diferentes aspectos de su desempeño.

Referencias Bibliograficas

** What is Supervised Learning?  | IBM. (s. f.). https://www.ibm.com/topics/supervised-learning

** MAE y RMSE: ¿qué métrica es mejor? (2020, 27 noviembre). ICHI.PRO. https://ichi.pro/es/mae-y-rmse-que-metrica-es-mejor-252933908062525

** Sanjuán, F. J. M. (2022, 24 noviembre). R cuadrado ajustado (Coeficiente de determinación ajustado). Economipedia. https://economipedia.com/definiciones/r-cuadrado-ajustado-coeficiente-de-determinacion-ajustado.html

Librerías

library(psych)
library(ggplot2)
library(dplyr)
library(GGally)
library(ggcorrplot)
library(dplyr)
library(car)
library(lmtest)
library(regclass)       
library(mctest)         
library(lmtest)        
library(spdep)          
library(sf)            
library(spData)         
library(spatialreg)     
library(caret)         
library(e1071)          
library(SparseM)       
library(Metrics)       
library(randomForest)  
library(jtools)         
library(xgboost)        
library(DiagrammeR)     
library(effects)        
library(rpart.plot)    
library(neuralnet)      
library(MASS)           
library(sp)

Análisis Exploratorio de los Datos (EDA)

df=read.csv("C:\\Users\\lesda\\OneDrive\\Documentos\\Concentracion IA\\Estadistica\\health_insurance ACT 1.csv")
summary(df)

##       age            sex                 bmi           children    
##  Min.   :18.00   Length:1338        Min.   :16.00   Min.   :0.000  
##  1st Qu.:27.00   Class :character   1st Qu.:26.30   1st Qu.:0.000  
##  Median :39.00   Mode  :character   Median :30.40   Median :1.000  
##  Mean   :39.21                      Mean   :30.67   Mean   :1.095  
##  3rd Qu.:51.00                      3rd Qu.:34.70   3rd Qu.:2.000  
##  Max.   :64.00                      Max.   :53.10   Max.   :5.000  
##     smoker             region             expenses    
##  Length:1338        Length:1338        Min.   : 1122  
##  Class :character   Class :character   1st Qu.: 4740  
##  Mode  :character   Mode  :character   Median : 9382  
##                                        Mean   :13270  
##                                        3rd Qu.:16640  
##                                        Max.   :63770

str(df)

## 'data.frame':    1338 obs. of  7 variables:
##  $ age     : int  19 18 28 33 32 31 46 37 37 60 ...
##  $ sex     : chr  "female" "male" "male" "male" ...
##  $ bmi     : num  27.9 33.8 33 22.7 28.9 25.7 33.4 27.7 29.8 25.8 ...
##  $ children: int  0 1 3 0 0 0 1 3 2 0 ...
##  $ smoker  : chr  "yes" "no" "no" "no" ...
##  $ region  : chr  "southwest" "southeast" "southeast" "northwest" ...
##  $ expenses: num  16885 1726 4449 21984 3867 ...

Medidas de tendencia central

describe(df)

##          vars    n     mean       sd  median  trimmed     mad     min      max
## age         1 1338    39.21    14.05   39.00    39.01   17.79   18.00    64.00
## sex*        2 1338     1.51     0.50    2.00     1.51    0.00    1.00     2.00
## bmi         3 1338    30.67     6.10   30.40    30.50    6.23   16.00    53.10
## children    4 1338     1.09     1.21    1.00     0.94    1.48    0.00     5.00
## smoker*     5 1338     1.20     0.40    1.00     1.13    0.00    1.00     2.00
## region*     6 1338     2.52     1.10    3.00     2.52    1.48    1.00     4.00
## expenses    7 1338 13270.42 12110.01 9382.03 11076.02 7440.81 1121.87 63770.43
##             range  skew kurtosis     se
## age         46.00  0.06    -1.25   0.38
## sex*         1.00 -0.02    -2.00   0.01
## bmi         37.10  0.28    -0.06   0.17
## children     5.00  0.94     0.19   0.03
## smoker*      1.00  1.46     0.14   0.01
## region*      3.00 -0.04    -1.33   0.03
## expenses 62648.56  1.51     1.59 331.07

moda <- function(x) {
  tabla_frecuencias <- table(x)
  modas <- names(tabla_frecuencias[tabla_frecuencias == max(tabla_frecuencias)])
  if (is.numeric(x)) {
    modas <- as.numeric(modas)
  }
  return(modas)
}
# Usar sapply para aplicar la función moda a cada columna
modas_por_columna <- sapply(df, moda)

# Mostrar las modas de cada columna
modas_por_columna

## $age
## [1] 18
## 
## $sex
## [1] "male"
## 
## $bmi
## [1] 27.6 33.3
## 
## $children
## [1] 0
## 
## $smoker
## [1] "no"
## 
## $region
## [1] "southeast"
## 
## $expenses
## [1] 1639.56

Medidas de dispersión

# Función para calcular estadísticas
calcular_estadisticas <- function(columna) {
  rango <- max(columna, na.rm = TRUE) - min(columna, na.rm = TRUE)
  varianza <- var(columna, na.rm = TRUE)
  desviacion_std <- sd(columna, na.rm = TRUE)
  rango_iqr <- IQR(columna, na.rm = TRUE)
  
  c(Rango = rango, Varianza = varianza, `Desviación Estándar` = desviacion_std, `Rango Intercuartílico` = rango_iqr)
}

# Aplicar la función a cada columna numérica de df
estadisticas_df <- sapply(df, function(x) if(is.numeric(x)) calcular_estadisticas(x) else NA)

# Mostrar las estadísticas
estadisticas_df

## $age
##                 Rango              Varianza   Desviación Estándar 
##              46.00000             197.40139              14.04996 
## Rango Intercuartílico 
##              24.00000 
## 
## $sex
## [1] NA
## 
## $bmi
##                 Rango              Varianza   Desviación Estándar 
##             37.100000             37.190265              6.098382 
## Rango Intercuartílico 
##              8.400000 
## 
## $children
##                 Rango              Varianza   Desviación Estándar 
##              5.000000              1.453213              1.205493 
## Rango Intercuartílico 
##              2.000000 
## 
## $smoker
## [1] NA
## 
## $region
## [1] NA
## 
## $expenses
##                 Rango              Varianza   Desviación Estándar 
##              62648.56          146652372.23              12110.01 
## Rango Intercuartílico 
##              11899.63

Identificación de patrones, tendencias

#Distribución de la edad

ggplot(df, aes(x=age)) + geom_histogram(binwidth=5, fill="pink", color="black") + labs(title="Distribución de Edad", x="Edad", y="Frecuencia")

#Distribución por sexo

ggplot(df, aes(x=sex)) + geom_bar(fill="magenta", color="black") + labs(title="Distribución por Sexo", x="Sexo", y="Conteo")

#Comparación entre BMI y edad por sexo
ggplot(df, aes(x=age, y=bmi)) + geom_point(aes(color=sex)) + labs(title="BMI vs Edad por Sexo", x="Edad", y="BMI") + theme_minimal()

#BMI por sexo

ggplot(df, aes(x=sex, y=bmi, fill=sex)) + geom_boxplot() + labs(title="Distribución de BMI por Sexo", x="Sexo", y="BMI")

#Comparación de fumadores y no fumadores
ggplot(df, aes(x=smoker)) + geom_bar(aes(fill=smoker)) + labs(title="Distribución de Fumadores", x="Fumador", y="Conteo")

#Cantidad de hijos
ggplot(df, aes(x=factor(children))) + geom_bar(fill="magenta", color="black") + labs(title="Distribución del Número de Hijos", x="Número de Hijos", y="Conteo")

#Distribución de gastos por región

 ggplot(df, aes(x=region, y=expenses, fill=region)) + geom_violin() + labs(title="Distribución de Gastos por Región", x="Región", y="Gastos")

#Relación entre BMI y número de hijos
ggplot(df, aes(x=factor(children), y=bmi)) + geom_jitter(aes(color=sex), width=0.2) + labs(title="BMI por Número de Hijos y Sexo", x="Número de Hijos", y="BMI")

Distribuciones de cada variable

  hist(df$age, 
     main = "Histograma de Edad", 
     xlab = "Edad", 
     ylab = "Frecuencia", 
     col = "orange", 
     border = "black")

hist(df$bmi, 
     main = "Histograma de Índice de masa corporal", 
     xlab = "Edad", 
     ylab = "Frecuencia", 
     col = "yellow", 
     border = "black")

hist(df$children, 
     main = "Histograma de número de hijos", 
     xlab = "Edad", 
     ylab = "Frecuencia", 
     col = "lightblue", 
     border = "black")

hist(df$expenses, 
     main = "Histograma de prima", 
     xlab = "Prima", 
     ylab = "Frecuencia", 
     col = "purple", 
     border = "black")

Relaciones entre variables independientes y la dependiente

# Relaciones entre variables independientes y la dependiente

ggplot(df, aes(x=age, y=expenses, color=sex)) +
  geom_point(alpha=0.5) +
  labs(title="Gastos vs. Edad, Diferenciados por Sexo", x="Edad", y="Gastos Médicos") +
  theme_minimal()

ggplot(df, aes(x=sex, y=expenses, color=sex)) +
  geom_point(position = position_jitterdodge(), alpha=0.5) +
  labs(title="Gastos Médicos Diferenciados por Sexo", x="Sexo", y="Gastos Médicos") +
  theme_minimal()

ggplot(df, aes(x=bmi, y=expenses, color=smoker)) +
  geom_point(alpha=0.5) +
  labs(title="Gastos vs. BMI, Diferenciados por Fumadores", x="BMI", y="Gastos Médicos") +
  theme_minimal()

ggplot(df, aes(x=factor(children), y=expenses)) +
  geom_point(position = position_jitter(width = 0.2), aes(color=factor(children)), alpha=0.5) +
  labs(title="Gastos vs. Número de Hijos", x="Número de Hijos", y="Gastos Médicos") +
  theme_minimal()

ggplot(df, aes(x=region, y=expenses, color=region)) +
  geom_point(position = position_jitter(width = 0.2), alpha=0.5) +
  labs(title="Gastos Médicos Diferenciados por Región", x="Región", y="Gastos Médicos") +
  theme_minimal()

ggplot(df, aes(x=region, y=expenses, fill=region)) + 
  geom_boxplot() +
  scale_fill_brewer(palette="Pastel1") +
  labs(title="Gastos por Región", x="Región", y="Gastos Médicos") +
  theme_minimal()

ggplot(df, aes(x=sex, y=expenses, fill=sex)) + 
  geom_boxplot() +
  scale_fill_brewer(palette="Set1") +
  labs(title="Gastos por Sexo", x="Sexo", y="Gastos Médicos") +
  theme_minimal()

ggplot(df, aes(x=smoker, y=expenses, fill=smoker)) + 
  geom_boxplot() +
  scale_fill_brewer(palette="Set2") +
  labs(title="Gastos por Fumador", x="Fumador", y="Gastos Médicos") +
  theme_minimal()

ggplot(df, aes(x=as.factor(children), y=expenses)) +
  geom_boxplot(aes(fill=as.factor(children))) +
  scale_fill_brewer(palette="Pastel1") + 
  labs(title="Gastos Médicos por Número de Hijos",
       x="Número de Hijos",
       y="Gastos Médicos") +
  theme_minimal()

Correlación entre variables

# Correlación entre variables
# Correlation
model.matrix(~0+., data=df) %>% 
  cor(use="pairwise.complete.obs") %>% 
  ggcorrplot(show.diag = F, type="lower", lab=TRUE, lab_size=2, 
             colors = c("pink", "white", "purple")) 

age : int 19 18 28 33 32 31 46 37 37 60 … $ sex : chr “female” “male” “male” “male” … $ bmi : num 27.9 33.8 33 22.7 28.9 25.7 33.4 27.7 29.8 25.8 … $ children: int 0 1 3 0 0 0 1 3 2 0 … $ smoker : chr “yes” “no” “no” “no” … $ region : chr “southwest” “southeast” “southeast” “northwest” … $ expenses: num 16885 1726 4449 21984 3867 …

Modelo de regresión lineal

ols_model <- lm(expenses ~ age + bmi + smoker + region, data = df)
summary(ols_model)

## 
## Call:
## lm(formula = expenses ~ age + bmi + smoker + region, data = df)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -11902.8  -3038.0   -997.5   1528.2  29408.0 
## 
## Coefficients:
##                  Estimate Std. Error t value Pr(>|t|)    
## (Intercept)     -11604.13     976.19 -11.887   <2e-16 ***
## age                258.62      11.93  21.679   <2e-16 ***
## bmi                340.09      28.67  11.862   <2e-16 ***
## smokeryes        23851.43     413.50  57.682   <2e-16 ***
## regionnorthwest   -303.33     477.84  -0.635   0.5257    
## regionsoutheast  -1039.16     480.48  -2.163   0.0307 *  
## regionsouthwest   -915.16     479.54  -1.908   0.0566 .  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 6085 on 1331 degrees of freedom
## Multiple R-squared:  0.7487, Adjusted R-squared:  0.7475 
## F-statistic: 660.8 on 6 and 1331 DF,  p-value: < 2.2e-16

# Para 'age'
df <- df[!(df$age < (quantile(df$age, 0.25) - 1.5*IQR(df$age)) | df$age > (quantile(df$age, 0.75) + 1.5*IQR(df$age))),]

# Para 'bmi'
df <- df[!(df$bmi < (quantile(df$bmi, 0.25) - 1.5*IQR(df$bmi)) | df$bmi > (quantile(df$bmi, 0.75) + 1.5*IQR(df$bmi))),]

df <- df[!(df$expenses < (quantile(df$expenses, 0.25) - 1.5 * IQR(df$expenses)) | df$expenses > (quantile(df$expenses, 0.75) + 1.5 * IQR(df$expenses))), ]


log_ols_model <- lm(log(expenses) ~ log(age) + log(bmi) + as.factor(smoker) + children  , data = df)
summary(log_ols_model)

## 
## Call:
## lm(formula = log(expenses) ~ log(age) + log(bmi) + as.factor(smoker) + 
##     children, data = df)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.73765 -0.18524 -0.08124  0.02947  2.24392 
## 
## Coefficients:
##                      Estimate Std. Error t value Pr(>|t|)    
## (Intercept)           3.40936    0.24249  14.060  < 2e-16 ***
## log(age)              1.34353    0.03264  41.163  < 2e-16 ***
## log(bmi)              0.13081    0.06617   1.977   0.0483 *  
## as.factor(smoker)yes  1.30216    0.04115  31.646  < 2e-16 ***
## children              0.08576    0.01043   8.220  5.3e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.4361 on 1186 degrees of freedom
## Multiple R-squared:  0.6993, Adjusted R-squared:  0.6983 
## F-statistic: 689.7 on 4 and 1186 DF,  p-value: < 2.2e-16

Evaluación del modelo

AIC(ols_model)  

## [1] 27123.48

AIC(log_ols_model) 

## [1] 1410.047

RMSE_ols_model     <- sqrt(mean(ols_model$residuals^2))
RMSE_log_ols_model <- sqrt(mean(log_ols_model$residuals^2))

RMSE_ols_model

## [1] 6068.772

RMSE_log_ols_model

## [1] 0.4351707

Pruebas de diagnóstico

Multicolinealidad

vif(log_ols_model) 

##          log(age)          log(bmi) as.factor(smoker)          children 
##          1.023662          1.084219          1.072588          1.007884

Heterocedasticidad

bptest(log_ols_model)

## 
##  studentized Breusch-Pagan test
## 
## data:  log_ols_model
## BP = 79.894, df = 4, p-value < 2.2e-16

Autocorrelación serial

NA

Autocorrelación espacial

NA

Normalidad de los residuales

shapiro.test(resid(log_ols_model))

## 
##  Shapiro-Wilk normality test
## 
## data:  resid(log_ols_model)
## W = 0.7389, p-value < 2.2e-16

#Gráfico de residuos
hist(log_ols_model$residuals,
col = "purple", 
border = "black")

Transformaciones a las variables

#Quitar outliers

df2=df
df2 <- df2[!(df$expenses < (quantile(df$expenses, 0.25) - 1.5 * IQR(df$expenses)) | df$expenses > (quantile(df$expenses, 0.75) + 1.5 * IQR(df$expenses))), ]
# Para 'age'
df2 <- df2[!(df$age < (quantile(df$age, 0.25) - 1.5*IQR(df$age)) | df$age > (quantile(df$age, 0.75) + 1.5*IQR(df$age))),]

# Para 'bmi'
df2 <- df2[!(df$bmi < (quantile(df$bmi, 0.25) - 1.5*IQR(df$bmi)) | df$bmi > (quantile(df$bmi, 0.75) + 1.5*IQR(df$bmi))),]

log_ols_model2 <- lm(log(expenses) ~ log(age) + log(bmi) + as.factor(smoker) + children  , data = df2)
summary(log_ols_model)

## 
## Call:
## lm(formula = log(expenses) ~ log(age) + log(bmi) + as.factor(smoker) + 
##     children, data = df)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.73765 -0.18524 -0.08124  0.02947  2.24392 
## 
## Coefficients:
##                      Estimate Std. Error t value Pr(>|t|)    
## (Intercept)           3.40936    0.24249  14.060  < 2e-16 ***
## log(age)              1.34353    0.03264  41.163  < 2e-16 ***
## log(bmi)              0.13081    0.06617   1.977   0.0483 *  
## as.factor(smoker)yes  1.30216    0.04115  31.646  < 2e-16 ***
## children              0.08576    0.01043   8.220  5.3e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.4361 on 1186 degrees of freedom
## Multiple R-squared:  0.6993, Adjusted R-squared:  0.6983 
## F-statistic: 689.7 on 4 and 1186 DF,  p-value: < 2.2e-16

Heterocedasticidad

bptest(log_ols_model2)

## 
##  studentized Breusch-Pagan test
## 
## data:  log_ols_model2
## BP = 75.24, df = 4, p-value = 1.773e-15

Normalidad de los residuales

shapiro.test(resid(log_ols_model2))

## 
##  Shapiro-Wilk normality test
## 
## data:  resid(log_ols_model2)
## W = 0.69358, p-value < 2.2e-16

Weight OLS

weights <- 1 / log(df2$expenses)

# Ajustar el modelo WLS
log_ols_wls_model <- lm(log(expenses) ~ log(age) + log(bmi) + as.factor(smoker) + children, data = df2, weights = weights)

# Ver el resumen del modelo
summary(log_ols_wls_model)

## 
## Call:
## lm(formula = log(expenses) ~ log(age) + log(bmi) + as.factor(smoker) + 
##     children, data = df2, weights = weights)
## 
## Weighted Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.22687 -0.04532 -0.01595  0.01621  0.74026 
## 
## Coefficients:
##                      Estimate Std. Error t value Pr(>|t|)    
## (Intercept)          3.351724   0.212186  15.796   <2e-16 ***
## log(age)             1.422077   0.029058  48.940   <2e-16 ***
## log(bmi)             0.049728   0.058220   0.854    0.393    
## as.factor(smoker)yes 1.279863   0.041505  30.837   <2e-16 ***
## children             0.091877   0.009305   9.874   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.1274 on 1126 degrees of freedom
##   (53 observations deleted due to missingness)
## Multiple R-squared:  0.7592, Adjusted R-squared:  0.7583 
## F-statistic: 887.4 on 4 and 1126 DF,  p-value: < 2.2e-16

Heterocedasticidad

bptest(log_ols_wls_model)

## 
##  studentized Breusch-Pagan test
## 
## data:  log_ols_wls_model
## BP = 11.771, df = 4, p-value = 0.01914

Normalidad de los residuales

shapiro.test(resid(log_ols_wls_model))

## 
##  Shapiro-Wilk normality test
## 
## data:  resid(log_ols_wls_model)
## W = 0.67331, p-value < 2.2e-16

hist(log_ols_wls_model$residuals,
     col = "pink", 
     border = "black")

Las transformaciones de las variables con la eliminación de los Outliers y el modelo Weight OLS lograron reducir significativamente la Heterocedasticidad, sin embargo se mantiene el problema de la anormalidad de los residuales.

#RMSE

RMSE_log_ols_wls_model <- sqrt(mean(log_ols_wls_model$residuals^2))
RMSE_log_ols_wls_model

## [1] 0.390953

Modelos de machine learning

Cross Validation

df2 <- na.omit(df2)  # Limpia el dataframe de NA
df2_alt <- df2 %>% dplyr::select(expenses, age,bmi,smoker,children)

df2$expenses= log(df2$expenses )
df2$age= log(df2$age )
df2$bmi= log(df2$bmi )

set.seed(123)
partition <- createDataPartition(y = df2$expenses, p = 0.7, list = FALSE)
train = df2[partition, ]
test = df2[-partition, ]
train2 = df2[partition, ]
test2 = df2[-partition, ]

XGBoost Regresión

# **NOTE**: The estimation regression method XGBoost are sensitive to specifying commands such as log() and I()^2 in the regression equation so we will do the data transformation and directly include it the dataset. 





# define explanatory variables (X's) and dependent variable (Y) in training set
train_x = data.matrix(train[, -7])
train_y = train[,7]

# define explanatory variables (X's) and dependent variable (Y) in testing set
test_x = data.matrix(train[, -7])
test_y = train[, 7]

# define final training and testing sets
xgb_train = xgb.DMatrix(data = train_x, label = train_y)
xgb_test  = xgb.DMatrix(data = test_x, label = test_y)

# Lets fit XGBoost regression model and display RMSE for both training and testing data at each round
watchlist = list(train=xgb_train, test=xgb_test)
model_xgb = xgb.train(data=xgb_train, max.depth=3, watchlist=watchlist, nrounds=70) # the more the number of rounds selected, the longer the time to display the results. 

## [1]  train-rmse:5.893701 test-rmse:5.893701 
## [2]  train-rmse:4.146337 test-rmse:4.146337 
## [3]  train-rmse:2.926488 test-rmse:2.926488 
## [4]  train-rmse:2.071960 test-rmse:2.071960 
## [5]  train-rmse:1.479140 test-rmse:1.479140 
## [6]  train-rmse:1.072780 test-rmse:1.072780 
## [7]  train-rmse:0.797338 test-rmse:0.797338 
## [8]  train-rmse:0.616903 test-rmse:0.616903 
## [9]  train-rmse:0.502248 test-rmse:0.502248 
## [10] train-rmse:0.433278 test-rmse:0.433278 
## [11] train-rmse:0.393131 test-rmse:0.393131 
## [12] train-rmse:0.370336 test-rmse:0.370336 
## [13] train-rmse:0.355950 test-rmse:0.355950 
## [14] train-rmse:0.348553 test-rmse:0.348553 
## [15] train-rmse:0.343831 test-rmse:0.343831 
## [16] train-rmse:0.338843 test-rmse:0.338843 
## [17] train-rmse:0.334898 test-rmse:0.334898 
## [18] train-rmse:0.331781 test-rmse:0.331781 
## [19] train-rmse:0.329272 test-rmse:0.329272 
## [20] train-rmse:0.328168 test-rmse:0.328168 
## [21] train-rmse:0.325709 test-rmse:0.325709 
## [22] train-rmse:0.322696 test-rmse:0.322696 
## [23] train-rmse:0.320312 test-rmse:0.320312 
## [24] train-rmse:0.318257 test-rmse:0.318257 
## [25] train-rmse:0.316940 test-rmse:0.316940 
## [26] train-rmse:0.315047 test-rmse:0.315047 
## [27] train-rmse:0.313476 test-rmse:0.313476 
## [28] train-rmse:0.309481 test-rmse:0.309481 
## [29] train-rmse:0.308494 test-rmse:0.308494 
## [30] train-rmse:0.308048 test-rmse:0.308048 
## [31] train-rmse:0.306446 test-rmse:0.306446 
## [32] train-rmse:0.302774 test-rmse:0.302774 
## [33] train-rmse:0.301287 test-rmse:0.301287 
## [34] train-rmse:0.300636 test-rmse:0.300636 
## [35] train-rmse:0.299763 test-rmse:0.299763 
## [36] train-rmse:0.296936 test-rmse:0.296936 
## [37] train-rmse:0.295638 test-rmse:0.295638 
## [38] train-rmse:0.295129 test-rmse:0.295129 
## [39] train-rmse:0.294046 test-rmse:0.294046 
## [40] train-rmse:0.293564 test-rmse:0.293564 
## [41] train-rmse:0.293201 test-rmse:0.293201 
## [42] train-rmse:0.290426 test-rmse:0.290426 
## [43] train-rmse:0.289572 test-rmse:0.289572 
## [44] train-rmse:0.288200 test-rmse:0.288200 
## [45] train-rmse:0.284922 test-rmse:0.284922 
## [46] train-rmse:0.284153 test-rmse:0.284153 
## [47] train-rmse:0.282600 test-rmse:0.282600 
## [48] train-rmse:0.280748 test-rmse:0.280748 
## [49] train-rmse:0.280357 test-rmse:0.280357 
## [50] train-rmse:0.279848 test-rmse:0.279848 
## [51] train-rmse:0.278802 test-rmse:0.278802 
## [52] train-rmse:0.277189 test-rmse:0.277189 
## [53] train-rmse:0.276742 test-rmse:0.276742 
## [54] train-rmse:0.275960 test-rmse:0.275960 
## [55] train-rmse:0.273845 test-rmse:0.273845 
## [56] train-rmse:0.272900 test-rmse:0.272900 
## [57] train-rmse:0.272006 test-rmse:0.272006 
## [58] train-rmse:0.270044 test-rmse:0.270044 
## [59] train-rmse:0.269023 test-rmse:0.269023 
## [60] train-rmse:0.268472 test-rmse:0.268472 
## [61] train-rmse:0.267644 test-rmse:0.267644 
## [62] train-rmse:0.263392 test-rmse:0.263392 
## [63] train-rmse:0.261552 test-rmse:0.261552 
## [64] train-rmse:0.261288 test-rmse:0.261288 
## [65] train-rmse:0.260638 test-rmse:0.260638 
## [66] train-rmse:0.258681 test-rmse:0.258681 
## [67] train-rmse:0.257493 test-rmse:0.257493 
## [68] train-rmse:0.257118 test-rmse:0.257118 
## [69] train-rmse:0.255843 test-rmse:0.255843 
## [70] train-rmse:0.253668 test-rmse:0.253668

# Looks like the lowest RMSE for both training and test dataset is achieved at 59 round. 
# Lets estimate our final regression model
reg_xgb = xgboost(data = xgb_train, max.depth = 3, nrounds = 59, verbose = 0) # setting verbose = 0 avoids to display the training and testing error for each round. 
prediction_xgb_test<-predict(reg_xgb, xgb_test)
RMSE_SVM <- rmse(prediction_xgb_test, test$expenses)

## Warning in actual - predicted: longitud de objeto mayor no es múltiplo de la
## longitud de uno menor

# Lets do some diagnostic check of regression residuals 
xgb_reg_residuals<-test$expenses - prediction_xgb_test

## Warning in test$expenses - prediction_xgb_test: longitud de objeto mayor no es
## múltiplo de la longitud de uno menor

plot(xgb_reg_residuals, xlab= "Dependent Variable", ylab = "Residuals", main = 'XGBoost Regression Residuals')
abline(0,0)

# Plot first 3 trees of model
xgb.plot.tree(model=reg_xgb, trees=0:2)

importance_matrix <- xgb.importance(model = reg_xgb)
xgb.plot.importance(importance_matrix, xlab = "Explanatory Variables X's Importance")

RMSE

# Cálculo del RMSE
RMSE_XGB <- rmse(test$expenses, prediction_xgb_test)

# Imprimir el RMSE
print(RMSE_XGB)

## [1] 1.005894

Árboles de decisión

decision_tree_model <- rpart(expenses ~ age +  bmi + children + smoker, data = train)

# summary(decision_tree_regression)
plot(decision_tree_model, compress = TRUE)
text(decision_tree_model, use.n = TRUE)

rpart.plot(decision_tree_model)

RMSE

### RMSE of DECISION TREE REGRESSION
decision_tree_prediction <- predict(decision_tree_model,test)
RMSE_tree <- rmse(decision_tree_prediction, test$expenses)
RMSE_tree

## [1] 0.376174

colnames(train)

## [1] "age"      "sex"      "bmi"      "children" "smoker"   "region"   "expenses"

Random Forest

rf_model <- randomForest(expenses ~  age + bmi + children + smoker, data= train, proximity=TRUE)
# random_forest<-randomForest(MEDV~.,data=train_alt,importance=TRUE, proximity=TRUE) 
print(rf_model) ### the train data set model accuracy is around 85%.

## 
## Call:
##  randomForest(formula = expenses ~ age + bmi + children + smoker,      data = train, proximity = TRUE) 
##                Type of random forest: regression
##                      Number of trees: 500
## No. of variables tried at each split: 1
## 
##           Mean of squared residuals: 0.1975525
##                     % Var explained: 65.4

# Prediction & Confusion Matrix – test data
rf_prediction <- predict(rf_model,test)

# confusionMatrix(rf_prediction_train_data, train$MEDV) # a confusion matrix is essentially a table that categorizes predictions against actual values.
RMSE_rf <- rmse(rf_prediction, test$expenses)


# Evalute Variables' Importance
# How to interpret varImpPlot()? The higher the value of mean decrease accuracy, the higher the importance of the variable in the model. 
# In other words, mean decrease accuracy represents how much removing each variable reduces the accuracy of the model.
varImpPlot(rf_model, n.var = 5, main = "Top 10 - Variable") # It displays a variable importance plot from the random forest model. 

importance(rf_model)                                        # It is worth mentioning that IncNodePurity by how much the model error increases by dropping each of the specified explanatory variables. 

##          IncNodePurity
## age          148.47066
## bmi           19.23523
## children      19.81647
## smoker        69.04689

                                                            # Briefly, varImpPlot() indicates each variable's importance in explaining the performance of the dependent variable (Y).

RMSE

RMSE_rf

## [1] 0.4324175

Redes Neuronales

df2 <- na.omit(df2)  # Limpia el dataframe de NA
df2_alt <- df2 %>% dplyr::select(expenses, age,bmi,smoker,children)

df2$expenses= log(df2$expenses )
df2$age= log(df2$age )
df2$bmi= log(df2$bmi )

set.seed(123)
partition <- createDataPartition(y = df2$expenses, p = 0.7, list = FALSE)
train = df2[partition, ]
test = df2[-partition, ]
train2 = df2[partition, ]
test2 = df2[-partition, ]

# Lets estimate a Neural Network Regression

train2 <- train2 %>%
  mutate(smoker = if_else(smoker == "yes", 1, 0))
test2 <- test2 %>%
  mutate(smoker = if_else(smoker == "yes", 1, 0))

nn_model <- neuralnet(expenses ~ age +  bmi + children + smoker, data = train2, hidden = c(5, 3), linear.output = TRUE) 

# Plot the neural network 
plot(nn_model)


nn_predictions <- neuralnet::compute(nn_model, test2)
nn_predictions_values <- nn_predictions$net.result


actual_values <- test2$expenses


# Calcular el RMSE
rmse_NN <- Metrics::rmse(actual_values, nn_predictions_values)

RMSE

# Imprimir el RMSE
print(rmse_NN)

## [1] 0.03888752

Evaluación y selección del modelo

RMSE de todos los modelos

# Crear el data.frame
rmse_df <- data.frame(
  Modelo = c("RMSE de log_ols_model", "RMSE de log_ols_wls_model", "RMSE de XGB", "RMSE de Tree", "RMSE de RF", "RMSE de NN"),
  RMSE = c(RMSE_log_ols_model, RMSE_log_ols_wls_model, RMSE_XGB, RMSE_tree, RMSE_rf, rmse_NN)
)

# Imprimir el data.frame creado
print(rmse_df)

##                      Modelo       RMSE
## 1     RMSE de log_ols_model 0.43517072
## 2 RMSE de log_ols_wls_model 0.39095301
## 3               RMSE de XGB 1.00589420
## 4              RMSE de Tree 0.37617398
## 5                RMSE de RF 0.43241751
## 6                RMSE de NN 0.03888752

# Crear el gráfico de barras
ggplot(rmse_df, aes(x = Modelo, y = RMSE, fill = Modelo)) +
  geom_bar(stat = "identity", show.legend = FALSE) + 
  theme_minimal() +
  theme(axis.text.x = element_text(angle = 45, hjust = 1)) +
  labs(title = "Comparación de RMSE por Modelo",
       x = "Modelo",
       y = "RMSE") +
  geom_text(aes(label = round(RMSE, 2)), vjust = -0.5, size = 3.5)

Conclusiones

Dado que los resultados obtenidos al evaluar la métrica RMSE indican que el modelo de redes neuronales tiene el menor valor, se decide seleccionarlo como el modelo que presenta el mejor desempeño en términos de precisión según esta métrica.

EDA

Tras nuestro análisis exhaustivo de los datos, identificamos que varias variables no seguían distribuciones normales, por lo que recurrimos a transformaciones logarítmicas para normalizarlas. También detectamos la presencia de outliers que afectaban negativamente la interpretación de nuestros modelos. Este análisis detallado nos brindó una comprensión más profunda de cada variable y su relación con la variable dependiente, lo que facilitó la creación de un modelo más preciso y con una mejor capacidad predictiva. Además, observamos una baja correlación entre las variables, lo que descartó preocupaciones sobre multicolinealidad y respaldó la suposición de independencia entre las variables explicativas.

Sin embargo, debido a las irregularidades en las distribuciones de las variables independientes y la naturaleza transversal de los datos, anticipamos posibles problemas de heteroscedasticidad. Estas preocupaciones se confirmaron durante la estimación de los modelos de regresión lineal. Para abordar estos desafíos, aplicamos transformaciones a las variables y utilizamos un modelo de regresión lineal ponderado, lo que mejoró la fiabilidad de los resultados y mitigó los problemas asociados con la heteroscedasticidad.

Modelo seleccionado

i. ¿Cuáles son las variables que contribuyen a explicar los cambios de la principal variable de estudio?

Según los coeficientes obtenidos del modelo de redes neuronales, podemos inferir que las variables que ejercen una mayor influencia en la predicción de la variable dependiente son age, smoke y bmi.

ii. ¿Cómo es el impacto de dichas variables explicativas sobre la variable dependiente?

Todas las variables mencionadas anteriormente tienen una influencia significativa en la variable dependiente y su efecto en la variable a predecir es positivo. Esto significa que a medida que estas variables aumentan, también aumenta la variable dependiente (expenses).

iii. ¿Los resultados estimados del modelo seleccionado son similares a los otros modelos estimados? ¿Cuáles son las diferencias?

Los resultados de todos los modelos son similares cuando se comparan con sus principales métricas de evaluación de resultados, en los modelos de regresión lineal, los AIC son similares y se mejoró esta mética con el modelo de regresión ponderado por los motivos de transformaciones y de la ponderación de las variables. Por otro lado, todos los modelos muestran un RMSE similar, estando todas en un rango similar de valores. Por último, las variables significativas en los modelos coinciden al colocar a age como la variable con más peso, seguido de smoke y bmi.