# Cargar paquetes
# install.packages(c("tidymodels", "tidyverse", "rpart")) # Instalar si es necesario
library(tidymodels)
library(tidyverse)
library(rpart) # Para el motor de árbol de decisión

# Configurar tema de ggplot
theme_set(theme_minimal(base_size = 14))

paquetes <- c("tidymodels", "tidyverse", "rpart")

# require(): intenta cargar el paquete y devuelve TRUE/FALSE
# Si no está instalado (FALSE), lo instala y carga
# character.only = TRUE: pkg es una variable; sin esto, R buscaría un paquete llamado "pkg"
for (pkg in paquetes) {
  if (!require(pkg, character.only = TRUE)) {
    install.packages(pkg, dependencies = TRUE)
    library(pkg, character.only = TRUE)
  }
}

# Cargar datos
datos <- read_csv("datos/indicadores_mundiales.csv")
# Lean el archivo directo desde la web:
# datos <- read_csv("https://raw.githubusercontent.com/danilofreire/introduccion-ia-ucu/main/clases/dia-01/datos/indicadores_mundiales.csv")

# Ver estructura
glimpse(datos) # Ver las variables y tipos

# Resumen rápido de las variables numéricas
datos |> select(gasto_educacion, acceso_internet, urbanizacion) |> summary()

# Medias por nivel de crecimiento (3 predictores ilustrativos)
datos |>
  group_by(crecimiento_alto) |>
  summarise(
    educacion_media     = mean(gasto_educacion),
    internet_medio      = mean(acceso_internet),
    urbanizacion_media  = mean(urbanizacion)
  )

# Frame listo para modelar (datos sigue intacto con pais y continente)
datos_modelo <- datos |>
  mutate(crecimiento_alto = factor(crecimiento_alto, levels = c("no", "si"))) |>
  select(-pais, -continente)

# Verificar
glimpse(datos_modelo)

# Crear nuevas variables sobre el frame que va al modelo
datos_modelo <- datos_modelo |>
  mutate(
    # Grupos de urbanización (rural, mixto, urbano)
    grupo_urbanizacion = cut(urbanizacion,
                             breaks = c(0, 50, 75, 100),
                             labels = c("rural", "mixto", "urbano"),
                             include.lowest = TRUE),

    # Intensidad fiscal social: cuánto del PIB va a salud + educación
    gasto_social = gasto_educacion + gasto_salud,

    # Indicador de alto desempleo (lo usaremos como ejemplo más adelante)
    alto_desempleo = factor(if_else(desempleo > 7, "alto", "bajo"))
  )

# Ver las nuevas variables
datos_modelo |> select(grupo_urbanizacion, gasto_social, alto_desempleo) |> head()

# Ejemplo: crear variable internet_grupos por terciles
datos_modelo <- datos_modelo |>
  mutate(
    internet_grupos = cut(acceso_internet,
                          breaks = quantile(acceso_internet, c(0, 1/3, 2/3, 1)),
                          labels = c("bajo", "medio", "alto"),
                          include.lowest = TRUE)
  )

# Fijar semilla
set.seed(2026)

# Dividir datos (usando el frame con la variable objetivo como factor)
datos_split <- initial_split(datos_modelo, prop = 0.75, strata = crecimiento_alto)
datos_train <- training(datos_split)
datos_test <- testing(datos_split)

# Crear folds para validación cruzada
folds <- vfold_cv(datos_train, v = 5, strata = crecimiento_alto)

cat("Train:", nrow(datos_train), "| Test:", nrow(datos_test))

# 1. Regresión logística (lineal)
modelo_logistico <- logistic_reg() |>
  set_engine("glm") |>
  set_mode("classification")

# 2. Árbol de decisión (no lineal, captura interacciones automáticamente)
modelo_arbol <- decision_tree() |>
  set_engine("rpart") |>
  set_mode("classification")

# Fórmula con las 8 variables numéricas originales
formula_basica <- crecimiento_alto ~ gasto_educacion + acceso_internet +
                                     urbanizacion + gasto_salud + inflacion +
                                     desempleo + inversion_extranjera +
                                     indice_gobierno_digital

# Las tres métricas juntas en un objeto que vamos a reusar todo el lab
mis_metricas <- metric_set(precision, recall, roc_auc)

eval_logistico <- fit_resamples(
  modelo_logistico, formula_basica,
  resamples = folds,
  metrics = mis_metricas,
  control = control_resamples(event_level = "second")   # "si" es la clase positiva
) |>
  collect_metrics() |>
  mutate(modelo = "Logístico")

eval_logistico

# El mismo bloque que antes, cambiando modelo_logistico → modelo_arbol
eval_arbol <- fit_resamples(
  modelo_arbol, formula_basica,
  resamples = folds,
  metrics = mis_metricas,
  control = control_resamples(event_level = "second")
) |>
  collect_metrics() |>
  mutate(modelo = "Árbol")

# Combinar y mostrar la comparación
resultados <- bind_rows(eval_logistico, eval_arbol)

# Una fila por modelo y métrica, ordenadas para comparar de a pares
resultados |>
  select(modelo, .metric, mean, std_err) |>
  arrange(.metric, modelo)

# Fórmula extendida con un feature derivado
formula_ext <- crecimiento_alto ~ acceso_internet + gasto_educacion + urbanizacion + gasto_salud + inflacion +
                                  alto_desempleo + inversion_extranjera + indice_gobierno_digital

# Logístico con formula_ext
eval_log_ext <- fit_resamples(modelo_logistico, formula_ext, resamples = folds,
  metrics = mis_metricas,
  control = control_resamples(event_level = "second")) |>
  collect_metrics() |> mutate(modelo = "Logístico (ext)")

# Comparar el logístico básico vs. el extendido
bind_rows(eval_logistico, eval_log_ext) |>
  select(modelo, .metric, mean, std_err) |>
  arrange(.metric, modelo)

# Entrenar el modelo elegido con todos los datos de train
ajuste_final <- modelo_logistico |>
  fit(formula_ext, data = datos_train)

# augment() junta predicciones (clases + probabilidades) con los datos originales
# en una sola línea. Equivale a hacer dos predict() + bind_cols().
pred_test <- ajuste_final |> augment(datos_test)

# Métricas finales: el mismo mis_metricas de la validación cruzada
pred_test |> mis_metricas(truth = crecimiento_alto, estimate = .pred_class, .pred_si, event_level = "second")

# 1. Dimensiones
dim(datos)

# 2. Distribución de la variable objetivo
datos |>
  count(crecimiento_alto) |>
  mutate(prop = round(n / sum(n), 3))

# 3. Valores faltantes por variable
colSums(is.na(datos))

# 4. Resumen de variables numéricas
datos |>
  select(where(is.numeric)) |>
  summary()

# Estadísticas detalladas
summary(datos$inflacion)

# Histograma (muestra la cola larga)
ggplot(datos, aes(x = inflacion)) +
  geom_histogram(bins = 30, fill = "#2d4563") +
  labs(x = "Inflación anual (%)", y = "Países")

# Densidad por grupo (¿separa crecimiento alto de bajo?)
ggplot(datos, aes(x = inflacion, fill = crecimiento_alto)) +
  geom_density(alpha = 0.5) +
  scale_fill_manual(values = c("#e74c3c", "#27ae60")) +
  labs(x = "Inflación anual (%)")

# Medias por grupo de las 4 variables no mostradas en el cuerpo
datos |>
  group_by(crecimiento_alto) |>
  summarise(
    salud         = mean(gasto_salud),
    inflacion     = mean(inflacion),
    desempleo     = mean(desempleo),
    ied           = mean(inversion_extranjera),
    gob_digital   = mean(indice_gobierno_digital)
  )

# Correlación entre todos los predictores numéricos (¿hay pares redundantes?)
datos |>
  select(where(is.numeric)) |>
  cor() |>
  round(2)

# Test formal para una variable concreta
t.test(indice_gobierno_digital ~ crecimiento_alto, data = datos)

datos |>
  select(crecimiento_alto, gasto_salud, desempleo, inversion_extranjera, indice_gobierno_digital) |>
  pivot_longer(-crecimiento_alto, names_to = "variable", values_to = "valor") |>
  ggplot(aes(x = crecimiento_alto, y = valor, fill = crecimiento_alto)) +
  geom_boxplot(alpha = 0.7) +
  facet_wrap(~variable, scales = "free_y") +
  scale_fill_manual(values = c("#e74c3c", "#27ae60")) +
  labs(title = "Las 4 variables no mostradas, por nivel de crecimiento", x = NULL, y = NULL) +
  theme(legend.position = "none")

ggplot(datos, aes(x = urbanizacion, y = acceso_internet, color = crecimiento_alto)) +
  geom_point(alpha = 0.6) +
  geom_smooth(method = "lm", se = FALSE) +
  scale_color_manual(values = c("#e74c3c", "#27ae60")) +
  labs(x = "Urbanización (%)", y = "Acceso a internet (%)")

# La correlación numérica complementa el gráfico
cor(datos$urbanizacion, datos$acceso_internet)

# Una posible solución para cada una de las cuatro opciones sugeridas
datos_modelo <- datos_modelo |>
  mutate(
    # 1. Gasto en educación alto (gasto_educacion)
    educacion_alta = if_else(gasto_educacion > median(gasto_educacion), "alta", "baja"),

    # 2. Acceso a internet por terciles (acceso_internet)
    internet_grupos = cut(acceso_internet,
                          breaks = quantile(acceso_internet, c(0, 1/3, 2/3, 1)),
                          labels = c("bajo", "medio", "alto"),
                          include.lowest = TRUE),

    # 3. Inflación alta (inflacion)
    inflacion_alta = if_else(inflacion > 15, "si", "no"),

    # 4. Apertura financiera: IED por encima de la mediana (inversion_extranjera)
    apertura_financiera = if_else(inversion_extranjera > median(inversion_extranjera), "alta", "baja")
  )

# Verificar las nuevas variables
datos_modelo |> count(educacion_alta)

datos_modelo |> count(internet_grupos)

datos_modelo |> count(inflacion_alta)

datos_modelo |> count(apertura_financiera)

# Una interacción no necesita una columna nueva: se agrega a la fórmula con ":".
# Probamos indice_gobierno_digital × acceso_internet (dos variables que se refuerzan)
formula_int <- crecimiento_alto ~ gasto_educacion + acceso_internet +
                                  urbanizacion + gasto_salud + inflacion +
                                  desempleo + inversion_extranjera +
                                  indice_gobierno_digital +
                                  indice_gobierno_digital:acceso_internet

# Evaluar el logístico (reusamos los mismos folds del cuerpo del lab)
eval_log_int <- fit_resamples(
  modelo_logistico, formula_int, resamples = folds,
  metrics = mis_metricas,
  control = control_resamples(event_level = "second")
) |>
  collect_metrics() |> mutate(modelo = "Logístico (interacción)")

# Evaluar el árbol con la misma fórmula
eval_arb_int <- fit_resamples(
  modelo_arbol, formula_int, resamples = folds,
  metrics = mis_metricas,
  control = control_resamples(event_level = "second")
) |>
  collect_metrics() |> mutate(modelo = "Árbol (interacción)")

bind_rows(eval_logistico, eval_log_int, eval_arbol, eval_arb_int) |>
  filter(.metric == "roc_auc") |>
  select(modelo, mean, std_err) |>
  arrange(desc(mean))

resultados |>
  ggplot(aes(x = modelo, y = mean, fill = modelo)) +
  geom_col(alpha = 0.8) +
  geom_errorbar(aes(ymin = mean - std_err, ymax = mean + std_err),
                width = 0.2) +
  facet_wrap(~.metric, scales = "free_y") +
  scale_fill_manual(values = c("#2d4563", "#27ae60")) +
  labs(title = "Comparación de modelos", y = "Valor", x = "") +
  theme(legend.position = "none")