Algoritmos de Aprendizaje Computacional

1 Introducción

1.1 Objetivos

El principal objetivo de este estudio es obtener un modelo ajustado al conjunto de datos de Metástasis de Cáncer de Mama para poder realizar predicciones efectivas sobre el estado PCR de un paciente dada una serie de variables o características sobre este. Consecuentemente, nace de manera natural un segundo objetivo principal, y es el de la determinación de las características más importantes para la determinación del estado del paciente, es decir, la realización de una selección de características. Todo esto se realizará con la finalidad de poner un modelo en producción, el cual pueda ser usado a través de la interfaz RShiny para la predicción efectiva del estado PCR de un paciente de Cáncer de Mama.

Objetivos secundarios de este proyecto incluyen el estudio de la estabilidad en la selección de características a lo largo de los diferentes métodos y modelos aplicados. Se realizará también una exploración de la estabilidad de los hiperparámetros de un mismo modelo a lo largo de diferentes semillas aleatorias establecidas.

Objetivo	Tipo
Obtener un modelo ajustado al conjunto de datos de Metástasis de Cáncer de Mama	Principal
Determinación de características más importantes en la predicción del estado PCR del paciente	Principal
Poner un modelo en producción, accesible a través de la interfaz RShiny	Principal
Estudiar la estabilidad en las características seleccionadas	Secundario
Explorar la estabilidad de los hiperparámetros de un mismo modelo a lo largo de diferentes semillas aleatorias establecidas.	Secundario

1.2 Trabajos Previos

2 Metodología

A lo largo del desarrollo de este proyecto se han utilizado una serie de datos y se han seguido unas directrices ordenadas. Todo esto queda redactado en este apartado.

2.1 Material

2.1.1 Conjunto de Datos Inicial

Se nos ha proporcionado un fichero en formato Comma-Separated Values (.csv) gracias al trabajo del Dr. José Manuel Jérez Aragonés. Las características de este dataset quedan descritas en la Figura 1.

Figura 1: Variables Originales del Conjunto

Se puede identificar la variable PCR como la variable target u objetivo de este estudio. Consecuentemente, se ha decidido eliminar la variable Muestra, ya que no parece guardar una relevancia clínica, al ser un identificador único para cada paciente, compuesto por un código identificativo del estudio en el que fue obtenido el dato, y una serie de números únicos para cada persona.

2.1.2 Pre-Procesamiento

El conjunto de datos del que se hará uso fue curado en un estudio anterior a este, habiéndose imputado los valores faltantes (NA) por la mediana o la moda, si estos provenían de una característica numérica, o de tipo factor, respectivamente (González 2024).

Por otra parte, para asegurar que se tiene total control sobre los datos con los que se introducen al modelo, se han creado dos versiones del conjunto de datos. La primera se trata de una versión en la que todas las variables son de tipo factor, y la variable edad ha sido estratificada en diferentes intervalos. Por otra parte, se ha generado un conjunto de datos numérico, el cual es el conjunto de datos de tipo factor codificado mediante técnicas como One-Hot Encoder, o manteniendo las variables con su mismo tipo, si estas eran numéricas originalmente. Se comentan a continuación las modificaciones específicas.

2.1.2.1 Conjunto de Datos Factor

Este conjunto de datos existe para asegurar la interpretabilidad de los resultados de los modelos. La única variable que puede suponer un problema es la variable Edad, ya que es una variable continua, sin embargo, este problema, como se ha comentado anteriormente, se resolverá estratificándola. La anchura de los bins se calculará siguiendo la regla de Sturges, en la que se toma el logaritmo en base 2 de la cantidad de muestras y el rango de la variable para poder computar la anchura más adecuada.

k = \frac{max(Edad) - min(Edad)}{\log_2(n) + 1}

# Convertir la variable 'Edad' a numérica
data_original$Edad <- as.numeric(data_original$Edad)
data_factor <- data_original

# Calculando el número de bins con la regla de Sturges
n <- nrow(data_factor)
k <- ceiling(log2(n) + 1)

# Calculando el rango y el tamaño de cada bin
minimo <- min(data_factor$Edad)
maximo <- max(data_factor$Edad)
ancho_bin <- (maximo - minimo) / k

# Estratificar la variable 'Edad' en grupos
data_factor$Edad_estratificada <- cut(data_factor$Edad, 
                                      breaks = seq(min(data_factor$Edad), max(data_factor$Edad), by = ancho_bin), 
                                      include.lowest = TRUE, right = FALSE)

# Asignar nombres de grupo a los niveles de Edad_estratificada
levels(data_factor$Edad_estratificada) <- paste0("Grupo", seq_along(levels(data_factor$Edad_estratificada)))

# Convertir todas las variables del dataframe a factor, excepto 'Edad'
data_factor <- data_factor %>%
  mutate_at(vars(-Edad), as.factor)

# Actualizar la variable Edad con los grupos
data_factor$Edad <- data_factor$Edad_estratificada
data_factor$Edad_estratificada <- NULL

knitr::kable(head(data_factor, 10))

Edad	REst	RPro	Her2	Estadio	NodAfec	Grado	Fenotipo	PCR
Grupo3	P	P	N	T2	N1	II	LumA	0
Grupo5	P	P	N	T3	N1	II	Normal	0
Grupo4	N	N	N	T3	N0	III	Basal	0
Grupo4	N	N	N	T2	N1	III	Basal	1
Grupo3	P	P	N	T3	N1	II	LumA	0
Grupo6	N	N	N	T4	N0	III	Basal	0
Grupo3	N	N	N	T2	N0	III	Basal	0
Grupo6	P	N	N	T3	N1	III	LumA	0
Grupo5	N	N	N	T3	N3	III	Basal	0
Grupo4	P	P	P	T4	N3	II	LumA	0

Esta estratificación se puede justificar con el hecho de que la interpretación de las predicciones del modelo puede ser más adecuada, ya que no estaríamos hablando de un valor concreto de edad que se repite a la hora de analizar los pacientes PCR-positivos, sino que tal vez se podría identificar una tendencia en los pacientes de un determinado rango de edad a sufrir una metástasis.

Grupo	Rango de Edad
Grupo 1	[24,29.1)
Grupo 2	[29.1,34.2)
Grupo 3	[34.2,39.3)
Grupo 4	[39.3,44.4)
Grupo 5	[44.4,49.5)
Grupo 6	[49.5,54.6)
Grupo 7	[54.6,59.7)
Grupo 8	[59.7,64.8)
Grupo 9	[64.8,69.9)
Grupo 10	[69.9,75)

2.1.2.2 Conjunto de Datos Numérico

El paradigma de este estudio engloba la búsqueda del algoritmo de aprendizaje computacional que mejor prediga el estado PCR de un paciente no visto, todo desde un punto de vista de la Ingeniería. Debido a esta falta de conocimiento clínico, se ha decidido no asumir un orden de las variables factor a la hora de encontrar su codificación numérica (por ejemplo, se ha decidido no asumir que un valor de Grado: II es más importante o severo que un Grado: I). Es port esto que se ha decidido usar la función dummyVars para realizar una codificación One-Hot a las variables. Este conjunto de datos existe para poder tener siempre presente los datos numéricos que se pasan a los modelos de ML, así como para que se use en algunas implementaciones concretas, sin embargo, la regla será usar el conjunto data_factor como entrada.

# Generación del conjunto de datos numérico
encode_data <- function(data, target_name) {
  # Convertir variables categóricas sin orden en dummies
  data <- dummyVars(~ ., data = data, fullRank = TRUE) %>% predict(data) %>% as.data.frame()
  # Asegurar que la variable objetivo esté al final
  target <- data[[target_name]]
  data[[target_name]] <- NULL
  data[[target_name]] <- target
  return(data)
}

data_numeric <- encode_data(data_original, "PCR")

knitr::kable(head(data_numeric, 10))

Edad	REstP	RProP	Her2P	EstadioT2	EstadioT3	EstadioT4	NodAfecN1	NodAfecN2	NodAfecN3	GradoII	GradoIII	FenotipoHer2	FenotipoLumA	FenotipoLumB	FenotipoNormal	PCR
37.8	1	1	0	1	0	0	1	0	0	1	0	0	1	0	0	0
45.8	1	1	0	0	1	0	1	0	0	1	0	0	0	0	1	0
40.7	0	0	0	0	1	0	0	0	0	0	1	0	0	0	0	0
40.8	0	0	0	1	0	0	1	0	0	0	1	0	0	0	0	1
35.5	1	1	0	0	1	0	1	0	0	1	0	0	1	0	0	0
52.2	0	0	0	0	0	1	0	0	0	0	1	0	0	0	0	0
38.2	0	0	0	1	0	0	0	0	0	0	1	0	0	0	0	0
54.2	1	0	0	0	1	0	1	0	0	0	1	0	1	0	0	0
46.6	0	0	0	0	1	0	0	0	1	0	1	0	0	0	0	0
40.8	1	1	1	0	0	1	0	0	1	1	0	0	1	0	0	0

2.1.2.3 Análisis Exploratorio

Se va a proceder a realizar un análisis de la distribución de los datos del conjunto data_factor. Este análisis pretende dar conclusiones sobre el sesgo y la varianza de los datos con la finalidad de tomar decisiones sobre la inclusión o deleción de variables del conjunto de datos.

blue = '#377eb8'
red = '#e41a1c'
plot_variable_distribution <- function(data, variable_name, target_name) {
  # Comprobar si las variables existen en el data.frame
  if (!(variable_name %in% names(data) && target_name %in% names(data))) {
    stop("Una o ambas variables especificadas no existen en el data.frame proporcionado.")
  }
  
  # Extraer la variable y la variable target del data.frame
  variable <- data[[variable_name]]
  target <- data[[target_name]]
  
  # Asegurar que la variable target es factor
  data[[target_name]] <- as.factor(data[[target_name]])
  
  # Determinar si la variable principal es numérica o factor
  if (is.numeric(variable)) {
    # Crear un histograma para variables numéricas
    p <- ggplot(data, aes_string(x = variable_name, fill = target_name)) +
      geom_histogram(position = "identity", alpha = 0.5, bins = 30) +
      scale_fill_manual(values = c("0" = blue, "1" = red)) +
      ggtitle(paste("Distribución de", variable_name, "por", target_name)) +
      xlab(variable_name) +
      ylab("Frecuencia")
  } else if (is.factor(variable)) {
    # Crear un gráfico de barras para factores
    p <- ggplot(data, aes_string(x = variable_name, fill = target_name)) +
      geom_bar(position = "stack") +
      scale_fill_manual(values = c("0" = blue, "1" = red)) +
      ggtitle(paste("Distribución de", variable_name, "por", target_name)) +
      xlab(variable_name) +
      ylab("Frecuencia") +
      theme(legend.position = "right")
  } else {
    stop("El tipo de la variable debe ser numérico o factor.")
  }
  
  # Imprimir el gráfico
  print(p)
}

Edad

plot_variable_distribution(data_factor, "Edad", "PCR")

Como se puede observar, la variable Edad parece seguir una distribución normal o gaussiana, en la que la mayor cantidad de los datos se encuentran en los grupos centrales. La mayoría de los datos están clasificados por PCR-Positivos (PCR=1), indicando un claro sesgo, sin embargo, se considera que hay suficientes datos en la variable como para mantenerla.

REst

plot_variable_distribution(data_factor, "REst", "PCR")

La variable REst muestra dos clases, para las cuales hay una diferencia de datos de, aproximadamente, 100 muestras. La clase N muestra un mayor equilibrio en cuanto a la relación de pacientes PCR-Positivos y PCR-Negativos que el valor P, para el cual, solo paroximadamente un sexto de los datos son PCR-Positivos.

RPro

plot_variable_distribution(data_factor, "RPro", "PCR")

La variable RPro muestra dos clases (N y P). Ambas clases tienen una distribución similar en términos de número de muestras totales. La clase N presenta una mayor proporción de pacientes PCR-Positivos (PCR=1) en comparación con la clase P.

Her2

plot_variable_distribution(data_factor, "Her2", "PCR")

La variable Her2 muestra dos clases (N y P). La clase N tiene una distribución mucho mayor de muestras en comparación con la clase P. En la clase N, hay una cantidad significativa de pacientes PCR-Positivos (PCR=1), mientras que en la clase P, los datos son escasos, lo que puede dificultar análisis concluyentes, es por esto que se ha decidido eliminar la variable Her2 de los conjuntos de datos.

data_numeric$Her2P <- NULL
data_factor$Her2 <- NULL
data_original$Her2 <- NULL

Estadio

plot_variable_distribution(data_factor, "Estadio", "PCR")

La variable Estadio presenta varias categorías (T0-T1, T2, T3, T4). La mayoría de las muestras se encuentran en la categoría T2, seguida de T3. Las categorías T0-T1 y T4 tienen menos muestras. En todas las categorías, los pacientes PCR-Negativos (PCR=0) son mayoría, pero la proporción de PCR-Positivos (PCR=1) varía entre las categorías, siendo notablemente menor en T4.

NodAfec

plot_variable_distribution(data_factor, "NodAfec", "PCR")

La variable NodAfec muestra varias categorías (N0, N1, N2, N3). La categoría N1 tiene la mayor cantidad de muestras, seguida de N0, con N2 y N3 teniendo menos muestras. En todas las categorías, los pacientes PCR-Negativos (PCR=0) son mayoría. La proporción de PCR-Positivos (PCR=1) es relativamente constante en todas las categorías, aunque menor en comparación con los PCR-Negativos.

Grado

plot_variable_distribution(data_factor, "Grado", "PCR")

La variable Grado tiene tres niveles (I, II, III). El nivel III tiene la mayor cantidad de muestras, seguido por II. El nivel I tiene significativamente menos muestras. En los niveles II y III, hay una proporción considerable de pacientes PCR-Positivos (PCR=1), especialmente en III.

Fenotipo

plot_variable_distribution(data_factor, "Fenotipo", "PCR")

La variable Fenotipo tiene varias categorías (Basal, Her2, LumA, LumB, Normal). La categoría Basal tiene la mayor cantidad de muestras, seguida de LumA. Las categorías Her2, LumB y Normal tienen menos muestras. En Basal, hay una proporción significativa de pacientes PCR-Positivos (PCR=1), mientras que en las otras categorías, la proporción es menor.

PCR

plot_variable_distribution(data_factor, "PCR", "PCR")

La variable PCR muestra la distribución de los pacientes en las dos clases (0 y 1). La mayoría de las muestras son PCR-Negativos (PCR=0), con una menor proporción de PCR-Positivos (PCR=1). Esta distribución sugiere un sesgo hacia los resultados negativos, lo cual se tendrá en cuenta en los análisis posteriores.

2.2 Métodos

Todos los objetivos de este proyecto quedarían anulados si no se sigue una metodología estricta, basada en la estimación honesta de parámetros y la rigurosidad para la evaluación del rendimiento de los modelos mediante la selección de una métrica que sea capaz de cuantificar la capacidad predictiva del modelo independientemente del umbral seleccionado para clasificar las muestras -ya que este depende de las necesidades específicas del equipo sanitario. Esta sección sirve como guía para mostrar cuál ha sido esta metodología adoptada, y así demostrar su eficacia y seguridad.

2.2.1 AUC y Rendimiento Aparente

Para poder estimar el rendimiento cualitativo de un modelo, se necesita tanto una métrica cuantificable que refleje efectivamente su capacidad predictiva, como un valor base (ground truth) con el que comparar. En este estudio se ha decidido utilizar el AUC (Area Under the Curve) como métrica representativa del rendimiento de un modelo, porque es una medida independiente del umbral utilizado para clasificar probabilidades calculadas por el modelo, esto es una característica decisiva, ya que, como se comentó anteriormente, el umbral podrá ser decidido por el clínico cuando el modelo sea puesto en producción, para que se ajuste a sus necesidades. Al utilizar el AUC como medida, podremos dar una estimación independiente del umbral a los clínicos, indicando que, elijan el umbral que elijan, el modelo rendirá de una manera preestablecida.

Por otra parte, se ha apostado por tomar el rendimiento aparente del modelo como la métrica base de comparación para un mismo modelo. El rendimiento aparente se calcula mediante el entrenamiento y evaluación del modelo con todo el conjunto de datos, obteniendo una estimación optimista de la capacidad predictiva del modelo, al ser evaluado con los datos que fue entrenado. Una característica que debe quedar clara es la de que esta métrica solo servirá para comparar diferentes métodos de selección de variables, búsqueda de hiperparámetros, o ajuste fino, de un mismo modelo, no entre modelos. Para realizar una comparación entre modelos se hará uso de la estimación honesta de parámetros, es decir, la doble validación cruzada 5x2.

2.2.2 Ajuste Fino: Doble Validación Cruzada 5x2

Doble Validación Cruzada

El método de doble validación cruzada 5x2 es una técnica robusta para la evaluación de modelos que busca obtener una estimación honesta del rendimiento predictivo de un modelo. Este método se basa en una estructura jerárquica de validación cruzada, compuesta por dos niveles conocidos como inner loop y outer loop.

2.2.2.1 Outer Loop

El outer loop se encarga de la evaluación externa del modelo y se basa en la división del conjunto de datos original en 5 pliegues. En cada iteración, uno de estos pliegues se reserva como conjunto de test, mientras que los cuatro pliegues restantes se utilizan para el entrenamiento del modelo. Este proceso se repite cinco veces, de manera que cada pliegue se utiliza una vez como conjunto de test. Es crucial destacar que en el outer loop no se permite ajustar parámetros ni realizar selección de variables utilizando el conjunto de test, garantizando así una evaluación imparcial y honesta del modelo.

2.2.2.2 Inner Loop

Dentro de cada iteración del outer loop, se implementa el inner loop, cuyo objetivo es el ajuste fino de los hiperparámetros del modelo y la selección de las variables más relevantes. El inner loop utiliza una nueva división del conjunto de entrenamiento en pliegues adicionales. Generalmente, se opta por una validación cruzada más sencilla dentro de este bucle interno. Aquí, uno de los pliegues del conjunto de entrenamiento actúa como validación, mientras que los restantes se utilizan para el ajuste del modelo. Este proceso se repite varias veces para identificar la combinación óptima de hiperparámetros y características que maximicen el rendimiento predictivo en el conjunto de validación.

2.2.2.3 Evaluación Final

Una vez completado el ajuste en el inner loop, se evalúa el modelo ajustado en el conjunto de test reservado en el outer loop. Este enfoque jerárquico asegura que el ajuste de los hiperparámetros y la selección de variables no influyan en la evaluación final del rendimiento del modelo. El rendimiento se cuantifica mediante métricas como el AUC (Área Bajo la Curva), proporcionando una estimación promedio del rendimiento del modelo a lo largo de las iteraciones del outer loop. Este método de doble validación cruzada 5x2 ofrece una estimación más confiable y menos optimista del rendimiento del modelo, adecuada para su implementación en escenarios clínicos donde la precisión es crucial.

2.2.3 Selección de Características

La selección de características es un paso fundamental en la construcción de modelos predictivos, ya que permite identificar las variables más relevantes para el problema en cuestión. Este proceso no solo mejora la interpretabilidad del modelo, cosa que permite a los investigadores y clínicos entender mejor qué variables están influyendo en las predicciones del modelo, sino que también reduce significativamente el tiempo computacional y el riesgo de sobreajuste, proporcionando un modelo más eficiente y robusto.

2.2.3.1 Metodología de Selección de Características

En este estudio, se han realizado cuatro esquemas de selección de características por cada ciclo de doble validación cruzada 5x2 (2CV 5X2). Estos esquemas incluyen tanto métodos de filtrado, como métodos wrapped, y métodos embedded, proporcionando una perspectiva completa sobre las variables más relevantes.

1. Todas las Variables del Conjunto: En este enfoque, se entrena el modelo utilizando todas las variables disponibles en el conjunto de datos, sin realizar ningún tipo de selección previa. Esto sirve como una línea base para comparar la efectividad de los métodos de selección de características.

2. Métodos de Filtrado: Análisis de Asociación: Este método de filtrado se basa en el análisis de asociación, donde se evalúa la relación individual de cada variable con la variable objetivo. Se seleccionan aquellas variables que muestran una fuerte asociación, lo que permite reducir el conjunto de características antes del entrenamiento del modelo. Para poder realizar este análisis de asociación se inspeccionará el mapa de calor proporcionado por la función plot_p_valores. Como se puede observar, las variables seleccionadas son Estadio, Fenotipo, Grado, REst y RPro. Por lo que estas variables constituirán el vector de características asociacion, usado para la ejecución de los métodos de filtrado en los modelos.

source("./aux_scripts/calculaPValor.R")
plot <- plot_p_valores(data_factor)
print(plot)

# Se construye el vector de variables factor y numéricas asociación a partir de los resultados del análisis

todas_variables_factor <- colnames(data_factor)[!colnames(data_factor) %in% "PCR"]
todas_variables_numeric <- colnames(data_numeric)[!colnames(data_numeric) %in% "PCR"]

select_associated_variables <- function(asociacion, total_variables_numeric) {
  selected_variables <- c()
  
  for (assoc in asociacion) {
    # Buscar variables que contengan el nombre en asociacion
    matching_vars <- total_variables_numeric[grepl(assoc, total_variables_numeric)]
    selected_variables <- c(selected_variables, matching_vars)
  }
  
  return(selected_variables)
}

asociacion <- c("Estadio", "Fenotipo", "Grado", "REst", "RPro")
asociacion_numeric <- select_associated_variables(asociacion, todas_variables_numeric)

3. Métodos Wrapped: StepAUC Backwards: En este método, se comienza con todas las variables incluidas en el modelo. Durante el proceso de entrenamiento en el inner loop, se prueban eliminando variables una por una. Si la eliminación de una variable mejora el AUC (Área Bajo la Curva), la variable se elimina permanentemente del conjunto. Este proceso se repite hasta que la eliminación de más variables no resulta en una mejora del AUC, determinando así el conjunto final de características seleccionadas.

4. Métodos de Embedded: Modelo Lasso: El modelo Lasso es un método de regularización que no solo ajusta el modelo, sino que también selecciona automáticamente las variables más importantes. Durante el inner loop, se ajusta un modelo Lasso y se seleccionan las variables con coeficientes diferentes de cero.

2.2.4 Modelos Aplicados

Con el objetivo de conseguir un modelo ajustado para la predicción efectiva del estado PCR de un paciente en mente, se han seleccionado una serie de algoritmos de Aprendizaje Máquina (Machine Learning, ML) con varios parámetros de customización. Estos modelos se eligen por su capacidad de manejar diferentes tipos de datos y su rendimiento en tareas de clasificación, proporcionando así una evaluación completa y robusta para el problema en cuestión.

2.2.4.1 Naive Bayes

Naive Bayes es un clasificador probabilístico basado en el teorema de Bayes, asumiendo una independencia condicional entre las características dadas la clase. En este estudio, se ha utilizado la implementación proporcionada en el paquete e1071. Formalmente, el clasificador Naive Bayes asigna una clase C_k a una instancia x = (x_1, x_2, \ldots, x_n) maximizando la probabilidad posterior P(C_k|x), que se calcula como:

P(C_k|x) = \frac{P(C_k) \cdot P(x|C_k)}{P(x)}

Dado que P(x) es constante para todas las clases, se simplifica a:

P(C_k|x) \propto P(C_k) \prod_{i=1}^{n} P(x_i|C_k)

En resumen, el clasificador Naive Bayes es una opción eficiente y efectiva para problemas de clasificación, sin embargo, al asumir la independencia de las probabilidades condicionales, puede no dar resultados muy buenos, sin embargo, proporcionará una base sólida para predicciones iniciales, así como una base para la comparación entre modelos.

2.2.4.2 Adaboost

Adaboost, o Adaptive Boosting, es un algoritmo de aprendizaje automático que combina múltiples clasificadores débiles para crear un clasificador fuerte. En este estudio, se ha utilizado el paquete mboost de R con la configuración family=Binomial(), adaptado para tareas de clasificación binaria. Adaboost trabaja iterativamente, ajustando los pesos de las muestras de entrenamiento para concentrarse en aquellas que fueron clasificadas incorrectamente en iteraciones anteriores.

De esta forma, Adaboost asigna pesos a cada muestra y ajusta un clasificador débil h_t(x) en cada iteración t. La salida final del clasificador se obtiene mediante una combinación ponderada de estos clasificadores débiles:

H(x) = \text{sign}\left( \sum_{t=1}^{T} \alpha_t h_t(x) \right)

donde \alpha_t es el peso asignado al clasificador h_t(x) basado en su precisión. En cada iteración, el peso de las muestras mal clasificadas se incrementa, enfocando así el entrenamiento subsecuente en estas muestras.

El principal hiperparámetro ajustado en este estudio ha sido el número de iteraciones T. Este parámetro determina cuántos clasificadores débiles se combinarán para formar el clasificador final. Aumentar el número de iteraciones puede mejorar la capacidad del modelo para corregir errores y, por ende, afinar la frontera de decisión, pero también aumenta el riesgo de sobreajuste si el modelo se vuelve demasiado complejo.

2.2.4.3 Árboles de Decisión

Un árbol de decisión divide iterativamente el conjunto de datos en subconjuntos más pequeños utilizando una serie de reglas de decisión basadas en las características del conjunto de datos. Cada división se realiza en un nodo del árbol y produce ramas que representan los posibles valores de las características, conduciendo a hojas que representan las predicciones finales. El objetivo es crear un modelo que prediga el valor de la variable objetivo aprendiendo reglas de decisión simples inferidas de los datos. En este estudio, se ha implementado este algoritmo utilizando el paquete rpart de R.

Los principales hiperparámetros ajustados en este estudio han sido:

• cp_values (complexity parameter): Este parámetro controla la poda del árbol de decisión. La poda es un proceso que elimina ramas del árbol que tienen poca importancia, con el fin de reducir su complejidad y prevenir el sobreajuste. Un valor más alto de cp resulta en una mayor poda, lo que produce un árbol más simple. Ajustar este parámetro afecta directamente la profundidad del árbol y su capacidad para capturar patrones en los datos. Seleccionar un cp óptimo es esencial para equilibrar el sesgo (si cp es muy alto) y la varianza (si cp es muy bajo).

• minsplit_values: Este parámetro define el número mínimo de observaciones que un nodo debe tener para ser considerado para una división. Valores más altos de minsplit resultan en árboles más pequeños y menos detallados, ya que se requieren más observaciones para realizar una división. Por otro lado, valores más bajos permiten una mayor cantidad de divisiones, resultando en árboles más complejos. Ajustar minsplit permite controlar la granularidad de las divisiones en el árbol, asegurando que el modelo capture los patrones relevantes sin sobreajustarse a los datos de entrenamiento.

2.2.4.4 Random Forest

Random Forest opera creando una colección de árboles de decisión, cada uno entrenado con un subconjunto aleatorio del conjunto de datos y un subconjunto aleatorio de características. Este enfoque introduce dos tipos de aleatoriedad: aleatoriedad en las muestras y aleatoriedad en las características, lo que ayuda a reducir la varianza del modelo y mitigar el sobreajuste.En este estudio, se ha utilizado el paquete randomForest de R.

Los principales hiperparámetros ajustados en este estudio han sido:

• ntree_values: Este parámetro define el número de árboles en el bosque. Un mayor número de árboles generalmente mejora la precisión del modelo, ya que la predicción final se basa en la agregación de un mayor número de árboles, reduciendo así la varianza. Sin embargo, un número muy alto de árboles puede aumentar el tiempo de entrenamiento sin proporcionar mejoras significativas en el rendimiento. Ajustar el valor de ntree permite encontrar un equilibrio entre la precisión y la eficiencia computacional del modelo.

• mtry_values: Este parámetro controla el número de características que se consideran para cada división en los árboles de decisión. Un valor más bajo de mtry puede reducir el riesgo de sobreajuste, ya que introduce más variabilidad entre los árboles individuales. Un valor más alto permite que cada árbol considere más características en cada división, lo que puede mejorar la precisión del modelo si las características adicionales son informativas. Ajustar mtry permite optimizar la capacidad del modelo para capturar las relaciones relevantes en los datos sin sobreajustarse.

En resumen, Random Forest es un potente algoritmo de clasificación que combina las predicciones de múltiples árboles de decisión para mejorar la precisión y la robustez del modelo. La cuidadosa selección de los parámetros ntree_values y mtry_values asegura un modelo equilibrado que generaliza bien en nuevos datos, proporcionando predicciones fiables del estado PCR de los pacientes.

2.2.4.5 K-Nearest Neighbours

K-Nearest Neighbours (KNN) clasifica una muestra basada en las clases de sus k vecinos más cercanos en el espacio de características. En este estudio, se ha utilizado el paquete class de R para implementar KNN.

KNN clasifica una muestra x basándose en la mayoría de las clases de sus k vecinos más cercanos, identificados utilizando una métrica de distancia. La distancia entre las muestras se calcula comúnmente usando la distancia euclidiana, pero se puede generalizar utilizando la distancia de Minkowski, definida como:

d(x_i, x_j) = \left( \sum_{m=1}^{M} |x_{im} - x_{jm}|^l \right)^{1/l}

donde l es el parámetro de distancia de Minkowski. Cuando l = 2, se corresponde con la distancia euclidiana, y cuando l = 1, se convierte en la distancia de Manhattan.

Los principales hiperparámetros ajustados en este estudio son:

• k (número de vecinos más cercanos): Este parámetro define el número de vecinos que se considerarán para determinar la clase de una nueva muestra. Valores pequeños de k hacen que el modelo sea más sensible al ruido en los datos, mientras que valores grandes pueden hacer que el modelo sea demasiado general y pierda precisión en detalles locales.

• l (parámetro de distancia de Minkowski): Este parámetro ajusta la métrica de distancia utilizada para identificar los vecinos más cercanos. Aunque el valor predeterminado es l = 2 (distancia euclidiana), ajustar l permite al modelo adaptarse mejor a la geometría del espacio de características. Por ejemplo, un l menor puede capturar mejor las relaciones lineales, mientras que un l mayor puede ser más adecuado para relaciones más complejas y no lineales.

Además, el parámetro probes = TRUE se utiliza para devolver las probabilidades de las clases predichas, proporcionando una medida de confianza en las predicciones del modelo.

En resumen, KNN es un algoritmo intuitivo y efectivo para la clasificación, especialmente adecuado para conjuntos de datos donde la estructura local de las muestras es crucial para la precisión de las predicciones. Ajustar los parámetros k y l permite optimizar el modelo para capturar las relaciones subyacentes en los datos, asegurando predicciones precisas y confiables del estado PCR de los pacientes.

2.2.4.6 Support Vector Machine

La Máquina de Vectores de Soporte, o, Support Vector Machine (SVM) es una técnica de clasificación que busca encontrar el hiperplano que mejor separa las clases en el espacio de características. En este estudio, se ha utilizado el paquete e1071 de R para implementar SVM.

SVM funciona proyectando los datos a un espacio de mayor dimensión donde se puede encontrar un hiperplano que maximiza el margen de separación entre las clases. El hiperplano se define por:

\mathbf{w} \cdot \mathbf{x} + b = 0

donde \mathbf{w} es el vector de pesos y b es el sesgo. El objetivo es encontrar \mathbf{w} y b que maximicen el margen, sujetando a las restricciones de clasificación correcta de las muestras. Formalmente, los principales hiperparámetros ajustados en este estudio han sido:

• Kernel: Este parámetro define la función de transformación que se utiliza para proyectar los datos en un espacio de mayor dimensión. Los kernels disponibles incluyen:

  • Lineal: K(\mathbf{x_i}, \mathbf{x_j}) = \mathbf{x_i} \cdot \mathbf{x_j}
  • Sigmoid: $ K(, ) = ( + r) $
  • Polinomial: $ K(, ) = ( + r)^d $
  • RBF (Radial Basis Function): $ K(, ) = (-|| - ||^2) $

  La elección del kernel afecta directamente la capacidad del SVM para capturar relaciones no lineales en los datos. Por ejemplo, un kernel lineal es adecuado para datos linealmente separables, mientras que un kernel RBF es más flexible y puede modelar relaciones complejas.

• Coste (C): Este parámetro controla el margen de tolerancia a errores en la clasificación. Un valor alto de C da como resultado un margen más estrecho y menos tolerancia a errores de clasificación, lo que puede llevar a un modelo que sobreajuste los datos de entrenamiento. Por otro lado, un valor bajo de C permite un margen más amplio y más errores de clasificación, lo que puede resultar en un modelo más general.

2.2.4.7 Artificial Neural Networks

Las Redes Neuronales Artificiales (ANN) son modelos de aprendizaje automático inspirados en la estructura y funcionamiento del cerebro humano. En este estudio, se ha utilizado la función nnet del paquete nnet de R para implementar un modelo de red neuronal artificial.

Las redes neuronales consisten en una capa de entrada, una o más capas ocultas, y una capa de salida. Cada neurona en estas capas realiza una operación lineal seguida de una función de activación no lineal. En este estudio, la arquitectura del modelo se definirá principalmente por el número de neuronas en la capa oculta (size) y el método de regularización (weight-decay).

• Tamaño de la Capa Oculta (size): Este parámetro especifica el número de neuronas en la capa oculta de la red neuronal. Un mayor número de neuronas (size) puede capturar patrones más complejos en los datos, pero también aumenta el riesgo de sobreajuste si el modelo se vuelve demasiado complejo.

• Parámetro de Regularización (decay): El parámetro weight-decay controla la magnitud de la penalización aplicada a los pesos del modelo durante el entrenamiento para prevenir el sobreajuste. Este método de regularización agrega un término a la función de pérdida que penaliza grandes valores de los pesos, manteniendo los pesos más pequeños y el modelo más sencillo. Ajustar el decay afecta la capacidad del modelo para generalizar, ya que valores muy altos pueden resultar en un modelo subajustado, mientras que valores muy bajos pueden no ser suficientes para evitar el sobreajuste.

3 Resultados

3.1 Naive Bayes

3.1.1 Aparente

evaluate_aparent_performance_model_naiveBayes <- function(data, target_var, model_func, vars = NULL, threshold, seed) {
  if (is.null(vars)) {
    vars <- setdiff(names(data), target_var)
  }
  set.seed(seed = seed)
  data[[target_var]] <- factor(data[[target_var]], levels = c(0, 1))
  formula <- as.formula(paste(target_var, "~", paste(vars, collapse = "+")))
  
  model <- model_func(formula, data)
  predictions <- predict(model, newdata = data, type = "raw")[, 2]
  
  actual_classes <- data[[target_var]]
  predicted_classes <- ifelse(predictions > threshold, 1, 0)
  confusion_matrix <- table(predicted_classes, actual_classes)
  
  tp <- ifelse("1" %in% rownames(confusion_matrix) && "1" %in% colnames(confusion_matrix), confusion_matrix["1", "1"], 0)
  tn <- ifelse("0" %in% rownames(confusion_matrix) && "0" %in% colnames(confusion_matrix), confusion_matrix["0", "0"], 0)
  fp <- ifelse("1" %in% rownames(confusion_matrix) && "0" %in% colnames(confusion_matrix), confusion_matrix["1", "0"], 0)
  fn <- ifelse("0" %in% rownames(confusion_matrix) && "1" %in% colnames(confusion_matrix), confusion_matrix["0", "1"], 0)
  
  accuracy <- (tp + tn) / (tp + tn + fp + fn)
  precision <- ifelse(tp + fp > 0, tp / (tp + fp), 0)
  recall <- ifelse(tp + fn > 0, tp / (tp + fn), 0)
  f1_score <- ifelse(precision + recall > 0, 2 * (precision * recall) / (precision + recall), 0)
  roc_obj <- pROC::roc(actual_classes, predictions)
  
  list(confusion_matrix = confusion_matrix, accuracy = accuracy, precision = precision, 
       recall = recall, f1_score = f1_score, roc_curve = roc_obj)
}
naiveBayes_model <- function(formula, data) {
  naiveBayes(formula = formula, data = data)
}
resultadosAparentesNaiveBayes <- evaluate_aparent_performance_model_naiveBayes(data = data_factor, 
                                                                               target_var = "PCR",
                                                                              model_func = naiveBayes_model,
                                                                              vars = ".",
                                                                              threshold = .35,
                                                                              seed = 90)

load("RData_Files_Algorithms/NaiveBayes_Aparente.RData")

3.1.2 Todas las Variables

double_cross_validation_naiveBayes <- function(data, target_name, outer, inner, variables, threshold, seed) {
  set.seed(seed)
  performance_metrics <- data.frame(Fold = integer(), TP = integer(), TN = integer(),
                                    FP = integer(), FN = integer(), AUC = numeric(), stringsAsFactors = FALSE)
  
  data[[target_name]] <- as.factor(data[[target_name]])
  outer_folds <- createFolds(data[[target_name]], k = outer, list = TRUE, returnTrain = FALSE)
  
  for (outer_index in seq_along(outer_folds)) {
    outer_train_data <- data[-outer_folds[[outer_index]], ]
    outer_test_data <- data[outer_folds[[outer_index]], ]
    
    cat(sprintf("Outer Fold %d\n", outer_index))
    cat("Training Data - Class 0: ", sum(outer_train_data[[target_name]] == "0"), "\n")
    cat("Training Data - Class 1: ", sum(outer_train_data[[target_name]] == "1"), "\n")
    cat("Testing Data - Class 0: ", sum(outer_test_data[[target_name]] == "0"), "\n")
    cat("Testing Data - Class 1: ", sum(outer_test_data[[target_name]] == "1"), "\n")
    
    if (sum(outer_train_data[[target_name]] == "0") == 0 || sum(outer_train_data[[target_name]] == "1") == 0) {
      cat("Skipping outer fold due to lack of class diversity\n")
      next
    }
    
    inner_folds <- createFolds(outer_train_data[[target_name]], k = inner, list = TRUE, returnTrain = FALSE)
    best_model <- NULL
    best_auc <- 0

    formula <- reformulate(variables, target_name)

    for (inner_index in seq_along(inner_folds)) {
      inner_train_data <- outer_train_data[-inner_folds[[inner_index]], ]
      inner_test_data <- outer_train_data[inner_folds[[inner_index]], ]

      if (sum(inner_train_data[[target_name]] == "0") == 0 || sum(inner_train_data[[target_name]] == "1") == 0) {
        cat("Skipping inner fold due to lack of class diversity\n")
        next
      }

      inner_auc <- numeric()

      result <- tryCatch({
        model <- naiveBayes(formula, data = inner_train_data)
        predictions <- predict(model, inner_test_data, type = "raw")[, 2]
        roc_curve <- roc(inner_test_data[[target_name]], predictions)
        inner_auc[inner_index] <- roc_curve$auc
      }, error = function(e) {
        cat("Error in inner fold ", inner_index, "\n")
        cat("Error message: ", e$message, "\n")
        inner_auc[inner_index] <- 0
      })

      if (mean(inner_auc, na.rm = TRUE) > best_auc) {
        best_auc <- mean(inner_auc, na.rm = TRUE)
        best_model <- model
      }
    }

    if (is.null(best_model)) {
      cat("No valid model found for outer fold ", outer_index, "\n")
      next
    }
    best_model <- naiveBayes(formula, data = outer_train_data)
    predictions <- predict(best_model, outer_test_data, type = "raw")[, 2]
    pred_class <- ifelse(predictions > threshold, "1", "0")
    confusion <- table(Actual = outer_test_data[[target_name]], Predicted = pred_class)

    if (nrow(confusion) < 2 || ncol(confusion) < 2) {
      TP <- ifelse(nrow(confusion) >= 2 && ncol(confusion) >= 2, confusion[2, 2], 0)
      TN <- ifelse(nrow(confusion) >= 1 && ncol(confusion) >= 1, confusion[1, 1], 0)
      FP <- ifelse(nrow(confusion) >= 1 && ncol(confusion) >= 2, confusion[1, 2], 0)
      FN <- ifelse(nrow(confusion) >= 2 && ncol(confusion) >= 1, confusion[2, 1], 0)
    } else {
      TP <- confusion[2, 2]
      TN <- confusion[1, 1]
      FP <- confusion[1, 2]
      FN <- confusion[2, 1]
    }
    auc <- roc(outer_test_data[[target_name]], predictions)$auc

    metrics <- setNames(c(outer_index, TP, TN, FP, FN, auc), names(performance_metrics))
    performance_metrics <- rbind(performance_metrics, as.data.frame(t(metrics)))
  }

  return(performance_metrics)
}
evaluate_with_seeds_NaiveBayes_Asociacion_Todas <- 
  function(evaluation_function, data, target_var, vars, threshold, seeds, outer, inner) {
  results_list <- list()

  for (seed in seeds) {
    set.seed(seed)
    cat("Evaluating with seed:", seed, "\n")

    results <- evaluation_function(data = data, 
                                   target_name = target_var, 
                                   variables = vars, 
                                   threshold = threshold, 
                                   outer = outer, 
                                   inner = inner, 
                                   seed = seed)
    results_list[[paste0("seed_", seed)]] <- results
  }

  return(results_list)
}

performance_naiveBayes_todasVariables <- evaluate_with_seeds_NaiveBayes_Asociacion_Todas(
  evaluation_function = double_cross_validation_naiveBayes,
  data = data_factor,
  outer = 5,
  inner = 2,
  target_var = "PCR",
  vars = todas_variables_factor,
  threshold = .35,
  seeds = seeds)

load("RData_Files_Algorithms/NaiveBayes_TodasVariables.RData")

3.1.3 Filtrado (Asociación)

performance_naiveBayes_Asociacion <- evaluate_with_seeds_NaiveBayes_Asociacion_Todas(
  evaluation_function = double_cross_validation_naiveBayes,
  data = data_factor,
  outer = 5,
  inner = 2,
  target_var = "PCR",
  vars = todas_variables_factor,
  threshold = .35,
  seeds = seeds)

load("RData_Files_Algorithms/NaiveBayes_Asociacion.RData")

3.1.4 Wrapped (StepAUC)

double_cross_validation_naiveBayes_stepAUC <- function(data, target_name, outer, inner, variables, threshold, seed) {
  set.seed(seed)
  performance_metrics <- data.frame(Fold = integer(), TP = integer(), TN = integer(),
                                    FP = integer(), FN = integer(), AUC = numeric(),
                                    Best_Variables = character(), stringsAsFactors = FALSE)
  
  inner_fold_metrics <- data.frame(Fold = integer(), Inner_Fold = integer(), 
                                   TP = integer(), TN = integer(), FP = integer(), 
                                   FN = integer(), AUC = numeric(), stringsAsFactors = FALSE)
  
  data[[target_name]] <- as.factor(data[[target_name]])
  outer_folds <- createFolds(data[[target_name]], k = outer, list = TRUE, returnTrain = FALSE)
  
  for (outer_index in seq_along(outer_folds)) {
    outer_train_data <- data[-outer_folds[[outer_index]], ]
    outer_test_data <- data[outer_folds[[outer_index]], ]
    
    cat(sprintf("Outer Fold %d\n", outer_index))
    cat("Training Data - Class 0: ", sum(outer_train_data[[target_name]] == "0"), "\n")
    cat("Training Data - Class 1: ", sum(outer_train_data[[target_name]] == "1"), "\n")
    cat("Testing Data - Class 0: ", sum(outer_test_data[[target_name]] == "0"), "\n")
    cat("Testing Data - Class 1: ", sum(outer_test_data[[target_name]] == "1"), "\n")
    
    inner_folds <- createFolds(outer_train_data[[target_name]], k = inner, list = TRUE, returnTrain = FALSE)
    best_model <- NULL
    best_auc <- 0
    best_vars <- variables
    
    stepAUC <- function(vars, outer_train_data, target_name, inner_folds) {
      current_vars <- vars
      best_inner_auc <- 0
      best_inner_model <- NULL
      best_inner_vars <- current_vars
      improved <- TRUE
      
      while (improved && length(current_vars) > 1) {
        improved <- FALSE
        best_auc_in_step <- best_inner_auc
        
        for (var in current_vars) {
          temp_vars <- setdiff(current_vars, var)
          inner_auc <- numeric()
          
          for (inner_index in seq_along(inner_folds)) {
            inner_train_data <- outer_train_data[-inner_folds[[inner_index]], ]
            inner_test_data <- outer_train_data[inner_folds[[inner_index]], ]
            
            result <- tryCatch({
              model <- naiveBayes(as.formula(paste(target_name, "~", paste(temp_vars, collapse = "+"))),
                                  data = inner_train_data)
              
              predictions <- predict(model, newdata = inner_test_data, type = "raw")[, 2]
              roc_curve <- roc(inner_test_data[[target_name]], predictions)
              inner_auc[inner_index] <- roc_curve$auc
              
            }, error = function(e) {
              cat("Error with variables:", paste(temp_vars, collapse = ", "), "\n")
              inner_auc[inner_index] <- NA  # Set to NA or some other indicator of failure
            })
          }
          
          mean_inner_auc <- mean(inner_auc, na.rm = TRUE)
          if (!is.na(mean_inner_auc) && mean_inner_auc > best_auc_in_step) {
            best_auc_in_step <- mean_inner_auc
            best_inner_model <- model
            best_inner_vars <- temp_vars
            improved <- TRUE
          }
        }
        
        if (improved) {
          current_vars <- best_inner_vars
          best_inner_auc <- best_auc_in_step
          cat("Vars improved: ", paste(current_vars, collapse=", "), "\n")
          cat("AUC: ", paste(best_inner_auc), "\n")
        }
      }
      
      list(model = best_inner_model, auc = best_inner_auc, vars = best_inner_vars)
    }
    
    best_result <- stepAUC(variables, outer_train_data, target_name, inner_folds)
    best_model <- best_result$model
    best_auc <- best_result$auc
    best_vars <- best_result$vars

    if (is.null(best_model)) {
      cat("Error: No valid model found for outer fold ", outer_index, "\n")
      next
    }

    # Evaluar el mejor modelo en el validation split del inner loop
    inner_auc_vals <- numeric()
    for (inner_index in seq_along(inner_folds)) {
      inner_train_data <- outer_train_data[-inner_folds[[inner_index]], ]
      inner_test_data <- outer_train_data[inner_folds[[inner_index]], ]
      
      predictions <- predict(best_model, newdata = inner_test_data, type = "raw")[, 2]
      roc_curve <- roc(inner_test_data[[target_name]], predictions)
      inner_auc_vals[inner_index] <- roc_curve$auc

      pred_class <- ifelse(predictions > threshold, "1", "0")
      confusion <- table(Actual = inner_test_data[[target_name]], Predicted = pred_class)
      TP <- ifelse(nrow(confusion) >= 2 && ncol(confusion) >= 2, confusion[2, 2], 0)
      TN <- ifelse(nrow(confusion) >= 1 && ncol(confusion) >= 1, confusion[1, 1], 0)
      FP <- ifelse(nrow(confusion) >= 1 && ncol(confusion) >= 2, confusion[1, 2], 0)
      FN <- ifelse(nrow(confusion) >= 2 && ncol(confusion) >= 1, confusion[2, 1], 0)
      auc <- roc_curve$auc
      
      inner_metrics <- setNames(c(outer_index, inner_index, TP, TN, FP, FN, auc), names(inner_fold_metrics))
      inner_fold_metrics <- rbind(inner_fold_metrics, as.data.frame(t(inner_metrics)))
    }
    
    best_inner_auc <- mean(inner_auc_vals, na.rm = TRUE)
    best_model <- naiveBayes(as.formula(paste(target_name, "~", paste(best_vars, collapse = "+"))),
                         data = outer_train_data)

    predictions <- predict(best_model, outer_test_data, type = "raw")[, 2]
    pred_class <- ifelse(predictions > threshold, "1", "0")
    confusion <- table(Actual = outer_test_data[[target_name]], Predicted = pred_class)
    
    if (nrow(confusion) < 2 || ncol(confusion) < 2) {
      TP <- ifelse(nrow(confusion) >= 2 && ncol(confusion) >= 2, confusion[2, 2], 0)
      TN <- ifelse(nrow(confusion) >= 1 && ncol(confusion) >= 1, confusion[1, 1], 0)
      FP <- ifelse(nrow(confusion) >= 1 && ncol(confusion) >= 2, confusion[1, 2], 0)
      FN <- ifelse(nrow(confusion) >= 2 && ncol(confusion) >= 1, confusion[2, 1], 0)
    } else {
      TP <- confusion[2, 2]
      TN <- confusion[1, 1]
      FP <- confusion[1, 2]
      FN <- confusion[2, 1]
    }
    auc <- roc(outer_test_data[[target_name]], predictions)$auc
    
    metrics <- setNames(c(outer_index, TP, TN, FP, FN, auc, 
                          paste(best_vars, collapse = ",")), 
                        names(performance_metrics))
    performance_metrics <- rbind(performance_metrics, as.data.frame(t(metrics)))
    
    cat("Performance:", paste(performance_metrics))
    cat(sprintf("Confusion Matrix for Fold %d:\n", outer_index))
    print(confusion)
    cat(sprintf("Metrics for Fold %d:\n", outer_index))
    cat(sprintf("TP: %d, TN: %d, FP: %d, FN: %d, AUC: %f, Best Variables: %s\n",
                TP, TN, FP, FN, auc, paste(best_vars, collapse = ",")))
  }
  
  list(performance_metrics = performance_metrics, inner_fold_metrics = inner_fold_metrics)
}
evaluate_with_seeds_NaiveBayes_stepAUC <- 
  function(evaluation_function, data, target_var, vars, threshold, seeds, outer, inner) {
  results_list <- list()

  for (seed in seeds) {
    set.seed(seed)
    cat("Evaluating with seed:", seed, "\n")

    results <- evaluation_function(data = data, 
                                   target_name = target_var, 
                                   variables = vars, 
                                   threshold = threshold, 
                                   outer = outer, 
                                   inner = inner, 
                                   seed = seed)
    results_list[[paste0("seed_", seed)]] <- results
  }

  return(results_list)
}

performance_naiveBayes_stepAUC <- evaluate_with_seeds_NaiveBayes_stepAUC(
  evaluation_function = double_cross_validation_naiveBayes_stepAUC,
  data = data_factor,
  outer = 5,
  inner = 2,
  vars = todas_variables_factor,
  target_var = "PCR",
  threshold = .35,
  seeds = seeds
)

load("RData_Files_Algorithms/NaiveBayes_StepAUC.RData")

3.1.5 Embedded (Lasso)

double_cross_validation_naiveBayes_lasso <- function(data, target_name, outer, inner, variables, threshold, seed) {
  set.seed(seed)
  performance_metrics <- data.frame(Fold = integer(), TP = integer(), TN = integer(),
                                    FP = integer(), FN = integer(), AUC = numeric(),
                                    Best_Variables = character(), stringsAsFactors = FALSE)
  
  inner_fold_metrics <- data.frame(Fold = integer(), Inner_Fold = integer(), 
                                   TP = integer(), TN = integer(), FP = integer(), 
                                   FN = integer(), AUC = numeric(), stringsAsFactors = FALSE)
  
  # Convertir las variables del dataset con model.matrix
  target <- as.factor(data[[target_name]])  # Aseguramos que la variable de respuesta sea factor
  features <- data[, setdiff(names(data), target_name)]
  data_matrix <- as.data.frame(model.matrix(~ ., data=features)[,-1])  # Eliminar el intercepto
  data_matrix[[target_name]] <- target
  
  outer_folds <- createFolds(data_matrix[[target_name]], k = outer)
  
  for (outer_index in seq_along(outer_folds)) {
    outer_train_data <- data_matrix[-outer_folds[[outer_index]], ]
    outer_test_data <- data_matrix[outer_folds[[outer_index]], ]
    
    cat(sprintf("Outer Fold %d\n", outer_index))
    cat("Training Data - Class 0: ", sum(outer_train_data[[target_name]] == 0), "\n")
    cat("Training Data - Class 1: ", sum(outer_train_data[[target_name]] == 1), "\n")
    cat("Testing Data - Class 0: ", sum(outer_test_data[[target_name]] == 0), "\n")
    cat("Testing Data - Class 1: ", sum(outer_test_data[[target_name]] == 1), "\n")
    
    inner_folds <- createFolds(outer_train_data[[target_name]], k = inner)
    best_model <- NULL
    best_auc <- 0
    best_vars <- colnames(outer_train_data)[!colnames(outer_train_data) %in% target_name]
    
    select_vars_lasso <- function(data, target_name) {
      target <- as.numeric(data[[target_name]]) - 1  # Convertimos a numérico para Lasso
      features <- data[, setdiff(names(data), target_name)]
      X <- model.matrix(~ ., data=features)[,-1]  # Eliminamos el intercepto
      y <- target
      
      # Ajustamos el modelo Lasso
      lasso_model <- cv.glmnet(X, y, alpha=1)
      lasso_coef <- coef(lasso_model, s = "lambda.min")
      
      # Seleccionamos las variables no nulas
      selected_vars <- rownames(lasso_coef)[lasso_coef[, 1] != 0]
      selected_vars <- setdiff(selected_vars, "(Intercept)")
      cat("Selected variables by Lasso: ", paste(selected_vars, collapse=", "), "\n")
      return(selected_vars)
    }
    
    stepAUC <- function(vars, data, target_name, inner_folds) {
      best_inner_auc <- 0
      best_inner_model <- NULL
      best_inner_vars <- vars
      improved <- TRUE
      
      while (improved && length(vars) > 1) {
        improved <- FALSE
        best_auc_in_step <- best_inner_auc
        
        auc_vals <- numeric()
        
        for (inner_index in seq_along(inner_folds)) {
          inner_train_data <- outer_train_data[-inner_folds[[inner_index]], ]
          inner_test_data <- outer_train_data[inner_folds[[inner_index]], ]
          
          selected_vars <- select_vars_lasso(inner_train_data, target_name)
          formula <- as.formula(paste(target_name, "~", paste(selected_vars, collapse = "+")))
          cat("Testing with formula: ", deparse(formula), "\n")
          
          result <- tryCatch({
            model <- naiveBayes(formula, data = inner_train_data)
            
            predictions <- predict(model, inner_test_data, type = "raw")[,2]
            roc_curve <- roc(inner_test_data[[target_name]], predictions)
            auc_vals[inner_index] <- roc_curve$auc
            cat("AUC for inner fold ", inner_index, ": ", roc_curve$auc, "\n")
            
          }, error = function(e) {
            cat("Error with variables: ", paste(selected_vars, collapse = ", "), "\n")
            cat("Error message: ", e$message, "\n")
            auc_vals[inner_index] <- NA  # Set to NA or some other indicator of failure
          })
        }
        
        mean_auc <- mean(auc_vals, na.rm = TRUE)
        cat("Mean AUC: ", mean_auc, "\n")
        if (!is.na(mean_auc) && mean_auc > best_auc_in_step) {
          best_auc_in_step <- mean_auc
          best_inner_model <- model
          best_inner_vars <- selected_vars
          improved <- TRUE
        }
        
        if (improved) {
          vars <- best_inner_vars
          best_inner_auc <- best_auc_in_step
          cat("Vars improved: ", paste(vars, collapse=", "), "\n")
          cat("AUC: ", best_inner_auc, "\n")
        }
      }
      
      list(model = best_inner_model, auc = best_inner_auc, vars = best_inner_vars)
    }
    
    best_result <- stepAUC(best_vars, outer_train_data, target_name, inner_folds)
    best_model <- best_result$model
    best_auc <- best_result$auc
    best_vars <- best_result$vars
    
    if (is.null(best_model)) {
      cat("Error: No valid model found for outer fold ", outer_index, "\n")
      next
    }

    # Evaluar el mejor modelo en el validation split del inner loop
    inner_auc_vals <- numeric()
    for (inner_index in seq_along(inner_folds)) {
      inner_train_data <- outer_train_data[-inner_folds[[inner_index]], ]
      inner_test_data <- outer_train_data[inner_folds[[inner_index]], ]
      
      predictions <- predict(best_model, inner_test_data, type = "raw")[,2]
      roc_curve <- roc(inner_test_data[[target_name]], predictions)
      inner_auc_vals[inner_index] <- roc_curve$auc
      cat("Inner fold ", inner_index, " AUC: ", roc_curve$auc, "\n")

      pred_class <- ifelse(predictions > threshold, 1, 0)
      confusion <- table(Actual = inner_test_data[[target_name]], Predicted = pred_class)
      TP <- ifelse(nrow(confusion) >= 2 && ncol(confusion) >= 2, confusion[2, 2], 0)
      TN <- ifelse(nrow(confusion) >= 1 && ncol(confusion) >= 1, confusion[1, 1], 0)
      FP <- ifelse(nrow(confusion) >= 1 && ncol(confusion) >= 2, confusion[1, 2], 0)
      FN <- ifelse(nrow(confusion) >= 2 && ncol(confusion) >= 1, confusion[2, 1], 0)
      auc <- roc_curve$auc
      
      inner_metrics <- setNames(c(outer_index, inner_index, TP, TN, FP, FN, auc), names(inner_fold_metrics))
      inner_fold_metrics <- rbind(inner_fold_metrics, as.data.frame(t(inner_metrics)))
    }
    
    best_inner_auc <- mean(inner_auc_vals, na.rm = TRUE)
    cat("Best inner AUC: ", best_inner_auc, "\n")
    best_model <- naiveBayes(as.formula(paste(target_name, "~", paste(best_vars, collapse = "+"))), data = outer_train_data)
    predictions <- predict(best_model, outer_test_data, type = "raw")[,2]
    pred_class <- ifelse(predictions > threshold, 1, 0)  # Asignar clases basadas en umbral
    confusion <- table(Actual = outer_test_data[[target_name]], Predicted = pred_class)
    
    if (nrow(confusion) < 2 || ncol(confusion) < 2) {
      TP <- ifelse(nrow(confusion) >= 2 && ncol(confusion) >= 2, confusion[2, 2], 0)
      TN <- ifelse(nrow(confusion) >= 1 && ncol(confusion) >= 1, confusion[1, 1], 0)
      FP <- ifelse(nrow(confusion) >= 1 && ncol(confusion) >= 2, confusion[1, 2], 0)
      FN <- ifelse(nrow(confusion) >= 2 && ncol(confusion) >= 1, confusion[2, 1], 0)
    } else {
      TP <- confusion[2, 2]
      TN <- confusion[1, 1]
      FP <- confusion[1, 2]
      FN <- confusion[2, 1]
    }
    auc <- roc(outer_test_data[[target_name]], predictions)$auc
    
    metrics <- setNames(c(outer_index, TP, TN, FP, FN, auc, 
                          paste(best_vars, collapse = ",")), 
                        names(performance_metrics))
    performance_metrics <- rbind(performance_metrics, as.data.frame(t(metrics)))
    
    cat("Performance:\n", paste(performance_metrics), "\n")
    cat(sprintf("Confusion Matrix for Fold %d:\n", outer_index))
    print(confusion)
    cat(sprintf("Metrics for Fold %d:\n", outer_index))
    cat(sprintf("TP: %d, TN: %d, FP: %d, FN: %d, AUC: %f, Best Variables: %s\n",
                TP, TN, FP, FN, auc, paste(best_vars, collapse = ",")))
  }
  
  list(performance_metrics = performance_metrics, inner_fold_metrics = inner_fold_metrics)
}
evaluate_with_seeds_NaiveBayes_Lasso <- 
  function(evaluation_function, data, target_var, vars, threshold, seeds, outer, inner) {
  results_list <- list()

  for (seed in seeds) {
    set.seed(seed)
    cat("Evaluating with seed:", seed, "\n")

    results <- evaluation_function(data = data, 
                                   target_name = target_var, 
                                   variables = vars, 
                                   threshold = threshold, 
                                   outer = outer, 
                                   inner = inner, 
                                   seed = seed)
    results_list[[paste0("seed_", seed)]] <- results
  }

  return(results_list)
  }

performance_naiveBayes_lasso <- evaluate_with_seeds_NaiveBayes_Lasso(
  evaluation_function = double_cross_validation_naiveBayes_lasso,
  data = data_factor,
  outer = 5,
  inner = 2,
  vars = todas_variables_factor,
  target_var = "PCR",
  threshold = .35,
  seeds = seeds
)

load("RData_Files_Algorithms/NaiveBayes_Lasso.RData")

3.2 Adaboost

3.2.1 Aparente

evaluate_aparent_performance_model_adaboost <- function(data, target_var, model_func, vars = NULL, iterations = 100, threshold) {
  if (is.null(vars)) {
    vars <- setdiff(names(data), target_var)
  }
  
  data[[target_var]] <- factor(data[[target_var]], levels = c(0, 1))
  formula <- as.formula(paste(target_var, "~", paste(vars, collapse = "+")))
  
  model <- model_func(formula, data, iterations)
  predictions <- predict(model, newdata = data, type = "response")
  
  actual_classes <- data[[target_var]]
  predicted_classes <- ifelse(predictions > threshold, 1, 0)
  confusion_matrix <- table(Predicted = predicted_classes, Actual = actual_classes)
  
  tp <- ifelse("1" %in% rownames(confusion_matrix) && "1" %in% colnames(confusion_matrix), confusion_matrix["1", "1"], 0)
  tn <- ifelse("0" %in% rownames(confusion_matrix) && "0" %in% colnames(confusion_matrix), confusion_matrix["0", "0"], 0)
  fp <- ifelse("1" %in% rownames(confusion_matrix) && "0" %in% colnames(confusion_matrix), confusion_matrix["1", "0"], 0)
  fn <- ifelse("0" %in% rownames(confusion_matrix) && "1" %in% colnames(confusion_matrix), confusion_matrix["0", "1"], 0)
  
  accuracy <- (tp + tn) / (tp + tn + fp + fn)
  precision <- ifelse(tp + fp > 0, tp / (tp + fp), 0)
  recall <- ifelse(tp + fn > 0, tp / (tp + fn), 0)
  f1_score <- ifelse(precision + recall > 0, 2 * (precision * recall) / (precision + recall), 0)
  roc_obj <- pROC::roc(actual_classes, predictions)
  
  list(confusion_matrix = confusion_matrix, accuracy = accuracy, precision = precision, 
       recall = recall, f1_score = f1_score, roc_curve = roc_obj)
}

adaboost_model <- function(formula, data, iterations) {
  model <- glmboost(formula, data = data, family = Binomial(), control = boost_control(mstop = iterations))
  model
}

resultadosAparentesAdaBoost <- evaluate_aparent_performance_model_adaboost(data = data_numeric, target_var = "PCR",
                                                                           model_func = adaboost_model,
                                                                           vars = ".",
                                                                           iterations = 100,
                                                                           threshold = 0.35)

load("RData_Files_Algorithms/AdaBoost_Aparente.RData")

3.2.2 Todas las Variables

double_cross_validation_adaboost <- function(data, target_name, outer, inner, iterations, variables, threshold, seed) {
  set.seed(seed)
  performance_metrics <- data.frame(Fold = integer(), TP = integer(), TN = integer(),
                                    FP = integer(), FN = integer(), AUC = numeric(), 
                                    Best_Iteration = integer(), stringsAsFactors = FALSE)
  
  data[[target_name]] <- as.factor(data[[target_name]])
  outer_folds <- createFolds(data[[target_name]], k = outer, list = TRUE, returnTrain = FALSE)
  
  for (outer_index in seq_along(outer_folds)) {
    outer_train_data <- data[-outer_folds[[outer_index]], ]
    outer_test_data <- data[outer_folds[[outer_index]], ]
    
    cat(sprintf("Outer Fold %d\n", outer_index))
    cat("Training Data - Class 0: ", sum(outer_train_data[[target_name]] == "0"), "\n")
    cat("Training Data - Class 1: ", sum(outer_train_data[[target_name]] == "1"), "\n")
    cat("Testing Data - Class 0: ", sum(outer_test_data[[target_name]] == "0"), "\n")
    cat("Testing Data - Class 1: ", sum(outer_test_data[[target_name]] == "1"), "\n")
    
    if (sum(outer_train_data[[target_name]] == "0") == 0 || sum(outer_train_data[[target_name]] == "1") == 0) {
      cat("Skipping outer fold due to lack of class diversity\n")
      next
    }
    
    inner_folds <- createFolds(outer_train_data[[target_name]], k = inner, list = TRUE, returnTrain = FALSE)
    best_model <- NULL
    best_auc <- 0
    best_iteration <- NULL

    grid <- expand.grid(iteration = iterations)
    formula <- reformulate(variables, target_name)

    for (params in 1:nrow(grid)) {
      inner_auc <- numeric()

      for (inner_index in seq_along(inner_folds)) {
        inner_train_data <- outer_train_data[-inner_folds[[inner_index]], ]
        inner_test_data <- outer_train_data[inner_folds[[inner_index]], ]

        if (sum(inner_train_data[[target_name]] == "0") == 0 || sum(inner_train_data[[target_name]] == "1") == 0) {
          cat("Skipping inner fold due to lack of class diversity\n")
          next
        }

        result <- tryCatch({
          model <- blackboost(formula, data = inner_train_data, family = Binomial(), control = boost_control(mstop = grid$iteration[params]))
          
          predictions <- predict(model, inner_test_data, type = "response")
          roc_curve <- roc(inner_test_data[[target_name]], predictions)
          inner_auc[inner_index] <- roc_curve$auc
        }, error = function(e) {
          cat("Error with parameters: iteration =", grid$iteration[params], "\n")
          cat("Error message: ", e$message, "\n")
          inner_auc[inner_index] <- 0
        })
      }

      if (mean(inner_auc, na.rm = TRUE) > best_auc) {
        best_auc <- mean(inner_auc, na.rm = TRUE)
        best_model <- model
        best_iteration <- grid$iteration[params]
      }
    }

    if (is.null(best_model)) {
      cat("No valid model found for outer fold ", outer_index, "\n")
      next
    }
    # Reentrenar el mejor modelo con todos los datos de outer_train_data
    best_model <- blackboost(formula, data = outer_train_data, family = Binomial(), control = boost_control(mstop = best_iteration))
    
    predictions <- predict(best_model, outer_test_data, type = "response")
    pred_class <- ifelse(predictions > threshold, "1", "0")
    confusion <- table(Actual = outer_test_data[[target_name]], Predicted = pred_class)

    if (nrow(confusion) < 2 || ncol(confusion) < 2) {
      TP <- ifelse(nrow(confusion) >= 2 && ncol(confusion) >= 2, confusion[2, 2], 0)
      TN <- ifelse(nrow(confusion) >= 1 && ncol(confusion) >= 1, confusion[1, 1], 0)
      FP <- ifelse(nrow(confusion) >= 1 && ncol(confusion) >= 2, confusion[1, 2], 0)
      FN <- ifelse(nrow(confusion) >= 2 && ncol(confusion) >= 1, confusion[2, 1], 0)
    } else {
      TP <- confusion[2, 2]
      TN <- confusion[1, 1]
      FP <- confusion[1, 2]
      FN <- confusion[2, 1]
    }
    auc <- roc(outer_test_data[[target_name]], predictions)$auc

    metrics <- setNames(c(outer_index, TP, TN, FP, FN, auc, best_iteration), names(performance_metrics))
    performance_metrics <- rbind(performance_metrics, as.data.frame(t(metrics)))
  }

  return(performance_metrics)
}

evaluate_with_seeds_AdaBoost_Asociacion_Todas <- 
  function(evaluation_function, data, target_var, vars, threshold, seeds, outer, inner, iterations) {
  results_list <- list()

  for (seed in seeds) {
    set.seed(seed)
    cat("Evaluating with seed:", seed, "\n")

    results <- evaluation_function(data = data, 
                                   target_name = target_var, 
                                   variables = vars, 
                                   threshold = threshold, 
                                   outer = outer, 
                                   inner = inner,
                                   iterations, 
                                   seed = seed)
    results_list[[paste0("seed_", seed)]] <- results
  }

  return(results_list)
}

performance_adaboost_todasvariables <- evaluate_with_seeds_AdaBoost_Asociacion_Todas(
  evaluation_function = double_cross_validation_adaboost,
  data = data_factor,
  outer = 5,
  inner = 2,
  target_var = "PCR",
  vars = todas_variables_factor,
  threshold = .35,
  seeds = seeds,
  iterations = c(50, 100, 150, 200)
)

load("RData_Files_Algorithms/AdaBoost_TodasVariables.RData")

3.2.3 Filtrado (Asociación)

performance_adaboost_Asociacion <- evaluate_with_seeds_AdaBoost_Asociacion_Todas(
  evaluation_function = double_cross_validation_adaboost,
  data = data_factor,
  outer = 5,
  inner = 2,
  target_var = "PCR",
  vars = asociacion,
  threshold = .35,
  seeds = seeds,
  iterations = c(50, 100, 150, 200)
)

load("RData_Files_Algorithms/AdaBoost_Asociacion.RData")

3.2.4 Wrapped (StepAUC)

double_cross_validation_adaboost_stepAUC <- function(data, target_name, outer, inner, iterations, variables, threshold, seed) {
  set.seed(seed)
  performance_metrics <- data.frame(Fold = integer(), TP = integer(), TN = integer(),
                                    FP = integer(), FN = integer(), AUC = numeric(),
                                    Best_Iterations = integer(), Best_Variables = character(), stringsAsFactors = FALSE)
  
  inner_fold_metrics <- data.frame(Fold = integer(), Inner_Fold = integer(), 
                                   TP = integer(), TN = integer(), FP = integer(), 
                                   FN = integer(), AUC = numeric(), stringsAsFactors = FALSE)
  
  data[[target_name]] <- as.factor(data[[target_name]])
  outer_folds <- createFolds(data[[target_name]], k = outer, list = TRUE, returnTrain = FALSE)
  
  for (outer_index in seq_along(outer_folds)) {
    outer_train_data <- data[-outer_folds[[outer_index]], ]
    outer_test_data <- data[outer_folds[[outer_index]], ]
    
    cat(sprintf("Outer Fold %d\n", outer_index))
    cat("Training Data - Class 0: ", sum(outer_train_data[[target_name]] == "0"), "\n")
    cat("Training Data - Class 1: ", sum(outer_train_data[[target_name]] == "1"), "\n")
    cat("Testing Data - Class 0: ", sum(outer_test_data[[target_name]] == "0"), "\n")
    cat("Testing Data - Class 1: ", sum(outer_test_data[[target_name]] == "1"), "\n")
    
    inner_folds <- createFolds(outer_train_data[[target_name]], k = inner, list = TRUE, returnTrain = FALSE)
    best_model <- NULL
    best_auc <- 0
    best_iterations <- NULL
    best_vars <- variables
    
    stepAUC <- function(vars, outer_train_data, target_name, iterations, inner_folds) {
      current_vars <- vars
      best_inner_auc <- 0
      best_inner_model <- NULL
      best_inner_iterations <- NULL
      best_inner_vars <- current_vars
      improved <- TRUE
      
      while (improved && length(current_vars) > 1) {
        improved <- FALSE
        best_auc_in_step <- best_inner_auc
        
        for (var in current_vars) {
          temp_vars <- setdiff(current_vars, var)
          inner_auc <- numeric()
          
          for (iter in iterations) {
            auc_vals <- numeric()
            
            for (inner_index in seq_along(inner_folds)) {
              inner_train_data <- outer_train_data[-inner_folds[[inner_index]], ]
              inner_test_data <- outer_train_data[inner_folds[[inner_index]], ]
              
              result <- tryCatch({
                model <- glmboost(as.formula(paste(target_name, "~", paste(temp_vars, collapse = "+"))),
                                  data = inner_train_data, family = Binomial(),
                                  control = boost_control(mstop = iter))
                
                predictions <- predict(model, newdata = inner_test_data, type = "response")
                roc_curve <- roc(inner_test_data[[target_name]], predictions)
                auc_vals[inner_index] <- roc_curve$auc
                
              }, error = function(e) {
                cat("Error with variables:", paste(temp_vars, collapse = ", "), "\n")
                auc_vals[inner_index] <- 0
              })
            }
            
            if (mean(auc_vals, na.rm = TRUE) > best_auc_in_step) {
              best_auc_in_step <- mean(auc_vals, na.rm = TRUE)
              best_inner_model <- model
              best_inner_iterations <- iter
              best_inner_vars <- temp_vars
              improved <- TRUE
            }
          }
        }
        
        if (improved) {
          current_vars <- best_inner_vars
          best_inner_auc <- best_auc_in_step
          cat("Vars improved: ", paste(current_vars, collapse=", "), "\n")
          cat("AUC: ", paste(best_inner_auc), "\n")
        }
      }
      
      list(model = best_inner_model, auc = best_inner_auc, 
           iterations = best_inner_iterations, 
           vars = best_inner_vars)
    }
    
    best_result <- stepAUC(variables, outer_train_data, target_name, iterations, inner_folds)
    best_model <- best_result$model
    best_auc <- best_result$auc
    best_iterations <- best_result$iterations
    best_vars <- best_result$vars

    # Evaluar el mejor modelo en el validation split del inner loop
    inner_auc_vals <- numeric()
    for (inner_index in seq_along(inner_folds)) {
      inner_train_data <- outer_train_data[-inner_folds[[inner_index]], ]
      inner_test_data <- outer_train_data[inner_folds[[inner_index]], ]
      
      predictions <- predict(best_model, newdata = inner_test_data, type = "response")
      roc_curve <- roc(inner_test_data[[target_name]], predictions)
      inner_auc_vals[inner_index] <- roc_curve$auc

      pred_class <- ifelse(predictions > threshold, "1", "0")
      confusion <- table(Actual = inner_test_data[[target_name]], Predicted = pred_class)
      TP <- ifelse(nrow(confusion) >= 2 && ncol(confusion) >= 2, confusion[2, 2], 0)
      TN <- ifelse(nrow(confusion) >= 1 && ncol(confusion) >= 1, confusion[1, 1], 0)
      FP <- ifelse(nrow(confusion) >= 1 && ncol(confusion) >= 2, confusion[1, 2], 0)
      FN <- ifelse(nrow(confusion) >= 2 && ncol(confusion) >= 1, confusion[2, 1], 0)
      auc <- roc_curve$auc
      
      inner_metrics <- setNames(c(outer_index, inner_index, TP, TN, FP, FN, auc), names(inner_fold_metrics))
      inner_fold_metrics <- rbind(inner_fold_metrics, as.data.frame(t(inner_metrics)))
    }
    
    best_inner_auc <- mean(inner_auc_vals, na.rm = TRUE)
    # Reentrenar el mejor modelo con todos los datos de outer_train_data
    best_model <- glmboost(as.formula(paste(target_name, "~", paste(best_vars, collapse = "+"))),
                           data = outer_train_data, family = Binomial(),
                           control = boost_control(mstop = best_iterations))
    predictions <- predict(best_model, outer_test_data, type = "response")
    pred_class <- ifelse(predictions > threshold, "1", "0")
    confusion <- table(Actual = outer_test_data[[target_name]], Predicted = pred_class)
    
    if (nrow(confusion) < 2 || ncol(confusion) < 2) {
      TP <- ifelse(nrow(confusion) >= 2 && ncol(confusion) >= 2, confusion[2, 2], 0)
      TN <- ifelse(nrow(confusion) >= 1 && ncol(confusion) >= 1, confusion[1, 1], 0)
      FP <- ifelse(nrow(confusion) >= 1 && ncol(confusion) >= 2, confusion[1, 2], 0)
      FN <- ifelse(nrow(confusion) >= 2 && ncol(confusion) >= 1, confusion[2, 1], 0)
    } else {
      TP <- confusion[2, 2]
      TN <- confusion[1, 1]
      FP <- confusion[1, 2]
      FN <- confusion[2, 1]
    }
    auc <- roc(outer_test_data[[target_name]], predictions)$auc
    
    metrics <- setNames(c(outer_index, TP, TN, FP, FN, auc, 
                          best_iterations, paste(best_vars, collapse = ",")), 
                        names(performance_metrics))
    performance_metrics <- rbind(performance_metrics, as.data.frame(t(metrics)))
    
    cat("Performance:", paste(performance_metrics))
    cat(sprintf("Confusion Matrix for Fold %d:\n", outer_index))
    print(confusion)
    cat(sprintf("Metrics for Fold %d:\n", outer_index))
    cat(sprintf("TP: %d, TN: %d, FP: %d, FN: %d, AUC: %f, Best Iterations: %d, Best Variables: %s\n",
                TP, TN, FP, FN, auc, best_iterations, paste(best_vars, collapse = ",")))
  }
  
  list(performance_metrics = performance_metrics, inner_fold_metrics = inner_fold_metrics)
}

evaluate_with_seeds_AdaBoost_stepAUC <- 
  function(evaluation_function, data, target_var, vars, threshold, seeds, outer, inner, iterations) {
  results_list <- list()

  for (seed in seeds) {
    set.seed(seed)
    cat("Evaluating with seed:", seed, "\n")

    results <- evaluation_function(data = data, 
                                   target_name = target_var, 
                                   variables = vars, 
                                   threshold = threshold, 
                                   outer = outer, 
                                   inner = inner, 
                                   seed = seed,
                                   iterations = iterations)
    results_list[[paste0("seed_", seed)]] <- results
  }

  return(results_list)
  }

performance_adaboost_stepAUC <- evaluate_with_seeds_AdaBoost_stepAUC(
  evaluation_function = double_cross_validation_adaboost_stepAUC,
  data = data_factor,
  outer = 5,
  inner = 2,
  target_var = "PCR",
  vars = todas_variables_factor,
  threshold = .35,
  iterations = c(50, 100, 150, 200),
  seeds = seeds
)

load("RData_Files_Algorithms/AdaBoost_StepAUC.RData")

3.2.5 Embedded (Lasso)

double_cross_validation_adaboost_lasso <- function(data, target_name, outer, inner, iterations, variables, threshold, seed) {
  set.seed(seed)
  performance_metrics <- data.frame(Fold = integer(), TP = integer(), TN = integer(),
                                    FP = integer(), FN = integer(), AUC = numeric(),
                                    Best_Iteration = integer(), 
                                    Best_Variables = character(), stringsAsFactors = FALSE)
  
  inner_fold_metrics <- data.frame(Fold = integer(), Inner_Fold = integer(), 
                                   TP = integer(), TN = integer(), FP = integer(), 
                                   FN = integer(), AUC = numeric(), stringsAsFactors = FALSE)
  
  # Convertir las variables del dataset con model.matrix
  target <- as.factor(data[[target_name]])  # Aseguramos que la variable de respuesta sea factor
  features <- data[, setdiff(names(data), target_name)]
  data_matrix <- as.data.frame(model.matrix(~ ., data=features)[,-1])  # Eliminar el intercepto
  data_matrix[[target_name]] <- target
  
  outer_folds <- createFolds(data_matrix[[target_name]], k = outer)
  
  for (outer_index in seq_along(outer_folds)) {
    outer_train_data <- data_matrix[-outer_folds[[outer_index]], ]
    outer_test_data <- data_matrix[outer_folds[[outer_index]], ]
    
    cat(sprintf("Outer Fold %d\n", outer_index))
    cat("Training Data - Class 0: ", sum(outer_train_data[[target_name]] == 0), "\n")
    cat("Training Data - Class 1: ", sum(outer_train_data[[target_name]] == 1), "\n")
    cat("Testing Data - Class 0: ", sum(outer_test_data[[target_name]] == 0), "\n")
    cat("Testing Data - Class 1: ", sum(outer_test_data[[target_name]] == 1), "\n")
    
    inner_folds <- createFolds(outer_train_data[[target_name]], k = inner)
    best_model <- NULL
    best_auc <- 0
    best_iteration <- NULL
    best_vars <- colnames(outer_train_data)[!colnames(outer_train_data) %in% target_name]
    
    grid <- expand.grid(iteration = iterations)
    
    select_vars_lasso <- function(data, target_name) {
      target <- as.numeric(data[[target_name]]) - 1  # Convertimos a numérico para Lasso
      features <- data[, setdiff(names(data), target_name)]
      X <- model.matrix(~ ., data=features)[,-1]  # Eliminamos el intercepto
      y <- target
      
      # Ajustamos el modelo Lasso
      lasso_model <- cv.glmnet(X, y, alpha=1)
      lasso_coef <- coef(lasso_model, s = "lambda.min")
      
      # Seleccionamos las variables no nulas
      selected_vars <- rownames(lasso_coef)[lasso_coef[, 1] != 0]
      selected_vars <- setdiff(selected_vars, "(Intercept)")
      #cat("Selected variables by Lasso: ", paste(selected_vars, collapse=", "), "\n")
      return(selected_vars)
    }
    
    stepAUC <- function(vars, data, target_name, grid, inner_folds) {
      best_inner_auc <- 0
      best_inner_model <- NULL
      best_inner_iteration <- NULL
      best_inner_vars <- vars
      improved <- TRUE
      
      while (improved && length(vars) > 1) {
        improved <- FALSE
        best_auc_in_step <- best_inner_auc
        
        for (params in 1:nrow(grid)) {
          auc_vals <- numeric()
          
          for (inner_index in seq_along(inner_folds)) {
            inner_train_data <- outer_train_data[-inner_folds[[inner_index]], ]
            inner_test_data <- outer_train_data[inner_folds[[inner_index]], ]
            
            selected_vars <- select_vars_lasso(inner_train_data, target_name)
            formula <- as.formula(paste(target_name, "~", paste(selected_vars, collapse = "+")))
            #cat("Testing with formula: ", deparse(formula), "\n")
            #cat("Grid parameters - iteration: ", grid$iteration[params], "\n")
            
            result <- tryCatch({
              model <- blackboost(formula, data = inner_train_data, family = Binomial(), 
                                  control = boost_control(mstop = grid$iteration[params]))
              
              predictions <- predict(model, inner_test_data, type = "response")
              roc_curve <- roc(inner_test_data[[target_name]], predictions)
              auc_vals[inner_index] <- roc_curve$auc
              #cat("AUC for inner fold ", inner_index, ": ", roc_curve$auc, "\n")
              
            }, error = function(e) {
              cat("Error with variables: ", paste(selected_vars, collapse = ", "), "\n")
              cat("Error message: ", e$message, "\n")
              auc_vals[inner_index] <- NA  # Set to NA or some other indicator of failure
            })
          }
          
          mean_auc <- mean(auc_vals, na.rm = TRUE)
          #cat("Mean AUC for params - iteration: ", grid$iteration[params], ": ", mean_auc, "\n")
          if (!is.na(mean_auc) && mean_auc > best_auc_in_step) {
            best_auc_in_step <- mean_auc
            best_inner_model <- model
            best_inner_iteration <- grid$iteration[params]
            best_inner_vars <- selected_vars
            improved <- TRUE
          }
        }
        
        if (improved) {
          vars <- best_inner_vars
          best_inner_auc <- best_auc_in_step
          cat("Vars improved: ", paste(vars, collapse=", "), "\n")
          cat("AUC: ", best_inner_auc, "\n")
        }
      }
      
      list(model = best_inner_model, auc = best_inner_auc, 
           iteration = best_inner_iteration, vars = best_inner_vars)
    }
    
    best_result <- stepAUC(best_vars, outer_train_data, target_name, grid, inner_folds)
    best_model <- best_result$model
    best_auc <- best_result$auc
    best_iteration <- best_result$iteration
    best_vars <- best_result$vars
    
    if (is.null(best_model)) {
      cat("Error: No valid model found for outer fold ", outer_index, "\n")
      next
    }

    # Evaluar el mejor modelo en el validation split del inner loop
    inner_auc_vals <- numeric()
    for (inner_index in seq_along(inner_folds)) {
      inner_train_data <- outer_train_data[-inner_folds[[inner_index]], ]
      inner_test_data <- outer_train_data[inner_folds[[inner_index]], ]
      
      predictions <- predict(best_model, inner_test_data, type = "response")
      roc_curve <- roc(inner_test_data[[target_name]], predictions)
      inner_auc_vals[inner_index] <- roc_curve$auc
      cat("Inner fold ", inner_index, " AUC: ", roc_curve$auc, "\n")

      pred_class <- ifelse(predictions > threshold, 1, 0)
      confusion <- table(Actual = inner_test_data[[target_name]], Predicted = pred_class)
      TP <- ifelse(nrow(confusion) >= 2 && ncol(confusion) >= 2, confusion[2, 2], 0)
      TN <- ifelse(nrow(confusion) >= 1 && ncol(confusion) >= 1, confusion[1, 1], 0)
      FP <- ifelse(nrow(confusion) >= 1 && ncol(confusion) >= 2, confusion[1, 2], 0)
      FN <- ifelse(nrow(confusion) >= 2 && ncol(confusion) >= 1, confusion[2, 1], 0)
      auc <- roc_curve$auc
      
      inner_metrics <- setNames(c(outer_index, inner_index, TP, TN, FP, FN, auc), names(inner_fold_metrics))
      inner_fold_metrics <- rbind(inner_fold_metrics, as.data.frame(t(inner_metrics)))
    }
    
    best_inner_auc <- mean(inner_auc_vals, na.rm = TRUE)
    cat("Best inner AUC: ", best_inner_auc, "\n")
    # Reentrenar el mejor modelo con todos los datos de outer_train_data
    best_model <- glmboost(as.formula(paste(target_name, "~", paste(best_vars, collapse = "+"))),
                           data = outer_train_data, family = Binomial(),
                           control = boost_control(mstop = best_iteration))
    predictions <- predict(best_model, outer_test_data, type = "response")
    pred_class <- ifelse(predictions > threshold, 1, 0)  # Asignar clases basadas en umbral
    confusion <- table(Actual = outer_test_data[[target_name]], Predicted = pred_class)
    
    if (nrow(confusion) < 2 || ncol(confusion) < 2) {
      TP <- ifelse(nrow(confusion) >= 2 && ncol(confusion) >= 2, confusion[2, 2], 0)
      TN <- ifelse(nrow(confusion) >= 1 && ncol(confusion) >= 1, confusion[1, 1], 0)
      FP <- ifelse(nrow(confusion) >= 1 && ncol(confusion) >= 2, confusion[1, 2], 0)
      FN <- ifelse(nrow(confusion) >= 2 && ncol(confusion) >= 1, confusion[2, 1], 0)
    } else {
      TP <- confusion[2, 2]
      TN <- confusion[1, 1]
      FP <- confusion[1, 2]
      FN <- confusion[2, 1]
    }
    auc <- roc(outer_test_data[[target_name]], predictions)$auc
    
    metrics <- setNames(c(outer_index, TP, TN, FP, FN, auc, 
                          best_iteration, paste(best_vars, collapse = ",")), 
                        names(performance_metrics))
    performance_metrics <- rbind(performance_metrics, as.data.frame(t(metrics)))
    
    cat("Performance:\n", paste(performance_metrics), "\n")
    cat(sprintf("Confusion Matrix for Fold %d:\n", outer_index))
    print(confusion)
    cat(sprintf("Metrics for Fold %d:\n", outer_index))
    cat(sprintf("TP: %d, TN: %d, FP: %d, FN: %d, AUC: %f, Best Iteration: %d, Best Variables: %s\n",
                TP, TN, FP, FN, auc, best_iteration, paste(best_vars, collapse = ",")))
  }
  
  list(performance_metrics = performance_metrics, inner_fold_metrics = inner_fold_metrics)
}
evaluate_with_seeds_AdaBoost_Lasso <- 
  function(evaluation_function, data, target_var, vars, threshold, seeds, outer, inner, iterations) {
  results_list <- list()

  for (seed in seeds) {
    set.seed(seed)
    cat("Evaluating with seed:", seed, "\n")

    results <- evaluation_function(data = data, 
                                   target_name = target_var, 
                                   variables = vars, 
                                   threshold = threshold, 
                                   outer = outer, 
                                   inner = inner, 
                                   seed = seed,
                                   iterations = iterations)
    results_list[[paste0("seed_", seed)]] <- results
  }

  return(results_list)
  }

performance_adaboost_lasso <- evaluate_with_seeds_AdaBoost_stepAUC(
  evaluation_function = double_cross_validation_adaboost_lasso,
  data = data_factor,
  outer = 5,
  inner = 2,
  target_var = "PCR",
  vars = todas_variables_factor,
  threshold = .35,
  iterations = c(50, 100, 150, 200),
  seeds = seeds
)

load("RData_Files_Algorithms/AdaBoost_Lasso.RData")

3.3 Random Forest

3.3.1 Aparente

evaluate_apparent_performance_model_randomForest <- function(data, target_var, model_func, vars = NULL, threshold = 0.5) {
  if (is.null(vars)) {
    vars <- setdiff(names(data), target_var)
  }
  
  data[[target_var]] <- factor(data[[target_var]], levels = c(0, 1))
  formula <- as.formula(paste(target_var, "~", paste(vars, collapse = "+")))
  
  model <- model_func(formula, data)
  predictions <- predict(model, newdata = data, type = "prob")[, 2]
  
  actual_classes <- data[[target_var]]
  predicted_classes <- ifelse(predictions > threshold, 1, 0)
  confusion_matrix <- table(predicted_classes, actual_classes)
  
  tp <- ifelse("1" %in% rownames(confusion_matrix) && "1" %in% colnames(confusion_matrix), confusion_matrix["1", "1"], 0)
  tn <- ifelse("0" %in% rownames(confusion_matrix) && "0" %in% colnames(confusion_matrix), confusion_matrix["0", "0"], 0)
  fp <- ifelse("1" %in% rownames(confusion_matrix) && "0" %in% colnames(confusion_matrix), confusion_matrix["1", "0"], 0)
  fn <- ifelse("0" %in% rownames(confusion_matrix) && "1" %in% colnames(confusion_matrix), confusion_matrix["0", "1"], 0)
  
  accuracy <- (tp + tn) / (tp + tn + fp + fn)
  precision <- ifelse(tp + fp > 0, tp / (tp + fp), 0)
  recall <- ifelse(tp + fn > 0, tp / (tp + fn), 0)
  f1_score <- ifelse(precision + recall > 0, 2 * (precision * recall) / (precision + recall), 0)
  roc_obj <- pROC::roc(actual_classes, predictions)
  
  list(confusion_matrix = confusion_matrix, accuracy = accuracy, precision = precision, 
       recall = recall, f1_score = f1_score, roc_curve = roc_obj)
}

rf_model <- function(formula, data) {
  set.seed(90)
  randomForest(formula = formula, data = data, ntree = 100, mtry = 2)
}


resultadosAparentesRF <- evaluate_apparent_performance_model_randomForest(data = data_factor, target_var = "PCR",
                                                             model_func = rf_model,
                                                             vars = ".",
                                                             threshold = .35)

load("RData_Files_Algorithms/RandomForest_Aparente.RData")

3.3.2 Todas las Variables

double_cross_validation_randomforest <- function(data, target_name, outer, inner, ntree_values, mtry_values, variables, threshold, seed) {
  set.seed(seed)
  performance_metrics <- data.frame(Fold = integer(), TP = integer(), TN = integer(),
                                    FP = integer(), FN = integer(), AUC = numeric(), 
                                    Best_ntree = integer(), Best_mtry = integer(), stringsAsFactors = FALSE)
  
  inner_fold_metrics <- data.frame(Fold = integer(), Inner_Fold = integer(), 
                                   TP = integer(), TN = integer(), FP = integer(), 
                                   FN = integer(), AUC = numeric(), stringsAsFactors = FALSE)
  
  data[[target_name]] <- as.factor(data[[target_name]])
  outer_folds <- createFolds(data[[target_name]], k = outer, list = TRUE, returnTrain = FALSE)
  
  for (outer_index in seq_along(outer_folds)) {
    outer_train_data <- data[-outer_folds[[outer_index]], ]
    outer_test_data <- data[outer_folds[[outer_index]], ]
    
    cat(sprintf("Outer Fold %d\n", outer_index))
    cat("Training Data - Class 0: ", sum(outer_train_data[[target_name]] == "0"), "\n")
    cat("Training Data - Class 1: ", sum(outer_train_data[[target_name]] == "1"), "\n")
    cat("Testing Data - Class 0: ", sum(outer_test_data[[target_name]] == "0"), "\n")
    cat("Testing Data - Class 1: ", sum(outer_test_data[[target_name]] == "1"), "\n")
    
    if (sum(outer_train_data[[target_name]] == "0") == 0 || sum(outer_train_data[[target_name]] == "1") == 0) {
      cat("Skipping outer fold due to lack of class diversity\n")
      next
    }
    
    inner_folds <- createFolds(outer_train_data[[target_name]], k = inner, list = TRUE, returnTrain = FALSE)
    best_model <- NULL
    best_auc <- 0
    best_ntree <- NULL
    best_mtry <- NULL

    grid <- expand.grid(ntree = ntree_values, mtry = mtry_values)
    formula <- reformulate(variables, target_name)

    for (params in 1:nrow(grid)) {
      inner_auc <- numeric()

      for (inner_index in seq_along(inner_folds)) {
        inner_train_data <- outer_train_data[-inner_folds[[inner_index]], ]
        inner_test_data <- outer_train_data[inner_folds[[inner_index]], ]

        if (sum(inner_train_data[[target_name]] == "0") == 0 || sum(inner_train_data[[target_name]] == "1") == 0) {
          cat("Skipping inner fold due to lack of class diversity\n")
          next
        }

        result <- tryCatch({
          model <- randomForest(formula, data = inner_train_data, ntree = grid$ntree[params], mtry = grid$mtry[params])
          
          predictions <- predict(model, inner_test_data, type = "prob")[, 2]
          roc_curve <- roc(inner_test_data[[target_name]], predictions)
          inner_auc[inner_index] <- roc_curve$auc
        }, error = function(e) {
          cat("Error with parameters: ntree =", grid$ntree[params], ", mtry =", grid$mtry[params], "\n")
          cat("Error message: ", e$message, "\n")
          inner_auc[inner_index] <- 0
        })
      }

      if (mean(inner_auc, na.rm = TRUE) > best_auc) {
        best_auc <- mean(inner_auc, na.rm = TRUE)
        best_model <- model
        best_ntree <- grid$ntree[params]
        best_mtry <- grid$mtry[params]
      }
    }

    if (is.null(best_model)) {
      cat("No valid model found for outer fold ", outer_index, "\n")
      next
    }

    # Evaluar el mejor modelo en el validation split del inner loop 
    inner_auc_vals <- numeric()
    for (inner_index in seq_along(inner_folds)) {
      inner_train_data <- outer_train_data[-inner_folds[[inner_index]], ]
      inner_test_data <- outer_train_data[inner_folds[[inner_index]], ]
      
      predictions <- predict(best_model, inner_test_data, type = "prob")[, 2]
      roc_curve <- roc(inner_test_data[[target_name]], predictions)
      inner_auc_vals[inner_index] <- roc_curve$auc

      pred_class <- ifelse(predictions > threshold, "1", "0")
      confusion <- table(Actual = inner_test_data[[target_name]], Predicted = pred_class)
      TP <- ifelse(nrow(confusion) >= 2 && ncol(confusion) >= 2, confusion[2, 2], 0)
      TN <- ifelse(nrow(confusion) >= 1 && ncol(confusion) >= 1, confusion[1, 1], 0)
      FP <- ifelse(nrow(confusion) >= 1 && ncol(confusion) >= 2, confusion[1, 2], 0)
      FN <- ifelse(nrow(confusion) >= 2 && ncol(confusion) >= 1, confusion[2, 1], 0)
      auc <- roc_curve$auc
      
      inner_metrics <- setNames(c(outer_index, inner_index, TP, TN, FP, FN, auc), names(inner_fold_metrics))
      inner_fold_metrics <- rbind(inner_fold_metrics, as.data.frame(t(inner_metrics)))
    }

    best_model <- randomForest(formula, data = outer_train_data, ntree = best_ntree, mtry = best_mtry)
    predictions <- predict(best_model, outer_test_data, type = "prob")[, 2]
    pred_class <- ifelse(predictions > threshold, "1", "0")
    confusion <- table(Actual = outer_test_data[[target_name]], Predicted = pred_class)

    if (nrow(confusion) < 2 || ncol(confusion) < 2) {
      TP <- ifelse(nrow(confusion) >= 2 && ncol(confusion) >= 2, confusion[2, 2], 0)
      TN <- ifelse(nrow(confusion) >= 1 && ncol(confusion) >= 1, confusion[1, 1], 0)
      FP <- ifelse(nrow(confusion) >= 1 && ncol(confusion) >= 2, confusion[1, 2], 0)
      FN <- ifelse(nrow(confusion) >= 2 && ncol(confusion) >= 1, confusion[2, 1], 0)
    } else {
      TP <- confusion[2, 2]
      TN <- confusion[1, 1]
      FP <- confusion[1, 2]
      FN <- confusion[2, 1]
    }
    auc <- roc(outer_test_data[[target_name]], predictions)$auc

    metrics <- setNames(c(outer_index, TP, TN, FP, FN, auc, best_ntree, best_mtry), names(performance_metrics))
    performance_metrics <- rbind(performance_metrics, as.data.frame(t(metrics)))
  }

  return(list(performance_metrics = performance_metrics, inner_fold_metrics = inner_fold_metrics))
}

evaluate_with_seeds_RandomForest_TodasVariables_Asociacion <- 
  function(evaluation_function, data, target_var, vars, threshold, seeds, outer, inner, ntree_values , mtry_values) {
  results_list <- list()

  for (seed in seeds) {
    set.seed(seed)
    cat("Evaluating with seed:", seed, "\n")

    results <- evaluation_function(data = data, 
                                   target_name = target_var, 
                                   variables = vars, 
                                   threshold = threshold, 
                                   outer = outer, 
                                   inner = inner, 
                                   seed = seed,
                                   ntree_values = ntree_values,
                                   mtry_values = mtry_values)
    results_list[[paste0("seed_", seed)]] <- results
  }

  return(results_list)
  }

performance_RandomForest_todasvariables <- evaluate_with_seeds_RandomForest_TodasVariables_Asociacion(
  evaluation_function = double_cross_validation_randomforest,
  data = data_factor,
  outer = 5,
  inner = 2,
  target_var = "PCR",
  ntree_values = c(100, 200, 300, 400, 500),###
  mtry_values  = c(floor(sqrt(length(todas_variables_factor))), 
                   floor(sqrt(length(todas_variables_factor))) - 1, 
                   floor(sqrt(length(todas_variables_factor))) + 1),
  vars = todas_variables_factor,
  threshold = .35,
  seeds = seeds
)

load("RData_Files_Algorithms/RandomForest_TodasVariables.RData")

3.3.3 Filtrado (Asociación)

performance_RandomForest_Aparente <- evaluate_with_seeds_RandomForest_TodasVariables_Asociacion(
  evaluation_function = double_cross_validation_randomforest,
  data = data_factor,
  outer = 5,
  inner = 2,
  target_var = "PCR",
  ntree_values = c(100, 200, 300, 400, 500),###
  mtry_values  = c(floor(sqrt(length(asociacion))), 
                   floor(sqrt(length(asociacion))) - 1, 
                   floor(sqrt(length(asociacion))) + 1),
  vars = asociacion,
  threshold = .35,
  seeds = seeds
)

load("RData_Files_Algorithms/RandomForest_Asociacion.RData")

3.3.4 Wrapped (StepAUC)

double_cross_validation_randomforest_stepAUC <- function(data, target_name, outer, inner, ntree_values, mtry_values, variables, threshold, seed) {
  set.seed(seed)
  performance_metrics <- data.frame(Fold = integer(), TP = integer(), TN = integer(),
                                    FP = integer(), FN = integer(), AUC = numeric(),
                                    Best_Ntree = integer(), Best_Mtry = integer(), 
                                    Best_Variables = character(), stringsAsFactors = FALSE)
  
  inner_fold_metrics <- data.frame(Fold = integer(), Inner_Fold = integer(), 
                                   TP = integer(), TN = integer(), FP = integer(), 
                                   FN = integer(), AUC = numeric(), stringsAsFactors = FALSE)
  
  data[[target_name]] <- as.factor(data[[target_name]])
  outer_folds <- createFolds(data[[target_name]], k = outer)
  
  for (outer_index in seq_along(outer_folds)) {
    outer_train_data <- data[-outer_folds[[outer_index]], ]
    outer_test_data <- data[outer_folds[[outer_index]], ]
    
    cat(sprintf("Outer Fold %d\n", outer_index))
    cat("Training Data - Class 0: ", sum(outer_train_data[[target_name]] == "0"), "\n")
    cat("Training Data - Class 1: ", sum(outer_train_data[[target_name]] == "1"), "\n")
    cat("Testing Data - Class 0: ", sum(outer_test_data[[target_name]] == "0"), "\n")
    cat("Testing Data - Class 1: ", sum(outer_test_data[[target_name]] == "1"), "\n")
    
    inner_folds <- createFolds(outer_train_data[[target_name]], k = inner)
    best_model <- NULL
    best_auc <- 0
    best_ntree <- NULL
    best_mtry <- NULL
    best_vars <- variables
    
    grid <- expand.grid(ntree = ntree_values, mtry = mtry_values)
    
    stepAUC <- function(vars, data, target_name, grid, inner_folds) {
      current_vars <- vars
      best_inner_auc <- 0
      best_inner_model <- NULL
      best_inner_ntree <- NULL
      best_inner_mtry <- NULL
      best_inner_vars <- current_vars
      improved <- TRUE
    
      while (improved && length(current_vars) > 1) {
        improved <- FALSE
        best_auc_in_step <- best_inner_auc
    
        for (var in current_vars) {
          temp_vars <- setdiff(current_vars, var)
          inner_auc <- numeric()
    
          for (params in 1:nrow(grid)) {
            auc_vals <- numeric()
    
            for (inner_index in seq_along(inner_folds)) {
              inner_train_data <- data[-inner_folds[[inner_index]], ]
              inner_test_data <- data[inner_folds[[inner_index]], ]
    
              result <- tryCatch({
                model <- randomForest(as.formula(paste(target_name, "~", paste(temp_vars, collapse = "+"))),
                                      data = inner_train_data, ntree = grid$ntree[params], mtry = grid$mtry[params])
    
                predictions <- predict(model, inner_test_data, type = "prob")[, 2]
                roc_curve <- roc(inner_test_data[[target_name]], predictions)
                auc_vals[inner_index] <- roc_curve$auc
    
              }, error = function(e) {
                cat("Error with variables:", paste(temp_vars, collapse = ", "), "\n")
                auc_vals[inner_index] <- 0
              })
            }
    
            if (mean(auc_vals, na.rm = TRUE) > best_auc_in_step) {
              best_auc_in_step <- mean(auc_vals, na.rm = TRUE)
              best_inner_model <- model
              best_inner_ntree <- grid$ntree[params]
              best_inner_mtry <- grid$mtry[params]
              best_inner_vars <- temp_vars
              improved <- TRUE
            }
          }
        }
    
        if (improved) {
          current_vars <- best_inner_vars
          best_inner_auc <- best_auc_in_step
          cat("Vars improved: ", paste(current_vars, collapse=", "), "\n")
          cat("AUC: ", paste(best_inner_auc), "\n")
        }
      }
    
      list(model = best_inner_model, auc = best_inner_auc, 
           ntree = best_inner_ntree, mtry = best_inner_mtry, 
           vars = best_inner_vars)
    }

    
    best_result <- stepAUC(variables, outer_train_data, target_name, grid, inner_folds)
    best_model <- best_result$model
    best_auc <- best_result$auc
    best_ntree <- best_result$ntree
    best_mtry <- best_result$mtry
    best_vars <- best_result$vars

    # Evaluar el mejor modelo en el validation split del inner loop
    inner_auc_vals <- numeric()
    for (inner_index in seq_along(inner_folds)) {
      inner_train_data <- outer_train_data[-inner_folds[[inner_index]], ]
      inner_test_data <- outer_train_data[inner_folds[[inner_index]], ]
    
      predictions <- predict(best_model, inner_test_data, type = "prob")[, 2]
      roc_curve <- roc(inner_test_data[[target_name]], predictions)
      inner_auc_vals[inner_index] <- roc_curve$auc
    
      pred_class <- ifelse(predictions > threshold, "1", "0")
      confusion <- table(Actual = inner_test_data[[target_name]], Predicted = pred_class)
      TP <- ifelse(nrow(confusion) >= 2 && ncol(confusion) >= 2, confusion[2, 2], 0)
      TN <- ifelse(nrow(confusion) >= 1 && ncol(confusion) >= 1, confusion[1, 1], 0)
      FP <- ifelse(nrow(confusion) >= 1 && ncol(confusion) >= 2, confusion[1, 2], 0)
      FN <- ifelse(nrow(confusion) >= 2 && ncol(confusion) >= 1, confusion[2, 1], 0)
      auc <- roc_curve$auc
    
      inner_metrics <- setNames(c(outer_index, inner_index, TP, TN, FP, FN, auc), names(inner_fold_metrics))
      inner_fold_metrics <- rbind(inner_fold_metrics, as.data.frame(t(inner_metrics)))
    }

    
    best_inner_auc <- mean(inner_auc_vals, na.rm = TRUE)
    # Reentrenar el mejor modelo con todos los datos de outer_train_data
    best_model <- randomForest(as.formula(paste(target_name, "~", paste(best_vars, collapse = "+"))),
                    data = outer_train_data, ntree = best_ntree, mtry = best_mtry)

    predictions <- predict(best_model, outer_test_data, type = "prob")[, 2]
    pred_class <- ifelse(predictions > threshold, "1", "0")
    confusion <- table(Actual = outer_test_data[[target_name]], Predicted = pred_class)
    
    if (nrow(confusion) < 2 || ncol(confusion) < 2) {
      TP <- ifelse(nrow(confusion) >= 2 && ncol(confusion) >= 2, confusion[2, 2], 0)
      TN <- ifelse(nrow(confusion) >= 1 && ncol(confusion) >= 1, confusion[1, 1], 0)
      FP <- ifelse(nrow(confusion) >= 1 && ncol(confusion) >= 2, confusion[1, 2], 0)
      FN <- ifelse(nrow(confusion) >= 2 && ncol(confusion) >= 1, confusion[2, 1], 0)
    } else {
      TP <- confusion[2, 2]
      TN <- confusion[1, 1]
      FP <- confusion[1, 2]
      FN <- confusion[2, 1]
    }
    auc <- roc(outer_test_data[[target_name]], predictions)$auc
    
    metrics <- setNames(c(outer_index, TP, TN, FP, FN, auc, 
                          best_ntree, best_mtry, paste(best_vars, collapse = ",")), 
                        names(performance_metrics))
    performance_metrics <- rbind(performance_metrics, as.data.frame(t(metrics)))
    
    cat("Performance:", paste(performance_metrics))
    cat(sprintf("Confusion Matrix for Fold %d:\n", outer_index))
    print(confusion)
    cat(sprintf("Metrics for Fold %d:\n", outer_index))
    cat(sprintf("TP: %d, TN: %d, FP: %d, FN: %d, AUC: %f, Best Ntree: %d, Best Mtry: %d, Best Variables: %s\n",
                TP, TN, FP, FN, auc, best_ntree, best_mtry, paste(best_vars, collapse = ",")))
  }
  
  return(list(performance_metrics = performance_metrics, inner_fold_metrics = inner_fold_metrics))
}

evaluate_with_seeds_RandomForest_StepAUC <- 
  function(evaluation_function, data, target_var, vars, threshold, seeds, outer, inner, ntree_values, mtry_values) {
  results_list <- list()

  for (seed in seeds) {
    set.seed(seed)
    cat("Evaluating with seed:", seed, "\n")

    results <- evaluation_function(data = data, 
                                   target_name = target_var, 
                                   variables = vars, 
                                   threshold = threshold, 
                                   outer = outer, 
                                   inner = inner, 
                                   seed = seed,
                                   ntree_values = ntree_values,
                                   mtry_values = mtry_values)
    results_list[[paste0("seed_", seed)]] <- results
  }

  return(results_list)
  }

performance_RandomForest_StepAUC <- evaluate_with_seeds_RandomForest_StepAUC(
  evaluation_function = double_cross_validation_randomforest_stepAUC,
  data = data_factor,
  outer = 5,
  inner = 2,
  target_var = "PCR",
  ntree_values = c(100, 200, 300, 400, 500),
  mtry_values  = c(floor(sqrt(length(todas_variables_factor))), 
                   floor(sqrt(length(todas_variables_factor))) - 1, 
                   floor(sqrt(length(todas_variables_factor))) + 1),
  vars = todas_variables_factor,
  threshold = .35,
  seeds = seeds
)

load("RData_Files_Algorithms/RandomForest_StepAUC.RData")

3.3.5 Embedded (Lasso)

Debido a dificultades del código, no se han podido obtener resultados del método de selección de variables embedded con el algoritmo Lasso para al modelo de Árbol de Decisión (DT).

3.4 Decision Tree

3.4.1 Aparente

evaluate_aparent_performance_model_rpart <- function(data, target_var, model_func, vars = NULL, threshold = 0.5) {
  if (is.null(vars)) {
    vars <- setdiff(names(data), target_var)
  }
  
  data[[target_var]] <- factor(data[[target_var]], levels = c(0, 1))
  formula <- as.formula(paste(target_var, "~", paste(vars, collapse = "+")))
  
  model <- model_func(formula, data)
  predictions <- predict(model, newdata = data, type = "prob")[, 2]
  
  actual_classes <- data[[target_var]]
  predicted_classes <- ifelse(predictions > threshold, 1, 0)
  confusion_matrix <- table(predicted_classes, actual_classes)
  
  tp <- ifelse("1" %in% rownames(confusion_matrix) && "1" %in% colnames(confusion_matrix), confusion_matrix["1", "1"], 0)
  tn <- ifelse("0" %in% rownames(confusion_matrix) && "0" %in% colnames(confusion_matrix), confusion_matrix["0", "0"], 0)
  fp <- ifelse("1" %in% rownames(confusion_matrix) && "0" %in% colnames(confusion_matrix), confusion_matrix["1", "0"], 0)
  fn <- ifelse("0" %in% rownames(confusion_matrix) && "1" %in% colnames(confusion_matrix), confusion_matrix["0", "1"], 0)
  
  accuracy <- (tp + tn) / (tp + tn + fp + fn)
  precision <- ifelse(tp + fp > 0, tp / (tp + fp), 0)
  recall <- ifelse(tp + fn > 0, tp / (tp + fn), 0)
  f1_score <- ifelse(precision + recall > 0, 2 * (precision * recall) / (precision + recall), 0)
  roc_obj <- pROC::roc(actual_classes, predictions)
  
  list(confusion_matrix = confusion_matrix, accuracy = accuracy, precision = precision, 
       recall = recall, f1_score = f1_score, roc_curve = roc_obj)
}


dt_model <- function(formula, data) {
  set.seed(90)
  rpart(formula = formula, data = data, method = "class", control = rpart.control(minsplit = 2, cp = 0))
}


resultadosAparentesDT <- evaluate_aparent_performance_model_rpart(data = data_factor, target_var = "PCR",
                                                             model_func = dt_model,
                                                             vars = ".",
                                                             threshold = .35)

load("RData_Files_Algorithms/DecisionTree_Aparente.RData")

3.4.2 Todas las Variables

double_cross_validation_rpart <- function(data, target_name, outer, inner, cp_values, minsplit_values, variables, threshold, seed) {
  set.seed(seed)
  performance_metrics <- data.frame(Fold = integer(), TP = integer(), TN = integer(),
                                    FP = integer(), FN = integer(), AUC = numeric(),
                                    Best_cp = numeric(), Best_minsplit = integer(), stringsAsFactors = FALSE)
  
  inner_fold_metrics <- data.frame(Fold = integer(), Inner_Fold = integer(), 
                                   TP = integer(), TN = integer(), FP = integer(), 
                                   FN = integer(), AUC = numeric(), stringsAsFactors = FALSE)
  
  data[[target_name]] <- as.factor(data[[target_name]])
  outer_folds <- createFolds(data[[target_name]], k = outer)
  
  for (outer_index in seq_along(outer_folds)) {
    outer_train_data <- data[-outer_folds[[outer_index]], ]
    outer_test_data <- data[outer_folds[[outer_index]], ]
    
    # Imprime la cantidad de instancias de cada clase en el outer train y test set
    cat(sprintf("Outer Fold %d\n", outer_index))
    cat("Training Data - Class 0: ", sum(outer_train_data[[target_name]] == 0), "\n")
    cat("Training Data - Class 1: ", sum(outer_train_data[[target_name]] == 1), "\n")
    cat("Testing Data - Class 0: ", sum(outer_test_data[[target_name]] == 0), "\n")
    cat("Testing Data - Class 1: ", sum(outer_test_data[[target_name]] == 1), "\n")
    
    inner_folds <- createFolds(outer_train_data[[target_name]], k = inner)
    best_model <- NULL
    best_auc <- 0
    best_cp <- NULL
    best_minsplit <- NULL
    
    grid <- expand.grid(cp = cp_values, minsplit = minsplit_values)
    
    for (params in 1:nrow(grid)) {
      inner_auc <- numeric()
      
      for (inner_index in seq_along(inner_folds)) {
        inner_train_data <- outer_train_data[-inner_folds[[inner_index]], ]
        inner_test_data <- outer_train_data[inner_folds[[inner_index]], ]
        
        model <- rpart(formula = as.formula(paste(target_name, "~", paste(variables, collapse = "+"))),
                       data = inner_train_data, control = rpart.control(cp = grid$cp[params], minsplit = grid$minsplit[params]))
        
        predictions <- predict(model, inner_test_data, type = "prob")[, 2]
        roc_curve <- roc(inner_test_data[[target_name]], predictions)
        inner_auc[inner_index] <- roc_curve$auc
      }
      
      if (mean(inner_auc) > best_auc) {
        best_auc = mean(inner_auc)
        best_model = model
        best_cp = grid$cp[params]
        best_minsplit = grid$minsplit[params]
      }
    }
    
    # Evaluar el mejor modelo en el validation split del inner loop 
    inner_auc_vals <- numeric()
    for (inner_index in seq_along(inner_folds)) {
      inner_train_data <- outer_train_data[-inner_folds[[inner_index]], ]
      inner_test_data <- outer_train_data[inner_folds[[inner_index]], ]
      
      predictions <- predict(best_model, inner_test_data, type = "prob")[, 2]
      roc_curve <- roc(inner_test_data[[target_name]], predictions)
      inner_auc_vals[inner_index] <- roc_curve$auc

      pred_class <- ifelse(predictions > threshold, 1, 0)
      confusion <- table(Actual = inner_test_data[[target_name]], Predicted = pred_class)
      TP <- ifelse(nrow(confusion) >= 2 && ncol(confusion) >= 2, confusion[2, 2], 0)
      TN <- ifelse(nrow(confusion) >= 1 && ncol(confusion) >= 1, confusion[1, 1], 0)
      FP <- ifelse(nrow(confusion) >= 1 && ncol(confusion) >= 2, confusion[1, 2], 0)
      FN <- ifelse(nrow(confusion) >= 2 && ncol(confusion) >= 1, confusion[2, 1], 0)
      auc <- roc_curve$auc
      
      inner_metrics <- setNames(c(outer_index, inner_index, TP, TN, FP, FN, auc), names(inner_fold_metrics))
      inner_fold_metrics <- rbind(inner_fold_metrics, as.data.frame(t(inner_metrics)))
    }
    
    best_model <- rpart(formula = as.formula(paste(target_name, "~", paste(variables, collapse = "+"))),
                       data = outer_train_data, control = rpart.control(cp = best_cp, minsplit = best_minsplit))
    
    predictions <- predict(best_model, outer_test_data, type = "prob")[, 2]
    pred_class <- ifelse(predictions > threshold, 1, 0)  # Asignar clases basadas en umbral
    confusion <- table(Actual = outer_test_data[[target_name]], Predicted = pred_class)
    
    if (nrow(confusion) < 2 || ncol(confusion) < 2) {
      TP <- ifelse(nrow(confusion) >= 2 && ncol(confusion) >= 2, confusion[2, 2], 0)
      TN <- ifelse(nrow(confusion) >= 1 && ncol(confusion) >= 1, confusion[1, 1], 0)
      FP <- ifelse(nrow(confusion) >= 1 && ncol(confusion) >= 2, confusion[1, 2], 0)
      FN <- ifelse(nrow(confusion) >= 2 && ncol(confusion) >= 1, confusion[2, 1], 0)
    } else {
      TP <- confusion[2, 2]
      TN <- confusion[1, 1]
      FP <- confusion[1, 2]
      FN <- confusion[2, 1]
    }
    auc = roc(outer_test_data[[target_name]], predictions)$auc
    
    metrics <- setNames(c(outer_index, TP, TN, FP, FN, auc, best_cp, best_minsplit), names(performance_metrics))
    performance_metrics <- rbind(performance_metrics, as.data.frame(t(metrics)))
  }
  
  return(list(performance_metrics = performance_metrics, inner_fold_metrics = inner_fold_metrics))
}


evaluate_with_seeds_DecisionTree_TodasVariables_Asociacion <- 
  function(evaluation_function, data, target_var, vars, threshold, seeds, outer, inner, cp_values, minsplit_values) {
  results_list <- list()

  for (seed in seeds) {
    set.seed(seed)
    cat("Evaluating with seed:", seed, "\n")

    results <- evaluation_function(data = data, 
                                   target_name = target_var, 
                                   variables = vars, 
                                   threshold = threshold, 
                                   outer = outer, 
                                   inner = inner, 
                                   seed = seed,
                                   cp_values = cp_values,
                                   minsplit_values = minsplit_values)
    results_list[[paste0("seed_", seed)]] <- results
  }

  return(results_list)
  }

performance_rpart_todasvariables <- evaluate_with_seeds_DecisionTree_TodasVariables_Asociacion(
  evaluation_function = double_cross_validation_rpart,
  data = data_factor,
  outer = 5,
  inner = 2,
  target_var = "PCR",
  cp_values = c(0.001, 0.01, 0.1),
  minsplit_values = c(2, 5, 10, 15, 20),
  vars = todas_variables_factor,
  threshold = .35,
  seeds = seeds
)

load("RData_Files_Algorithms/DecisionTree_TodasVariables.RData")

3.4.3 Filtrado (Asociación)

performance_rpart_Asociacion <- evaluate_with_seeds_DecisionTree_TodasVariables_Asociacion(
  evaluation_function = double_cross_validation_rpart,
  data = data_factor,
  outer = 5,
  inner = 2,
  target_var = "PCR",
  cp_values = c(0.001, 0.01, 0.1),
  minsplit_values = c(2, 5, 10, 15, 20),
  vars = asociacion,
  threshold = .35,
  seeds = seeds
)

load("RData_Files_Algorithms/DecisionTree_Asociacion.RData")

3.4.4 Wrapped (StepAUC)

double_cross_validation_rpart_stepAUC <- function(data, target_name, outer, inner, cps, minsplits, variables, threshold, seed) {
  set.seed(seed)
  performance_metrics <- data.frame(Fold = integer(), TP = integer(), TN = integer(),
                                    FP = integer(), FN = integer(), AUC = numeric(),
                                    Best_CP = numeric(), Best_Minsplit = integer(), 
                                    Best_Variables = character(), stringsAsFactors = FALSE)
  
  inner_fold_metrics <- data.frame(Fold = integer(), Inner_Fold = integer(), 
                                   TP = integer(), TN = integer(), FP = integer(), 
                                   FN = integer(), AUC = numeric(), stringsAsFactors = FALSE)
  
  data[[target_name]] <- as.factor(data[[target_name]])
  outer_folds <- createFolds(data[[target_name]], k = outer)
  
  for (outer_index in seq_along(outer_folds)) {
    outer_train_data <- data[-outer_folds[[outer_index]], ]
    outer_test_data <- data[outer_folds[[outer_index]], ]
    
    cat(sprintf("Outer Fold %d\n", outer_index))
    cat("Training Data - Class 0: ", sum(outer_train_data[[target_name]] == "0"), "\n")
    cat("Training Data - Class 1: ", sum(outer_train_data[[target_name]] == "1"), "\n")
    cat("Testing Data - Class 0: ", sum(outer_test_data[[target_name]] == "0"), "\n")
    cat("Testing Data - Class 1: ", sum(outer_test_data[[target_name]] == "1"), "\n")
    
    inner_folds <- createFolds(outer_train_data[[target_name]], k = inner)
    best_model <- NULL
    best_auc <- 0
    best_cp <- NULL
    best_minsplit <- NULL
    best_vars <- variables
    
    grid <- expand.grid(cp = cps, minsplit = minsplits)
    
    stepAUC <- function(vars, data, target_name, grid, inner_folds) {
      current_vars <- vars
      best_inner_auc <- 0
      best_inner_model <- NULL
      best_inner_cp <- NULL
      best_inner_minsplit <- NULL
      best_inner_vars <- current_vars
      improved <- TRUE
    
      while (improved && length(current_vars) > 1) {
        improved <- FALSE
        best_auc_in_step <- best_inner_auc
    
        for (var in current_vars) {
          temp_vars <- setdiff(current_vars, var)
          inner_auc <- numeric()
    
          for (params in 1:nrow(grid)) {
            auc_vals <- numeric()
    
            for (inner_index in seq_along(inner_folds)) {
              inner_train_data <- data[-inner_folds[[inner_index]], ]
              inner_test_data <- data[inner_folds[[inner_index]], ]
    
              result <- tryCatch({
                model <- rpart(as.formula(paste(target_name, "~", paste(temp_vars, collapse = "+"))),
                               data = inner_train_data, control = rpart.control(cp = grid$cp[params], minsplit = grid$minsplit[params]))
    
                predictions <- predict(model, inner_test_data, type = "prob")[, 2]
                roc_curve <- roc(inner_test_data[[target_name]], predictions)
                auc_vals[inner_index] <- roc_curve$auc
    
              }, error = function(e) {
                cat("Error with variables:", paste(temp_vars, collapse = ", "), "\n")
                auc_vals[inner_index] <- 0
              })
            }
    
            if (mean(auc_vals, na.rm = TRUE) > best_auc_in_step) {
              best_auc_in_step <- mean(auc_vals, na.rm = TRUE)
              best_inner_model <- model
              best_inner_cp <- grid$cp[params]
              best_inner_minsplit <- grid$minsplit[params]
              best_inner_vars <- temp_vars
              improved <- TRUE
            }
          }
        }
    
        if (improved) {
          current_vars <- best_inner_vars
          best_inner_auc <- best_auc_in_step
          cat("Vars improved: ", paste(current_vars, collapse=", "), "\n")
          cat("AUC: ", paste(best_inner_auc), "\n")
        }
      }
    
      list(model = best_inner_model, auc = best_inner_auc, 
           cp = best_inner_cp, minsplit = best_inner_minsplit, 
           vars = best_inner_vars)
    }

    
    best_result <- stepAUC(variables, outer_train_data, target_name, grid, inner_folds)
    best_model <- best_result$model
    best_auc <- best_result$auc
    best_cp <- best_result$cp
    best_minsplit <- best_result$minsplit
    best_vars <- best_result$vars

    # Evaluar el mejor modelo en el validation split del inner loop
    inner_auc_vals <- numeric()
    for (inner_index in seq_along(inner_folds)) {
      inner_train_data <- outer_train_data[-inner_folds[[inner_index]], ]
      inner_test_data <- outer_train_data[inner_folds[[inner_index]], ]
    
      predictions <- predict(best_model, inner_test_data, type = "prob")[, 2]
      roc_curve <- roc(inner_test_data[[target_name]], predictions)
      inner_auc_vals[inner_index] <- roc_curve$auc
    
      pred_class <- ifelse(predictions > threshold, "1", "0")
      confusion <- table(Actual = inner_test_data[[target_name]], Predicted = pred_class)
      TP <- ifelse(nrow(confusion) >= 2 && ncol(confusion) >= 2, confusion[2, 2], 0)
      TN <- ifelse(nrow(confusion) >= 1 && ncol(confusion) >= 1, confusion[1, 1], 0)
      FP <- ifelse(nrow(confusion) >= 1 && ncol(confusion) >= 2, confusion[1, 2], 0)
      FN <- ifelse(nrow(confusion) >= 2 && ncol(confusion) >= 1, confusion[2, 1], 0)
      auc <- roc_curve$auc
    
      inner_metrics <- setNames(c(outer_index, inner_index, TP, TN, FP, FN, auc), names(inner_fold_metrics))
      inner_fold_metrics <- rbind(inner_fold_metrics, as.data.frame(t(inner_metrics)))
    }

    
    best_inner_auc <- mean(inner_auc_vals, na.rm = TRUE)
    # Reentrenar el mejor modelo con todos los datos de outer_train_data
    best_model <- rpart(as.formula(paste(target_name, "~", paste(best_vars, collapse = "+"))),
                    data = outer_train_data, control = rpart.control(cp = best_cp, minsplit = best_minsplit))

    predictions <- predict(best_model, outer_test_data, type = "prob")[, 2]
    pred_class <- ifelse(predictions > threshold, "1", "0")
    confusion <- table(Actual = outer_test_data[[target_name]], Predicted = pred_class)
    
    if (nrow(confusion) < 2 || ncol(confusion) < 2) {
      TP <- ifelse(nrow(confusion) >= 2 && ncol(confusion) >= 2, confusion[2, 2], 0)
      TN <- ifelse(nrow(confusion) >= 1 && ncol(confusion) >= 1, confusion[1, 1], 0)
      FP <- ifelse(nrow(confusion) >= 1 && ncol(confusion) >= 2, confusion[1, 2], 0)
      FN <- ifelse(nrow(confusion) >= 2 && ncol(confusion) >= 1, confusion[2, 1], 0)
    } else {
      TP <- confusion[2, 2]
      TN <- confusion[1, 1]
      FP <- confusion[1, 2]
      FN <- confusion[2, 1]
    }
    auc <- roc(outer_test_data[[target_name]], predictions)$auc
    
    metrics <- setNames(c(outer_index, TP, TN, FP, FN, auc, 
                          best_cp, best_minsplit, paste(best_vars, collapse = ",")), 
                        names(performance_metrics))
    performance_metrics <- rbind(performance_metrics, as.data.frame(t(metrics)))
    
    cat("Performance:", paste(performance_metrics))
    cat(sprintf("Confusion Matrix for Fold %d:\n", outer_index))
    print(confusion)
    cat(sprintf("Metrics for Fold %d:\n", outer_index))
    cat(sprintf("TP: %d, TN: %d, FP: %d, FN: %d, AUC: %f, Best CP: %f, Best Minsplit: %d, Best Variables: %s\n",
                TP, TN, FP, FN, auc, best_cp, best_minsplit, paste(best_vars, collapse = ",")))
  }
  
  return(list(performance_metrics = performance_metrics, inner_fold_metrics = inner_fold_metrics))
}

evaluate_with_seeds_DecisionTree_StepAUC <- 
  function(evaluation_function, data, target_var, vars, threshold, seeds, outer, inner, cp_values, minsplit_values) {
  results_list <- list()

  for (seed in seeds) {
    set.seed(seed)
    cat("Evaluating with seed:", seed, "\n")

    results <- evaluation_function(data = data, 
                                   target_name = target_var, 
                                   variables = vars, 
                                   threshold = threshold, 
                                   outer = outer, 
                                   inner = inner, 
                                   seed = seed,
                                   cps = cp_values,
                                   minsplits = minsplit_values)
    results_list[[paste0("seed_", seed)]] <- results
  }

  return(results_list)
  }

performance_rpart_StepAUC <- evaluate_with_seeds_DecisionTree_StepAUC(
  evaluation_function = double_cross_validation_rpart_stepAUC,
  data = data_factor,
  outer = 5,
  inner = 2,
  target_var = "PCR",
  cp_values = c(0.001, 0.01, 0.1),
  minsplit_values = c(2, 5, 10, 15, 20),
  vars = todas_variables_factor,
  threshold = .35,
  seeds = seeds
)
performance_rpart_StepAUC

load("RData_Files_Algorithms/DecisionTree_StepAUC.RData")

3.4.5 Embedded (Lasso)

double_cross_validation_rpart_Lasso <- function(data, target_name, outer, inner, cps, minsplits, variables, threshold, seed) {
  set.seed(seed)
  performance_metrics <- data.frame(Fold = integer(), TP = integer(), TN = integer(),
                                    FP = integer(), FN = integer(), AUC = numeric(),
                                    Best_cp = numeric(), Best_minsplit = numeric(), 
                                    Best_Variables = character(), stringsAsFactors = FALSE)
  
  inner_fold_metrics <- data.frame(Fold = integer(), Inner_Fold = integer(), 
                                   TP = integer(), TN = integer(), FP = integer(), 
                                   FN = integer(), AUC = numeric(), stringsAsFactors = FALSE)
  
  # Convertir las variables del dataset con model.matrix
  target <- as.numeric(data[[target_name]]) - 1
  features <- data[, setdiff(names(data), target_name)]
  data_matrix <- as.data.frame(model.matrix(~ ., data=features)[,-1])  # Eliminar el intercepto
  data_matrix[[target_name]] <- target
  
  outer_folds <- createFolds(data_matrix[[target_name]], k = outer)
  
  for (outer_index in seq_along(outer_folds)) {
    outer_train_data <- data_matrix[-outer_folds[[outer_index]], ]
    outer_test_data <- data_matrix[outer_folds[[outer_index]], ]
    
    cat(sprintf("Outer Fold %d\n", outer_index))
    cat("Training Data - Class 0: ", sum(outer_train_data[[target_name]] == 0), "\n")
    cat("Training Data - Class 1: ", sum(outer_train_data[[target_name]] == 1), "\n")
    cat("Testing Data - Class 0: ", sum(outer_test_data[[target_name]] == 0), "\n")
    cat("Testing Data - Class 1: ", sum(outer_test_data[[target_name]] == 1), "\n")
    
    inner_folds <- createFolds(outer_train_data[[target_name]], k = inner)
    best_model <- NULL
    best_auc <- 0
    best_cp <- NULL
    best_minsplit <- NULL
    best_vars <- colnames(outer_train_data)[!colnames(outer_train_data) %in% target_name]
    
    grid <- expand.grid(cp = cps, minsplit = minsplits)
    
    select_vars_lasso <- function(data, target_name) {
      target <- as.numeric(data[[target_name]])
      features <- data[, setdiff(names(data), target_name)]
      X <- model.matrix(~ ., data=features)[,-1]  # Eliminamos el intercepto
      y <- target
      
      # Ajustamos el modelo Lasso
      lasso_model <- cv.glmnet(X, y, alpha=1)
      lasso_coef <- coef(lasso_model, s = "lambda.min")
      
      # Seleccionamos las variables no nulas
      selected_vars <- rownames(lasso_coef)[lasso_coef[, 1] != 0]
      selected_vars <- setdiff(selected_vars, "(Intercept)")
      #cat("Selected variables by Lasso: ", paste(selected_vars, collapse=", "), "\n")
      return(selected_vars)
    }
    
    stepAUC <- function(vars, data, target_name, grid, inner_folds) {
      best_inner_auc <- 0
      best_inner_model <- NULL
      best_inner_cp <- NULL
      best_inner_minsplit <- NULL
      best_inner_vars <- vars
      improved <- TRUE
      
      while (improved && length(vars) > 1) {
        improved <- FALSE
        best_auc_in_step <- best_inner_auc
        
        for (params in 1:nrow(grid)) {
          auc_vals <- numeric()
          
          for (inner_index in seq_along(inner_folds)) {
            inner_train_data <- outer_train_data[-inner_folds[[inner_index]], ]
            inner_test_data <- outer_train_data[inner_folds[[inner_index]], ]
            
            selected_vars <- select_vars_lasso(inner_train_data, target_name)
            formula <- as.formula(paste(target_name, "~", paste(selected_vars, collapse = "+")))
            #cat("Testing with formula: ", deparse(formula), "\n")
            #cat("Grid parameters - cp: ", grid$cp[params], ", minsplit: ", grid$minsplit[params], "\n")
            
            result <- tryCatch({
              model <- rpart(formula, data = inner_train_data, method = "class", 
                             control = rpart.control(cp = grid$cp[params], minsplit = grid$minsplit[params]))
              
              predictions <- predict(model, inner_test_data, type = "prob")[,2]
              roc_curve <- roc(inner_test_data[[target_name]], predictions)
              auc_vals[inner_index] <- roc_curve$auc
              #cat("AUC for inner fold ", inner_index, ": ", roc_curve$auc, "\n")
              
            }, error = function(e) {
              cat("Error with variables: ", paste(selected_vars, collapse = ", "), "\n")
              cat("Error message: ", e$message, "\n")
              auc_vals[inner_index] <- NA  # Set to NA or some other indicator of failure
            })
          }
          
          mean_auc <- mean(auc_vals, na.rm = TRUE)
          #cat("Mean AUC for params - cp: ", grid$cp[params], ", minsplit: ", grid$minsplit[params], ": ", mean_auc, "\n")
          if (!is.na(mean_auc) && mean_auc > best_auc_in_step) {
            best_auc_in_step <- mean_auc
            best_inner_model <- model
            best_inner_cp <- grid$cp[params]
            best_inner_minsplit <- grid$minsplit[params]
            best_inner_vars <- selected_vars
            improved <- TRUE
          }
        }
        
        if (improved) {
          vars <- best_inner_vars
          best_inner_auc <- best_auc_in_step
          cat("Vars improved: ", paste(vars, collapse=", "), "\n")
          cat("AUC: ", best_inner_auc, "\n")
        }
      }
      
      list(model = best_inner_model, auc = best_inner_auc, 
           cp = best_inner_cp, minsplit = best_inner_minsplit, 
           vars = best_inner_vars)
    }
    
    best_result <- stepAUC(best_vars, outer_train_data, target_name, grid, inner_folds)
    best_model <- best_result$model
    best_auc <- best_result$auc
    best_cp <- best_result$cp
    best_minsplit <- best_result$minsplit
    best_vars <- best_result$vars
    
    if (is.null(best_model)) {
      cat("Error: No valid model found for outer fold ", outer_index, "\n")
      next
    }

    # Evaluar el mejor modelo en el validation split del inner loop
    inner_auc_vals <- numeric()
    for (inner_index in seq_along(inner_folds)) {
      inner_train_data <- outer_train_data[-inner_folds[[inner_index]], ]
      inner_test_data <- outer_train_data[inner_folds[[inner_index]], ]
      
      predictions <- predict(best_model, inner_test_data, type = "prob")[,2]
      roc_curve <- roc(inner_test_data[[target_name]], predictions)
      inner_auc_vals[inner_index] <- roc_curve$auc
      #cat("Inner fold ", inner_index, " AUC: ", roc_curve$auc, "\n")

      pred_class <- ifelse(predictions > threshold, 1, 0)
      confusion <- table(Actual = inner_test_data[[target_name]], Predicted = pred_class)
      TP <- ifelse(nrow(confusion) >= 2 && ncol(confusion) >= 2, confusion[2, 2], 0)
      TN <- ifelse(nrow(confusion) >= 1 && ncol(confusion) >= 1, confusion[1, 1], 0)
      FP <- ifelse(nrow(confusion) >= 1 && ncol(confusion) >= 2, confusion[1, 2], 0)
      FN <- ifelse(nrow(confusion) >= 2 && ncol(confusion) >= 1, confusion[2, 1], 0)
      auc <- roc_curve$auc
      
      inner_metrics <- setNames(c(outer_index, inner_index, TP, TN, FP, FN, auc), names(inner_fold_metrics))
      inner_fold_metrics <- rbind(inner_fold_metrics, as.data.frame(t(inner_metrics)))
    }
    
    best_inner_auc <- mean(inner_auc_vals, na.rm = TRUE)
    cat("Best inner AUC: ", best_inner_auc, "\n")
    best_model <- rpart(as.formula(paste(target_name, "~", paste(best_vars, collapse = "+"))),
                        data = outer_train_data, method = "class", 
                        control = rpart.control(cp = best_cp, minsplit = best_minsplit))
    predictions <- predict(best_model, outer_test_data, type = "prob")[,2]
    pred_class <- ifelse(predictions > threshold, 1, 0)  # Asignar clases basadas en umbral
    confusion <- table(Actual = outer_test_data[[target_name]], Predicted = pred_class)
    
    if (nrow(confusion) < 2 || ncol(confusion) < 2) {
      TP <- ifelse(nrow(confusion) >= 2 && ncol(confusion) >= 2, confusion[2, 2], 0)
      TN <- ifelse(nrow(confusion) >= 1 && ncol(confusion) >= 1, confusion[1, 1], 0)
      FP <- ifelse(nrow(confusion) >= 1 && ncol(confusion) >= 2, confusion[1, 2], 0)
      FN <- ifelse(nrow(confusion) >= 2 && ncol(confusion) >= 1, confusion[2, 1], 0)
    } else {
      TP <- confusion[2, 2]
      TN <- confusion[1, 1]
      FP <- confusion[1, 2]
      FN <- confusion[2, 1]
    }
    auc <- roc(outer_test_data[[target_name]], predictions)$auc
    
    metrics <- setNames(c(outer_index, TP, TN, FP, FN, auc, 
                          best_cp, best_minsplit, paste(best_vars, collapse = ",")), 
                        names(performance_metrics))
    performance_metrics <- rbind(performance_metrics, as.data.frame(t(metrics)))
    
    cat("Performance:\n", paste(performance_metrics), "\n")
    cat(sprintf("Confusion Matrix for Fold %d:\n", outer_index))
    print(confusion)
    cat(sprintf("Metrics for Fold %d:\n", outer_index))
    cat(sprintf("TP: %d, TN: %d, FP: %d, FN: %d, AUC: %f, Best cp: %f, Best minsplit: %d, Best Variables: %s\n",
                TP, TN, FP, FN, auc, best_cp, best_minsplit, paste(best_vars, collapse = ",")))
  }
  
  list(performance_metrics = performance_metrics, inner_fold_metrics = inner_fold_metrics)
}


evaluate_with_seeds_DecisionTree_lasso <- 
  function(evaluation_function, data, target_var, vars, threshold, seeds, outer, inner, cp_values, minsplit_values) {
  results_list <- list()

  for (seed in seeds) {
    set.seed(seed)
    cat("Evaluating with seed:", seed, "\n")

    results <- evaluation_function(data = data, 
                                   target_name = target_var, 
                                   variables = vars, 
                                   threshold = threshold, 
                                   outer = outer, 
                                   inner = inner, 
                                   seed = seed,
                                   cps = cp_values,
                                   minsplits = minsplit_values)
    results_list[[paste0("seed_", seed)]] <- results
  }

  return(results_list)
  }

performance_rpart_lasso <- evaluate_with_seeds_DecisionTree_lasso(
  evaluation_function = double_cross_validation_rpart_Lasso,
  data = data_factor,
  outer = 5,
  inner = 2,
  target_var = "PCR",
  cp_values = c(0.001, 0.01, 0.1),
  minsplit_values = c(2, 5, 10, 15, 20),
  vars = todas_variables_factor,
  threshold = .35,
  seeds = seeds
)

load("RData_Files_Algorithms/DecisionTree_Lasso.RData")

3.5 K-Nearest Neighbours

3.5.1 Aparente

evaluate_aparent_performance_model_knn <- function(data, target_var, model_func, vars = NULL, threshold = 0.5) {
  if (is.null(vars)) {
    vars <- setdiff(names(data), target_var)
  }
  
  data[[target_var]] <- factor(data[[target_var]], levels = c(0, 1))
  formula <- as.formula(paste(target_var, "~", paste(vars, collapse = "+")))
  
  model <- model_func(formula, data)
  train_x <- model$train_x
  train_y <- model$train_y
  test_x <- model.matrix(formula, data)[, -1]
  predictions <- knn(train = train_x, test = test_x, cl = train_y, k = model$k, prob = TRUE)
  predictions <- attr(predictions, "prob")
  
  actual_classes <- data[[target_var]]
  predicted_classes <- ifelse(predictions > threshold, 1, 0)
  confusion_matrix <- table(predicted_classes, actual_classes)
  
  tp <- ifelse("1" %in% rownames(confusion_matrix) && "1" %in% colnames(confusion_matrix), confusion_matrix["1", "1"], 0)
  tn <- ifelse("0" %in% rownames(confusion_matrix) && "0" %in% colnames(confusion_matrix), confusion_matrix["0", "0"], 0)
  fp <- ifelse("1" %in% rownames(confusion_matrix) && "0" %in% colnames(confusion_matrix), confusion_matrix["1", "0"], 0)
  fn <- ifelse("0" %in% rownames(confusion_matrix) && "1" %in% colnames(confusion_matrix), confusion_matrix["0", "1"], 0)
  
  accuracy <- (tp + tn) / (tp + tn + fp + fn)
  precision <- ifelse(tp + fp > 0, tp / (tp + fp), 0)
  recall <- ifelse(tp + fn > 0, tp / (tp + fn), 0)
  f1_score <- ifelse(precision + recall > 0, 2 * (precision * recall) / (precision + recall), 0)
  roc_obj <- pROC::roc(actual_classes, predictions)
  
  list(confusion_matrix = confusion_matrix, accuracy = accuracy, precision = precision, 
       recall = recall, f1_score = f1_score, roc_curve = roc_obj)
}


knn_model <- function(formula, data) {
  train_x <- model.matrix(formula, data)[, -1]
  train_y <- data[[all.vars(formula)[1]]]
  list(train_x = train_x, train_y = train_y, k = 5)  
}

resultadosAparentesKNN <- evaluate_aparent_performance_model_knn(data = data_numeric, target_var = "PCR",
                                                             model_func = knn_model,
                                                             vars = ".",
                                                             threshold = .35)

load("RData_Files_Algorithms/KNN_Aparente.RData")

3.5.2 Todas las Variables

double_cross_validation_knn <- function(data, target_name, outer, inner, ks, ls, variables, threshold, seed) {
  set.seed(seed)
  performance_metrics <- data.frame(Fold = integer(), TP = integer(), TN = integer(),
                                    FP = integer(), FN = integer(), AUC = numeric(),
                                    Best_k = integer(), Best_l = numeric(), stringsAsFactors = FALSE)
  
  data[[target_name]] <- as.factor(data[[target_name]])
  outer_folds <- createFolds(data[[target_name]], k = outer, list = TRUE, returnTrain = FALSE)
  
  for (outer_index in seq_along(outer_folds)) {
    outer_train_data <- data[-outer_folds[[outer_index]], ]
    outer_test_data <- data[outer_folds[[outer_index]], ]
    
    cat(sprintf("Outer Fold %d\n", outer_index))
    cat("Training Data - Class 0: ", sum(outer_train_data[[target_name]] == "0"), "\n")
    cat("Training Data - Class 1: ", sum(outer_train_data[[target_name]] == "1"), "\n")
    cat("Testing Data - Class 0: ", sum(outer_test_data[[target_name]] == "0"), "\n")
    cat("Testing Data - Class 1: ", sum(outer_test_data[[target_name]] == "1"), "\n")
    
    if (sum(outer_train_data[[target_name]] == "0") == 0 || sum(outer_train_data[[target_name]] == "1") == 0) {
      cat("Skipping outer fold due to lack of class diversity\n")
      next
    }
    
    inner_folds <- createFolds(outer_train_data[[target_name]], k = inner, list = TRUE, returnTrain = FALSE)
    best_model <- NULL
    best_auc <- 0
    best_k <- NULL
    best_l <- NULL
    
    grid <- expand.grid(k = ks, l = ls)
    formula <- reformulate(variables, target_name)

    for (params in 1:nrow(grid)) {
      inner_auc <- numeric()

      for (inner_index in seq_along(inner_folds)) {
        inner_train_data <- outer_train_data[-inner_folds[[inner_index]], ]
        inner_test_data <- outer_train_data[inner_folds[[inner_index]], ]

        if (sum(inner_train_data[[target_name]] == "0") == 0 || sum(inner_train_data[[target_name]] == "1") == 0) {
          cat("Skipping inner fold due to lack of class diversity\n")
          next
        }

        train_x <- inner_train_data[, variables, drop = FALSE]
        test_x <- inner_test_data[, variables, drop = FALSE]
        train_y <- inner_train_data[[target_name]]
        test_y <- inner_test_data[[target_name]]
        
        result <- tryCatch({
          predictions <- knn(train = train_x, test = test_x, cl = train_y, k = grid$k[params], prob = TRUE)
          probabilities <- attr(predictions, "prob")
          probabilities <- ifelse(predictions == "1", probabilities, 1 - probabilities)
          roc_curve <- roc(as.numeric(as.character(test_y)), probabilities)
          inner_auc[inner_index] <- roc_curve$auc
        }, error = function(e) {
          cat("Error with parameters: k =", grid$k[params], ", l =", grid$l[params], "\n")
          cat("Error message: ", e$message, "\n")
          inner_auc[inner_index] <- 0
        })
      }

      if (mean(inner_auc, na.rm = TRUE) > best_auc) {
        best_auc <- mean(inner_auc, na.rm = TRUE)
        best_model <- list(train_x = train_x, train_y = train_y, k = grid$k[params], l = grid$l[params])
        best_k <- grid$k[params]
        best_l <- grid$l[params]
      }
    }

    if (is.null(best_model)) {
      cat("No valid model found for outer fold ", outer_index, "\n")
      next
    }

    # Reentrenar el mejor modelo con todos los datos de outer_train_data
    train_x <- outer_train_data[, variables, drop = FALSE]
    train_y <- outer_train_data[[target_name]]
    
    # Evaluar el mejor modelo en el conjunto de prueba del outer fold
    test_x <- outer_test_data[, variables, drop = FALSE]
    test_y <- outer_test_data[[target_name]]
    predictions <- knn(train = train_x, test = test_x, cl = train_y, k = best_k, prob = TRUE)
    probabilities <- attr(predictions, "prob")
    probabilities <- ifelse(predictions == "1", probabilities, 1 - probabilities)
    pred_class <- ifelse(probabilities > threshold, "1", "0")
    confusion <- table(Actual = test_y, Predicted = pred_class)
    
    if (nrow(confusion) < 2 || ncol(confusion) < 2) {
      TP <- ifelse(nrow(confusion) >= 2 && ncol(confusion) >= 2, confusion[2, 2], 0)
      TN <- ifelse(nrow(confusion) >= 1 && ncol(confusion) >= 1, confusion[1, 1], 0)
      FP <- ifelse(nrow(confusion) >= 1 && ncol(confusion) >= 2, confusion[1, 2], 0)
      FN <- ifelse(nrow(confusion) >= 2 && ncol(confusion) >= 1, confusion[2, 1], 0)
    } else {
      TP <- confusion[2, 2]
      TN <- confusion[1, 1]
      FP <- confusion[1, 2]
      FN <- confusion[2, 1]
    }
    auc <- roc(test_y, probabilities)$auc
    
    metrics <- setNames(c(outer_index, TP, TN, FP, FN, auc, best_k, best_l), names(performance_metrics))
    performance_metrics <- rbind(performance_metrics, as.data.frame(t(metrics)))
    
    cat("Performance:", paste(performance_metrics))
    cat(sprintf("Confusion Matrix for Fold %d:\n", outer_index))
    print(confusion)
    cat(sprintf("Metrics for Fold %d:\n", outer_index))
    cat(sprintf("TP: %d, TN: %d, FP: %d, FN: %d, AUC: %f, Best k: %d, Best l: %f\n",
                TP, TN, FP, FN, auc, best_k, best_l))
  }
  
  return(performance_metrics)
}
evaluate_with_seeds_KNN_TodasVariables_Asociacion <- 
  function(evaluation_function, data, target_var, vars, threshold, seeds, outer, inner, ks, ls) {
  results_list <- list()

  for (seed in seeds) {
    set.seed(seed)
    cat("Evaluating with seed:", seed, "\n")

    results <- evaluation_function(data = data, 
                                   target_name = target_var, 
                                   variables = vars, 
                                   threshold = threshold, 
                                   outer = outer, 
                                   inner = inner, 
                                   seed = seed,
                                   ks = ks,
                                   ls = ls)
    results_list[[paste0("seed_", seed)]] <- results
  }

  return(results_list)
  }

performance_knn_todasvariables <- evaluate_with_seeds_KNN_TodasVariables_Asociacion(
  evaluation_function = double_cross_validation_knn,
  data = data_numeric,
  outer = 5,
  inner = 2,
  target_var = "PCR",
  ks = c(10, 16, 20, 25, 30, 32),
  ls = c(0, 1, 2),
  vars = todas_variables_numeric,
  threshold = .35,
  seeds = seeds
)

load("RData_Files_Algorithms/KNN_TodasVariables.RData")

3.5.3 Filtrado (Asociación)

performance_knn_Asociacion <- evaluate_with_seeds_KNN_TodasVariables_Asociacion(
  evaluation_function = double_cross_validation_knn,
  data = data_numeric,
  outer = 5,
  inner = 2,
  target_var = "PCR",
  ks = c(10, 16, 20, 25, 30, 32),
  ls = c(0, 1, 2),
  vars = asociacion_numeric,
  threshold = .35,
  seeds = seeds
)

load("RData_Files_Algorithms/KNN_Asociacion.RData")

3.5.4 Wrapped (StepAUC)

double_cross_validation_knn_stepAUC <- function(data, target_name, outer, inner, ks, ls, variables, threshold, seed) {
  set.seed(seed)
  performance_metrics <- data.frame(Fold = integer(), TP = integer(), TN = integer(),
                                    FP = integer(), FN = integer(), AUC = numeric(),
                                    Best_k = integer(), Best_l = numeric(), Best_Variables = character(), stringsAsFactors = FALSE)
  
  inner_fold_metrics <- data.frame(Fold = integer(), Inner_Fold = integer(), 
                                   TP = integer(), TN = integer(), FP = integer(), 
                                   FN = integer(), AUC = numeric(), stringsAsFactors = FALSE)
  
  data[[target_name]] <- as.factor(data[[target_name]])
  outer_folds <- createFolds(data[[target_name]], k = outer, list = TRUE, returnTrain = FALSE)
  
  for (outer_index in seq_along(outer_folds)) {
    outer_train_data <- data[-outer_folds[[outer_index]], ]
    outer_test_data <- data[outer_folds[[outer_index]], ]
    
    cat(sprintf("Outer Fold %d\n", outer_index))
    cat("Training Data - Class 0: ", sum(outer_train_data[[target_name]] == "0"), "\n")
    cat("Training Data - Class 1: ", sum(outer_train_data[[target_name]] == "1"), "\n")
    cat("Testing Data - Class 0: ", sum(outer_test_data[[target_name]] == "0"), "\n")
    cat("Testing Data - Class 1: ", sum(outer_test_data[[target_name]] == "1"), "\n")
    
    inner_folds <- createFolds(outer_train_data[[target_name]], k = inner, list = TRUE, returnTrain = FALSE)
    best_model <- NULL
    best_auc <- 0
    best_k <- NULL
    best_l <- NULL
    best_vars <- variables
    
    grid <- expand.grid(k = ks, l = ls)
    
    stepAUC <- function(vars, outer_train_data, target_name, grid, inner_folds) {
      current_vars <- vars
      best_inner_auc <- 0
      best_inner_model <- NULL
      best_inner_k <- NULL
      best_inner_l <- NULL
      best_inner_vars <- current_vars
      improved <- TRUE
      
      while (improved && length(current_vars) > 1) {
        improved <- FALSE
        best_auc_in_step <- best_inner_auc
        
        for (var in current_vars) {
          temp_vars <- setdiff(current_vars, var)
          inner_auc <- numeric()
          
          for (params in 1:nrow(grid)) {
            auc_vals <- numeric()
            
            for (inner_index in seq_along(inner_folds)) {
              inner_train_data <- outer_train_data[-inner_folds[[inner_index]], ]
              inner_test_data <- outer_train_data[inner_folds[[inner_index]], ]
              
              train_x <- inner_train_data[, temp_vars, drop = FALSE]
              test_x <- inner_test_data[, temp_vars, drop = FALSE]
              train_y <- inner_train_data[[target_name]]
              test_y <- inner_test_data[[target_name]]
              
              result <- tryCatch({
                predictions <- knn(train = train_x, test = test_x, cl = train_y, k = grid$k[params], l = grid$l[params], prob = TRUE)
                probabilities <- attr(predictions, "prob")
                roc_curve <- roc(as.numeric(as.character(test_y)), probabilities)
                auc_vals[inner_index] <- roc_curve$auc
                
              }, error = function(e) {
                cat("Error with variables:", paste(temp_vars, collapse = ", "), "\n")
                cat("Error message: ", e$message, "\n")
                auc_vals[inner_index] <- 0  # Set to NA or some other indicator of failure
              })
            }
            
            if (mean(auc_vals, na.rm = TRUE) > best_auc_in_step) {
              best_auc_in_step <- mean(auc_vals, na.rm = TRUE)
              best_inner_model <- list(train_x = train_x, train_y = train_y, k = grid$k[params], l = grid$l[params])
              best_inner_k <- grid$k[params]
              best_inner_l <- grid$l[params]
              best_inner_vars <- temp_vars
              improved <- TRUE
            }
          }
        }
        
        if (improved) {
          current_vars <- best_inner_vars
          best_inner_auc <- best_auc_in_step
          cat("Vars improved: ", paste(current_vars, collapse=", "), "\n")
          cat("AUC: ", paste(best_inner_auc), "\n")
        }
      }
      
      list(model = best_inner_model, auc = best_inner_auc, 
           k = best_inner_k, l = best_inner_l, 
           vars = best_inner_vars)
      
    }
    
    # Después de seleccionar el mejor modelo y parámetros en los inner folds
    best_result <- stepAUC(variables, outer_train_data, target_name, grid, inner_folds)
    best_model <- best_result$model
    best_auc <- best_result$auc
    best_k <- best_result$k
    best_l <- best_result$l
    best_vars <- best_result$vars

    # Evaluar el mejor modelo en el validation split del inner loop
    inner_auc_vals <- numeric()
    for (inner_index in seq_along(inner_folds)) {
      inner_train_data <- outer_train_data[-inner_folds[[inner_index]], ]
      inner_test_data <- outer_train_data[inner_folds[[inner_index]], ]
      
      train_x <- best_model$train_x
      test_x <- inner_test_data[, best_vars, drop = FALSE]
      train_y <- best_model$train_y
      test_y <- inner_test_data[[target_name]]
      
      predictions <- knn(train = train_x, test = test_x, cl = train_y, k = best_k, l = best_l, prob = TRUE)
      probabilities <- attr(predictions, "prob")
      roc_curve <- roc(as.numeric(as.character(test_y)), probabilities)
      inner_auc_vals[inner_index] <- roc_curve$auc

      pred_class <- ifelse(probabilities > threshold, "1", "0")
      confusion <- table(Actual = test_y, Predicted = pred_class)
      TP <- ifelse(nrow(confusion) >= 2 && ncol(confusion) >= 2, confusion[2, 2], 0)
      TN <- ifelse(nrow(confusion) >= 1 && ncol(confusion) >= 1, confusion[1, 1], 0)
      FP <- ifelse(nrow(confusion) >= 1 && ncol(confusion) >= 2, confusion[1, 2], 0)
      FN <- ifelse(nrow(confusion) >= 2 && ncol(confusion) >= 1, confusion[2, 1], 0)
      auc <- roc_curve$auc
      
      inner_metrics <- setNames(c(outer_index, inner_index, TP, TN, FP, FN, auc), names(inner_fold_metrics))
      inner_fold_metrics <- rbind(inner_fold_metrics, as.data.frame(t(inner_metrics)))
    }
    
    best_inner_auc <- mean(inner_auc_vals, na.rm = TRUE)
    
    # Reentrenar el mejor modelo con todos los datos de outer_train_data
    train_x <- outer_train_data[, best_vars, drop = FALSE]
    train_y <- outer_train_data[[target_name]]
    best_model <- list(train_x = train_x, train_y = train_y, k = best_k, l = best_l)
    
    # Evaluar el mejor modelo en el conjunto de prueba del outer fold
    test_x <- outer_test_data[, best_vars, drop = FALSE]
    test_y <- outer_test_data[[target_name]]
    predictions <- knn(train = best_model$train_x, test = test_x, cl = best_model$train_y, k = best_k, l = best_l, prob = TRUE)
    probabilities <- attr(predictions, "prob")
    pred_class <- ifelse(probabilities > threshold, "1", "0")
    confusion <- table(Actual = test_y, Predicted = pred_class)
    
    if (nrow(confusion) < 2 || ncol(confusion) < 2) {
      TP <- ifelse(nrow(confusion) >= 2 && ncol(confusion) >= 2, confusion[2, 2], 0)
      TN <- ifelse(nrow(confusion) >= 1 && ncol(confusion) >= 1, confusion[1, 1], 0)
      FP <- ifelse(nrow(confusion) >= 1 && ncol(confusion) >= 2, confusion[1, 2], 0)
      FN <- ifelse(nrow(confusion) >= 2 && ncol(confusion) >= 1, confusion[2, 1], 0)
    } else {
      TP <- confusion[2, 2]
      TN <- confusion[1, 1]
      FP <- confusion[1, 2]
      FN <- confusion[2, 1]
    }
    auc <- roc(test_y, probabilities)$auc
    
    metrics <- setNames(c(outer_index, TP, TN, FP, FN, auc, 
                          best_k, best_l, paste(best_vars, collapse = ",")), 
                        names(performance_metrics))
    performance_metrics <- rbind(performance_metrics, as.data.frame(t(metrics)))
    
    cat("Performance:", paste(performance_metrics))
    cat(sprintf("Confusion Matrix for Fold %d:\n", outer_index))
    print(confusion)
    cat(sprintf("Metrics for Fold %d:\n", outer_index))
    cat(sprintf("TP: %d, TN: %d, FP: %d, FN: %d, AUC: %f, Best k: %d, Best l: %f, Best Variables: %s\n",
                TP, TN, FP, FN, auc, best_k, best_l, paste(best_vars, collapse = ",")))
  }
  
  list(performance_metrics = performance_metrics, inner_fold_metrics = inner_fold_metrics)
}
evaluate_with_seeds_KNN_stepAUC <- 
  function(evaluation_function, data, target_var, vars, threshold, seeds, outer, inner, ks, ls) {
  results_list <- list()

  for (seed in seeds) {
    set.seed(seed)
    cat("Evaluating with seed:", seed, "\n")

    results <- evaluation_function(data = data, 
                                   target_name = target_var, 
                                   variables = vars, 
                                   threshold = threshold, 
                                   outer = outer, 
                                   inner = inner, 
                                   seed = seed,
                                   ks = ks,
                                   ls = ls)
    results_list[[paste0("seed_", seed)]] <- results
  }

  return(results_list)
  }

performance_knn_stepauc <- evaluate_with_seeds_KNN_stepAUC(
  evaluation_function = double_cross_validation_knn_stepAUC,
  data = data_numeric,
  outer = 5,
  inner = 2,
  target_var = "PCR",
  ks = c(10, 16, 20, 25, 30, 32),
  ls = c(0, 1, 2),
  vars = todas_variables_numeric,
  threshold = .35,
  seeds = seeds
)

load("RData_Files_Algorithms/KNN_StepAUC.RData")

3.5.5 Embedded (Lasso)

double_cross_validation_knn_lasso <- function(data, target_name, outer, inner, ks, ls, variables, threshold, seed) {
  set.seed(seed)
  performance_metrics <- data.frame(Fold = integer(), TP = integer(), TN = integer(),
                                    FP = integer(), FN = integer(), AUC = numeric(),
                                    Best_k = integer(), Best_l = integer(), 
                                    Best_Variables = character(), stringsAsFactors = FALSE)
  
  inner_fold_metrics <- data.frame(Fold = integer(), Inner_Fold = integer(), 
                                   TP = integer(), TN = integer(), FP = integer(), 
                                   FN = integer(), AUC = numeric(), stringsAsFactors = FALSE)
  
  # Convertir las variables del dataset con model.matrix
  target <- as.numeric(data[[target_name]]) - 1
  features <- data[, setdiff(names(data), target_name)]
  data_matrix <- as.data.frame(model.matrix(~ ., data=features)[,-1])  # Eliminar el intercepto
  data_matrix[[target_name]] <- target
  
  outer_folds <- createFolds(data_matrix[[target_name]], k = outer)
  
  for (outer_index in seq_along(outer_folds)) {
    outer_train_data <- data_matrix[-outer_folds[[outer_index]], ]
    outer_test_data <- data_matrix[outer_folds[[outer_index]], ]
    
    cat(sprintf("Outer Fold %d\n", outer_index))
    cat("Training Data - Class 0: ", sum(outer_train_data[[target_name]] == 0), "\n")
    cat("Training Data - Class 1: ", sum(outer_train_data[[target_name]] == 1), "\n")
    cat("Testing Data - Class 0: ", sum(outer_test_data[[target_name]] == 0), "\n")
    cat("Testing Data - Class 1: ", sum(outer_test_data[[target_name]] == 1), "\n")
    
    inner_folds <- createFolds(outer_train_data[[target_name]], k = inner)
    best_model <- NULL
    best_auc <- 0
    best_k <- NULL
    best_l <- NULL
    best_vars <- colnames(outer_train_data)[!colnames(outer_train_data) %in% target_name]
    
    grid <- expand.grid(k = ks, l = ls)
    
    select_vars_lasso <- function(data, target_name) {
      target <- as.numeric(data[[target_name]])
      features <- data[, setdiff(names(data), target_name)]
      X <- model.matrix(~ ., data=features)[,-1]  # Eliminamos el intercepto
      y <- target
      
      # Ajustamos el modelo Lasso
      lasso_model <- cv.glmnet(X, y, alpha=1)
      lasso_coef <- coef(lasso_model, s = "lambda.min")
      
      # Seleccionamos las variables no nulas
      selected_vars <- rownames(lasso_coef)[lasso_coef[, 1] != 0]
      selected_vars <- setdiff(selected_vars, "(Intercept)")
      #cat("Selected variables by Lasso: ", paste(selected_vars, collapse=", "), "\n")
      return(selected_vars)
    }
    
    stepAUC <- function(vars, data, target_name, grid, inner_folds) {
      best_inner_auc <- 0
      best_inner_model <- NULL
      best_inner_k <- NULL
      best_inner_l <- NULL
      best_inner_vars <- vars
      improved <- TRUE
      
      while (improved && length(vars) > 1) {
        improved <- FALSE
        best_auc_in_step <- best_inner_auc
        
        for (params in 1:nrow(grid)) {
          auc_vals <- numeric()
          
          for (inner_index in seq_along(inner_folds)) {
            inner_train_data <- outer_train_data[-inner_folds[[inner_index]], ]
            inner_test_data <- outer_train_data[inner_folds[[inner_index]], ]
            
            selected_vars <- select_vars_lasso(inner_train_data, target_name)
            formula <- as.formula(paste(target_name, "~", paste(selected_vars, collapse = "+")))
            #cat("Testing with formula: ", deparse(formula), "\n")
            #cat("Grid parameters - k: ", grid$k[params], ", l: ", grid$l[params], "\n")
            
            result <- tryCatch({
              train_x <- inner_train_data[, selected_vars, drop = FALSE]
              test_x <- inner_test_data[, selected_vars, drop = FALSE]
              train_y <- inner_train_data[[target_name]]
              test_y <- inner_test_data[[target_name]]
              
              predictions <- knn(train_x, test_x, train_y, k = grid$k[params], l = grid$l[params], prob = TRUE)
              prob_predictions <- attr(predictions, "prob")
              prob_predictions <- ifelse(predictions == "1", prob_predictions, 1 - prob_predictions)
              roc_curve <- roc(inner_test_data[[target_name]], prob_predictions)
              auc_vals[inner_index] <- roc_curve$auc
              #cat("AUC for inner fold ", inner_index, ": ", roc_curve$auc, "\n")
              
            }, error = function(e) {
              cat("Error with variables: ", paste(selected_vars, collapse = ", "), "\n")
              cat("Error message: ", e$message, "\n")
              auc_vals[inner_index] <- NA  # Set to NA or some other indicator of failure
            })
          }
          
          mean_auc <- mean(auc_vals, na.rm = TRUE)
          #cat("Mean AUC for params - k: ", grid$k[params], ", l: ", grid$l[params], ": ", mean_auc, "\n")
          if (!is.na(mean_auc) && mean_auc > best_auc_in_step) {
            best_auc_in_step <- mean_auc
            best_inner_model <- list(k = grid$k[params], l = grid$l[params])
            best_inner_k <- grid$k[params]
            best_inner_l <- grid$l[params]
            best_inner_vars <- selected_vars
            improved <- TRUE
          }
        }
        
        if (improved) {
          vars <- best_inner_vars
          best_inner_auc <- best_auc_in_step
          cat("Vars improved: ", paste(vars, collapse=", "), "\n")
          cat("AUC: ", best_inner_auc, "\n")
        }
      }
      
      list(model = best_inner_model, auc = best_inner_auc, 
           k = best_inner_k, l = best_inner_l, 
           vars = best_inner_vars)
    }
    
    best_result <- stepAUC(best_vars, outer_train_data, target_name, grid, inner_folds)
    best_model <- best_result$model
    best_auc <- best_result$auc
    best_k <- best_result$k
    best_l <- best_result$l
    best_vars <- best_result$vars
    
    if (is.null(best_model)) {
      cat("Error: No valid model found for outer fold ", outer_index, "\n")
      next
    }

    # Evaluar el mejor modelo en el validation split del inner loop
    inner_auc_vals <- numeric()
    for (inner_index in seq_along(inner_folds)) {
      inner_train_data <- outer_train_data[-inner_folds[[inner_index]], ]
      inner_test_data <- outer_train_data[inner_folds[[inner_index]], ]
      
      train_x <- inner_train_data[, best_vars, drop = FALSE]
      test_x <- inner_test_data[, best_vars, drop = FALSE]
      train_y <- inner_train_data[[target_name]]
      test_y <- inner_test_data[[target_name]]
      
      predictions <- knn(train_x, test_x, train_y, k = best_k, l = best_l, prob = TRUE)
      prob_predictions <- attr(predictions, "prob")
      prob_predictions <- ifelse(predictions == "1", prob_predictions, 1 - prob_predictions)
      roc_curve <- roc(inner_test_data[[target_name]], prob_predictions)
      inner_auc_vals[inner_index] <- roc_curve$auc
      cat("Inner fold ", inner_index, " AUC: ", roc_curve$auc, "\n")

      pred_class <- ifelse(prob_predictions > threshold, 1, 0)
      confusion <- table(Actual = inner_test_data[[target_name]], Predicted = pred_class)
      TP <- ifelse(nrow(confusion) >= 2 && ncol(confusion) >= 2, confusion[2, 2], 0)
      TN <- ifelse(nrow(confusion) >= 1 && ncol(confusion) >= 1, confusion[1, 1], 0)
      FP <- ifelse(nrow(confusion) >= 1 && ncol(confusion) >= 2, confusion[1, 2], 0)
      FN <- ifelse(nrow(confusion) >= 2 && ncol(confusion) >= 1, confusion[2, 1], 0)
      auc <- roc_curve$auc
      
      inner_metrics <- setNames(c(outer_index, inner_index, TP, TN, FP, FN, auc), names(inner_fold_metrics))
      inner_fold_metrics <- rbind(inner_fold_metrics, as.data.frame(t(inner_metrics)))
    }
    
    best_inner_auc <- mean(inner_auc_vals, na.rm = TRUE)
    cat("Best inner AUC: ", best_inner_auc, "\n")
    
    # Reentrenar el mejor modelo con todos los datos de outer_train_data
    train_x <- outer_train_data[, best_vars, drop = FALSE]
    train_y <- outer_train_data[[target_name]]
    best_model <- list(train_x = train_x, train_y = train_y, k = best_k, l = best_l)
    
    # Evaluar el mejor modelo en el conjunto de prueba del outer fold
    test_x <- outer_test_data[, best_vars, drop = FALSE]
    test_y <- outer_test_data[[target_name]]
    predictions <- knn(train = best_model$train_x, test = test_x, cl = best_model$train_y, k = best_k, l = best_l, prob = TRUE)
    prob_predictions <- attr(predictions, "prob")
    prob_predictions <- ifelse(predictions == "1", prob_predictions, 1 - prob_predictions)
    pred_class <- ifelse(prob_predictions > threshold, 1, 0)
    confusion <- table(Actual = test_y, Predicted = pred_class)
    
    if (nrow(confusion) < 2 || ncol(confusion) < 2) {
      TP <- ifelse(nrow(confusion) >= 2 && ncol(confusion) >= 2, confusion[2, 2], 0)
      TN <- ifelse(nrow(confusion) >= 1 && ncol(confusion) >= 1, confusion[1, 1], 0)
      FP <- ifelse(nrow(confusion) >= 1 && ncol(confusion) >= 2, confusion[1, 2], 0)
      FN <- ifelse(nrow(confusion) >= 2 && ncol(confusion) >= 1, confusion[2, 1], 0)
    } else {
      TP <- confusion[2, 2]
      TN <- confusion[1, 1]
      FP <- confusion[1, 2]
      FN <- confusion[2, 1]
    }
    auc <- roc(outer_test_data[[target_name]], prob_predictions)$auc
    
    metrics <- setNames(c(outer_index, TP, TN, FP, FN, auc, 
                          best_k, best_l, paste(best_vars, collapse = ",")), 
                        names(performance_metrics))
    performance_metrics <- rbind(performance_metrics, as.data.frame(t(metrics)))
    
    cat("Performance:\n", paste(performance_metrics), "\n")
    cat(sprintf("Confusion Matrix for Fold %d:\n", outer_index))
    print(confusion)
    cat(sprintf("Metrics for Fold %d:\n", outer_index))
    cat(sprintf("TP: %d, TN: %d, FP: %d, FN: %d, AUC: %f, Best k: %d, Best l: %d, Best Variables: %s\n",
                TP, TN, FP, FN, auc, best_k, best_l, paste(best_vars, collapse = ",")))
  }
  
  list(performance_metrics = performance_metrics, inner_fold_metrics = inner_fold_metrics)
}
evaluate_with_seeds_KNN_lasso <- 
  function(evaluation_function, data, target_var, vars, threshold, seeds, outer, inner, ks, ls) {
  results_list <- list()

  for (seed in seeds) {
    set.seed(seed)
    cat("Evaluating with seed:", seed, "\n")

    results <- evaluation_function(data = data, 
                                   target_name = target_var, 
                                   variables = vars, 
                                   threshold = threshold, 
                                   outer = outer, 
                                   inner = inner, 
                                   seed = seed,
                                   ks = ks,
                                   ls = ls)
    results_list[[paste0("seed_", seed)]] <- results
  }

  return(results_list)
  }

performance_knn_lasso <- evaluate_with_seeds_KNN_lasso(
  evaluation_function = double_cross_validation_knn_lasso,
  data = data_numeric,
  outer = 5,
  inner = 2,
  target_var = "PCR",
  ks = c(10, 16, 20, 25, 30, 32),
  ls = c(0, 1, 2),
  vars = todas_variables_numeric,
  threshold = .35,
  seeds = seeds
)

load("RData_Files_Algorithms/KNN_Lasso.RData")

3.6 Support Vector Machines

3.6.1 Aparente

evaluate_aparent_performance_model_svm <- function(data, target_var, model_func, vars = NULL, threshold = 0.5) {
  if (is.null(vars)) {
    vars <- setdiff(names(data), target_var)
  }
  
  data[[target_var]] <- factor(data[[target_var]], levels = c(0, 1))
  formula <- as.formula(paste(target_var, "~", paste(vars, collapse = "+")))
  
  model <- model_func(formula, data)
  predictions <- predict(model, newdata = data, probability = TRUE)
  predictions <- attr(predictions, "probabilities")[, 2]
  
  actual_classes <- data[[target_var]]
  predicted_classes <- ifelse(predictions > threshold, 1, 0)
  confusion_matrix <- table(predicted_classes, actual_classes)
  
  tp <- ifelse("1" %in% rownames(confusion_matrix) && "1" %in% colnames(confusion_matrix), confusion_matrix["1", "1"], 0)
  tn <- ifelse("0" %in% rownames(confusion_matrix) && "0" %in% colnames(confusion_matrix), confusion_matrix["0", "0"], 0)
  fp <- ifelse("1" %in% rownames(confusion_matrix) && "0" %in% colnames(confusion_matrix), confusion_matrix["1", "0"], 0)
  fn <- ifelse("0" %in% rownames(confusion_matrix) && "1" %in% colnames(confusion_matrix), confusion_matrix["0", "1"], 0)
  
  accuracy <- (tp + tn) / (tp + tn + fp + fn)
  precision <- ifelse(tp + fp > 0, tp / (tp + fp), 0)
  recall <- ifelse(tp + fn > 0, tp / (tp + fn), 0)
  f1_score <- ifelse(precision + recall > 0, 2 * (precision * recall) / (precision + recall), 0)
  roc_obj <- pROC::roc(actual_classes, predictions)
  
  list(confusion_matrix = confusion_matrix, accuracy = accuracy, precision = precision, 
       recall = recall, f1_score = f1_score, roc_curve = roc_obj)
}


svm_model <- function(formula, data) {
  set.seed(90)
  x <- model.matrix(formula, data)
  y <- data[[all.vars(formula)[1]]]
  y <- factor(y, levels = c(0, 1))
  svm(x = x, y = y, kernel = "linear", cost = 100, probability = TRUE, maxiter=300, random_state=90)
}
# svm(x = x, y = y, kernel = "poly", cost = 5, gamma = .5, probability = TRUE, maxiter=150)
# svm(x = x, y = y, kernel = "linear", cost = 50, probability = TRUE, maxiter=300, random_state=90)


resultadosAparentesSVM <- evaluate_aparent_performance_model_svm(data = data_numeric, target_var = "PCR",
                                                             model_func = svm_model,
                                                             vars = ".",
                                                             threshold = .35)

load("RData_Files_Algorithms/SVM_Aparente.RData")

3.6.2 Todas las Variables

double_cross_validation_svm <- function(data, target_name, outer, inner, kernels, costs, variables, threshold, seed) {
  set.seed(seed)
  performance_metrics <- data.frame(Fold = integer(), TP = integer(), TN = integer(),
                                    FP = integer(), FN = integer(), AUC = numeric(),
                                    Best_Kernel = character(), Best_Cost = numeric(), stringsAsFactors = FALSE)

  data[[target_name]] <- as.factor(data[[target_name]])
  outer_folds <- createFolds(data[[target_name]], k = outer, list = TRUE, returnTrain = FALSE)
  
  for (outer_index in seq_along(outer_folds)) {
    outer_train_data <- data[-outer_folds[[outer_index]], ]
    outer_test_data <- data[outer_folds[[outer_index]], ]

    cat(sprintf("Outer Fold %d\n", outer_index))
    cat("Training Data - Class 0: ", sum(outer_train_data[[target_name]] == "0"), "\n")
    cat("Training Data - Class 1: ", sum(outer_train_data[[target_name]] == "1"), "\n")
    cat("Testing Data - Class 0: ", sum(outer_test_data[[target_name]] == "0"), "\n")
    cat("Testing Data - Class 1: ", sum(outer_test_data[[target_name]] == "1"), "\n")
    
    if (sum(outer_train_data[[target_name]] == "0") == 0 || sum(outer_train_data[[target_name]] == "1") == 0) {
      cat("Skipping outer fold due to lack of class diversity\n")
      next
    }

    inner_folds <- createFolds(outer_train_data[[target_name]], k = inner, list = TRUE, returnTrain = FALSE)
    best_model <- NULL
    best_auc <- 0
    best_kernel <- NULL
    best_cost <- NULL

    grid <- expand.grid(kernel = kernels, cost = costs)

    for (params in 1:nrow(grid)) {
      inner_auc <- numeric()

      for (inner_index in seq_along(inner_folds)) {
        inner_train_data <- outer_train_data[-inner_folds[[inner_index]], ]
        inner_test_data <- outer_train_data[inner_folds[[inner_index]], ]

        if (sum(inner_train_data[[target_name]] == "0") == 0 || sum(inner_train_data[[target_name]] == "1") == 0) {
          cat("Skipping inner fold due to lack of class diversity\n")
          next
        }

        result <- tryCatch({
          model <- svm(reformulate(variables, target_name), 
                       data = inner_train_data, kernel = grid$kernel[params], cost = grid$cost[params], 
                       probability = TRUE)

          pred <- predict(model, inner_test_data, probability = TRUE)
          prob <- attr(pred, "probabilities")[, 2]
          roc_curve <- pROC::roc(inner_test_data[[target_name]], prob)
          inner_auc[inner_index] <- roc_curve$auc
        }, error = function(e) {
          cat("Error with parameters: kernel =", grid$kernel[params], ", cost =", grid$cost[params], "\n")
          cat("Error message: ", e$message, "\n")
          inner_auc[inner_index] <- 0
        })
      }

      if (mean(inner_auc, na.rm = TRUE) > best_auc) {
        best_auc <- mean(inner_auc, na.rm = TRUE)
        best_model <- model
        best_kernel <- grid$kernel[params]
        best_cost <- grid$cost[params]
      }
    }
    best_model <- svm(reformulate(variables, target_name), 
                      data = outer_train_data, kernel = best_kernel, cost = best_cost, 
                      probability = TRUE)
    predictions <- predict(best_model, outer_test_data, probability = TRUE)
    prob <- attr(predictions, "probabilities")[, 2]
    pred_class <- ifelse(prob > threshold, 1, 0)  # Asignar clases basadas en umbral
    confusion <- table(outer_test_data[[target_name]], pred_class)
    
    if (nrow(confusion) < 2 || ncol(confusion) < 2) {
      TP <- ifelse(nrow(confusion) >= 2 && ncol(confusion) >= 2, confusion[2, 2], 0)
      TN <- ifelse(nrow(confusion) >= 1 && ncol(confusion) >= 1, confusion[1, 1], 0)
      FP <- ifelse(nrow(confusion) >= 1 && ncol(confusion) >= 2, confusion[1, 2], 0)
      FN <- ifelse(nrow(confusion) >= 2 && ncol(confusion) >= 1, confusion[2, 1], 0)
    } else {
      TP <- confusion[2, 2]
      TN <- confusion[1, 1]
      FP <- confusion[1, 2]
      FN <- confusion[2, 1]
    }
    auc <- pROC::roc(outer_test_data[[target_name]], prob)$auc

    # Crear un vector con los nombres correctos y convertir a dataframe
    metrics <- setNames(c(outer_index, TP, TN, FP, FN, auc, best_kernel, best_cost), names(performance_metrics))
    performance_metrics <- rbind(performance_metrics, as.data.frame(t(metrics)))
  }
  
  return(performance_metrics)
}

evaluate_with_seeds_SVM_TodasVariables_Asociacion <- 
  function(evaluation_function, data, target_var, vars, threshold, seeds, outer, inner, kernels, costs) {
  results_list <- list()

  for (seed in seeds) {
    set.seed(seed)
    cat("Evaluating with seed:", seed, "\n")

    results <- evaluation_function(data = data, 
                                   target_name = target_var, 
                                   variables = vars, 
                                   threshold = threshold, 
                                   outer = outer, 
                                   inner = inner, 
                                   seed = seed,
                                   kernels = kernels,
                                   costs = costs)
    results_list[[paste0("seed_", seed)]] <- results
  }

  return(results_list)
  }

performance_SVM_todasvariables <- evaluate_with_seeds_SVM_TodasVariables_Asociacion(
  evaluation_function = double_cross_validation_svm,
  data = data_factor,
  outer = 5,
  inner = 2,
  target_var = "PCR",
  kernels = c("linear", "poly"),
  costs = c(50, 100, 200),
  vars = todas_variables_factor,
  threshold = .35,
  seeds = seeds
)

load("RData_Files_Algorithms/SVM_TodasVariables.RData")

3.6.3 Filtrado (Asociación)

performance_svm_asociacion <- evaluate_with_seeds_SVM_TodasVariables_Asociacion(
  evaluation_function = double_cross_validation_svm,
  data = data_factor,
  outer = 5,
  inner = 2,
  target_var = "PCR",
  kernels = c("linear", "poly"),
  costs = c(50, 100, 200),
  vars = asociacion,
  threshold = .35,
  seeds = seeds
)

load("RData_Files_Algorithms/SVM_Asociacion.RData")

3.6.4 Wrapped (StepAUC)

double_cross_validation_svm_stepAUC <- function(data, target_name, outer, inner, kernels, costs, variables, threshold, seed) {
  set.seed(seed)
  performance_metrics <- data.frame(Fold = integer(), TP = integer(), TN = integer(),
                                    FP = integer(), FN = integer(), AUC = numeric(),
                                    Best_Kernel = character(), Best_Cost = numeric(), Best_Variables = character(), stringsAsFactors = FALSE)
  
  inner_fold_metrics <- data.frame(Fold = integer(), Inner_Fold = integer(), 
                                   TP = integer(), TN = integer(), FP = integer(), 
                                   FN = integer(), AUC = numeric(), stringsAsFactors = FALSE)
  
  data[[target_name]] <- as.factor(data[[target_name]])
  outer_folds <- createFolds(data[[target_name]], k = outer, list = TRUE, returnTrain = FALSE)
  
  for (outer_index in seq_along(outer_folds)) {
    outer_train_data <- data[-outer_folds[[outer_index]], ]
    outer_test_data <- data[outer_folds[[outer_index]], ]
    
    cat(sprintf("Outer Fold %d\n", outer_index))
    cat("Training Data - Class 0: ", sum(outer_train_data[[target_name]] == 0), "\n")
    cat("Training Data - Class 1: ", sum(outer_train_data[[target_name]] == 1), "\n")
    cat("Training Data - Total: ", nrow(outer_train_data), "\n")
    cat("Testing Data - Class 0: ", sum(outer_test_data[[target_name]] == 0), "\n")
    cat("Testing Data - Class 1: ", sum(outer_test_data[[target_name]] == 1), "\n")
    cat("Testing Data - Total: ", nrow(outer_test_data), "\n")
    
    inner_folds <- createFolds(outer_train_data[[target_name]], k = inner, list = TRUE, returnTrain = FALSE)
    best_model <- NULL
    best_auc <- 0
    best_kernel <- NULL
    best_cost <- NULL
    best_vars <- variables
    
    grid <- expand.grid(kernel = kernels, cost = costs)
    
    stepAUC <- function(vars, outer_train_data, target_name, grid, inner_folds) {
      current_vars <- vars
      best_inner_auc <- 0
      best_inner_model <- NULL
      best_inner_kernel <- NULL
      best_inner_cost <- NULL
      best_inner_vars <- current_vars
      improved <- TRUE
      
      while (improved && length(current_vars) > 1) {
        improved <- FALSE
        best_auc_in_step <- best_inner_auc
        
        for (var in current_vars) {
          temp_vars <- setdiff(current_vars, var)
          inner_auc <- numeric()
          
          for (params in 1:nrow(grid)) {
            auc_vals <- numeric()
            
            for (inner_index in seq_along(inner_folds)) {
              inner_train_data <- outer_train_data[-inner_folds[[inner_index]], ]
              inner_test_data <- outer_train_data[inner_folds[[inner_index]], ]
              
              result <- tryCatch({
                model <- svm(as.formula(paste(target_name, "~", paste(temp_vars, collapse = "+"))),
                             data = inner_train_data, kernel = grid$kernel[params], cost = grid$cost[params], 
                             probability = TRUE)
                
                predictions <- predict(model, newdata = inner_test_data, probability = TRUE)
                predictions <- attr(predictions, "probabilities")[, 2]
                roc_curve <- roc(inner_test_data[[target_name]], predictions)
                auc_vals[inner_index] <- roc_curve$auc
                
              }, error = function(e) {
                cat("Error with variables:", paste(temp_vars, collapse = ", "), "\n")
                cat("Error message: ", e$message, "\n")
                auc_vals[inner_index] <- 0  # Set to NA or some other indicator of failure
              })
            }
            
            if (mean(auc_vals, na.rm = TRUE) > best_auc_in_step) {
              best_auc_in_step <- mean(auc_vals, na.rm = TRUE)
              best_inner_model <- model
              best_inner_kernel <- grid$kernel[params]
              best_inner_cost <- grid$cost[params]
              best_inner_vars <- temp_vars
              improved <- TRUE
            }
          }
        }
        
        if (improved) {
          current_vars <- best_inner_vars
          best_inner_auc <- best_auc_in_step
          cat("Vars improved: ", paste(current_vars, collapse=", "), "\n")
          cat("AUC: ", paste(best_inner_auc), "\n")
        }
      }
      
      list(model = best_inner_model, auc = best_inner_auc, 
           kernel = best_inner_kernel, cost = best_inner_cost, 
           vars = best_inner_vars)
      
    }
    
    best_result <- stepAUC(variables, outer_train_data, target_name, grid, inner_folds)
    best_model <- best_result$model
    best_auc <- best_result$auc
    best_kernel <- best_result$kernel
    best_cost <- best_result$cost
    best_vars <- best_result$vars

    # Evaluar el mejor modelo en el validation split del inner loop 
    inner_auc_vals <- numeric()
    for (inner_index in seq_along(inner_folds)) {
      inner_train_data <- outer_train_data[-inner_folds[[inner_index]], ]
      inner_test_data <- outer_train_data[inner_folds[[inner_index]], ]
      
      predictions <- predict(best_model, newdata = inner_test_data, probability = TRUE)
      predictions <- attr(predictions, "probabilities")[, 2]
      roc_curve <- roc(inner_test_data[[target_name]], predictions)
      inner_auc_vals[inner_index] <- roc_curve$auc

      pred_class <- ifelse(predictions > threshold, 1, 0)
      confusion <- table(Actual = inner_test_data[[target_name]], Predicted = pred_class)
      TP <- ifelse(nrow(confusion) >= 2 && ncol(confusion) >= 2, confusion[2, 2], 0)
      TN <- ifelse(nrow(confusion) >= 1 && ncol(confusion) >= 1, confusion[1, 1], 0)
      FP <- ifelse(nrow(confusion) >= 1 && ncol(confusion) >= 2, confusion[1, 2], 0)
      FN <- ifelse(nrow(confusion) >= 2 && ncol(confusion) >= 1, confusion[2, 1], 0)
      auc <- roc_curve$auc
      
      inner_metrics <- setNames(c(outer_index, inner_index, TP, TN, FP, FN, auc), names(inner_fold_metrics))
      inner_fold_metrics <- rbind(inner_fold_metrics, as.data.frame(t(inner_metrics)))
    }
    
    best_inner_auc <- mean(inner_auc_vals, na.rm = TRUE)
    # Reentrenar el mejor modelo con todos los datos de outer_train_data
    best_model <- svm(as.formula(paste(target_name, "~", paste(best_vars, collapse = "+"))),
                      data = outer_train_data, kernel = best_kernel, cost = best_cost, 
                      probability = TRUE)
    predictions <- predict(best_model, newdata = outer_test_data, probability = TRUE)
    predictions <- attr(predictions, "probabilities")[, 2]
    pred_class <- ifelse(predictions > threshold, 1, 0)  # Asignar clases basadas en umbral
    confusion <- table(Actual = outer_test_data[[target_name]], Predicted = pred_class)
    
    if (nrow(confusion) < 2 || ncol(confusion) < 2) {
      TP <- ifelse(nrow(confusion) >= 2 && ncol(confusion) >= 2, confusion[2, 2], 0)
      TN <- ifelse(nrow(confusion) >= 1 && ncol(confusion) >= 1, confusion[1, 1], 0)
      FP <- ifelse(nrow(confusion) >= 1 && ncol(confusion) >= 2, confusion[1, 2], 0)
      FN <- ifelse(nrow(confusion) >= 2 && ncol(confusion) >= 1, confusion[2, 1], 0)
    } else {
      TP <- confusion[2, 2]
      TN <- confusion[1, 1]
      FP <- confusion[1, 2]
      FN <- confusion[2, 1]
    }
    auc <- roc(outer_test_data[[target_name]], predictions)$auc
    
    metrics <- setNames(c(outer_index, TP, TN, FP, FN, auc, 
                          best_kernel, best_cost, paste(best_vars, collapse = ",")), 
                        names(performance_metrics))
    performance_metrics <- rbind(performance_metrics, as.data.frame(t(metrics)))
    
    cat("Performance:", paste(performance_metrics))
    cat(sprintf("Confusion Matrix for Fold %d:\n", outer_index))
    print(confusion)
    cat(sprintf("Metrics for Fold %d:\n", outer_index))
    cat(sprintf("TP: %d, TN: %d, FP: %d, FN: %d, AUC: %f, Best Kernel: %s, Best Cost: %f, Best Variables: %s\n",
                TP, TN, FP, FN, auc, best_kernel, best_cost, paste(best_vars, collapse = ",")))
  }
  
  list(performance_metrics = performance_metrics, inner_fold_metrics = inner_fold_metrics)
}

evaluate_with_seeds_SVM_stepAUC <- 
  function(evaluation_function, data, target_var, vars, threshold, seeds, outer, inner, kernels, costs) {
  results_list <- list()

  for (seed in seeds) {
    set.seed(seed)
    cat("Evaluating with seed:", seed, "\n")

    results <- evaluation_function(data = data, 
                                   target_name = target_var, 
                                   variables = vars, 
                                   threshold = threshold, 
                                   outer = outer, 
                                   inner = inner, 
                                   seed = seed,
                                   kernels = kernels,
                                   costs = costs)
    results_list[[paste0("seed_", seed)]] <- results
  }

  return(results_list)
  }

performance_SVM_stepauc <- evaluate_with_seeds_SVM_stepAUC(
  evaluation_function = double_cross_validation_svm_stepAUC,
  data = data_numeric,
  outer = 5,
  inner = 2,
  target_var = "PCR",
  kernels = c("linear", "poly"),
  costs = c(50, 100, 200),
  vars = todas_variables_numeric,
  threshold = .35,
  seeds = seeds
)

load("RData_Files_Algorithms/SVM_StepAUC.RData")

3.6.5 Embedded (Lasso)

double_cross_validation_svm_lasso <- function(data, target_name, outer, inner, kernels, costs, variables, threshold, seed) {
  set.seed(seed)
  performance_metrics <- data.frame(Fold = integer(), TP = integer(), TN = integer(),
                                    FP = integer(), FN = integer(), AUC = numeric(),
                                    Best_Kernel = character(), Best_Cost = numeric(),
                                    Best_Variables = character(), stringsAsFactors = FALSE)
  
  inner_fold_metrics <- data.frame(Fold = integer(), Inner_Fold = integer(), 
                                   TP = integer(), TN = integer(), FP = integer(), 
                                   FN = integer(), AUC = numeric(), stringsAsFactors = FALSE)
  
  # Convertir las variables del dataset con model.matrix
  target <- as.factor(data[[target_name]])  # Aseguramos que la variable de respuesta sea factor
  features <- data[, setdiff(names(data), target_name)]
  data_matrix <- as.data.frame(model.matrix(~ ., data=features)[,-1])  # Eliminar el intercepto
  data_matrix[[target_name]] <- target
  
  outer_folds <- createFolds(data_matrix[[target_name]], k = outer)
  
  for (outer_index in seq_along(outer_folds)) {
    outer_train_data <- data_matrix[-outer_folds[[outer_index]], ]
    outer_test_data <- data_matrix[outer_folds[[outer_index]], ]
    
    cat(sprintf("Outer Fold %d\n", outer_index))
    cat("Training Data - Class 0: ", sum(outer_train_data[[target_name]] == 0), "\n")
    cat("Training Data - Class 1: ", sum(outer_train_data[[target_name]] == 1), "\n")
    cat("Testing Data - Class 0: ", sum(outer_test_data[[target_name]] == 0), "\n")
    cat("Testing Data - Class 1: ", sum(outer_test_data[[target_name]] == 1), "\n")
    
    inner_folds <- createFolds(outer_train_data[[target_name]], k = inner)
    best_model <- NULL
    best_auc <- 0
    best_kernel <- NULL
    best_cost <- NULL
    best_vars <- colnames(outer_train_data)[!colnames(outer_train_data) %in% target_name]
    
    grid <- expand.grid(kernel = kernels, cost = costs)
    
    select_vars_lasso <- function(data, target_name) {
      target <- as.numeric(data[[target_name]]) - 1  # Convertimos a numérico para Lasso
      features <- data[, setdiff(names(data), target_name)]
      X <- model.matrix(~ ., data=features)[,-1]  # Eliminamos el intercepto
      y <- target
      
      # Ajustamos el modelo Lasso
      lasso_model <- cv.glmnet(X, y, alpha=1)
      lasso_coef <- coef(lasso_model, s = "lambda.min")
      
      # Seleccionamos las variables no nulas
      selected_vars <- rownames(lasso_coef)[lasso_coef[, 1] != 0]
      selected_vars <- setdiff(selected_vars, "(Intercept)")
      cat("Selected variables by Lasso: ", paste(selected_vars, collapse=", "), "\n")
      return(selected_vars)
    }
    
    stepAUC <- function(vars, data, target_name, grid, inner_folds) {
      best_inner_auc <- 0
      best_inner_model <- NULL
      best_inner_kernel <- NULL
      best_inner_cost <- NULL
      best_inner_vars <- vars
      improved <- TRUE
      
      while (improved && length(vars) > 1) {
        improved <- FALSE
        best_auc_in_step <- best_inner_auc
        
        for (params in 1:nrow(grid)) {
          auc_vals <- numeric()
          
          for (inner_index in seq_along(inner_folds)) {
            inner_train_data <- outer_train_data[-inner_folds[[inner_index]], ]
            inner_test_data <- outer_train_data[inner_folds[[inner_index]], ]
            
            selected_vars <- select_vars_lasso(inner_train_data, target_name)
            formula <- as.formula(paste(target_name, "~", paste(selected_vars, collapse = "+")))
            cat("Testing with formula: ", deparse(formula), "\n")
            cat("Grid parameters - kernel: ", grid$kernel[params], ", cost: ", grid$cost[params], "\n")
            
            result <- tryCatch({
              model <- svm(formula, data = inner_train_data, kernel = grid$kernel[params], cost = grid$cost[params], probability = TRUE)
              
              predictions <- predict(model, inner_test_data, probability = TRUE)
              predictions_prob <- attr(predictions, "probabilities")[, 2]
              roc_curve <- roc(inner_test_data[[target_name]], predictions_prob)
              auc_vals[inner_index] <- roc_curve$auc
              cat("AUC for inner fold ", inner_index, ": ", roc_curve$auc, "\n")
              
            }, error = function(e) {
              cat("Error with variables: ", paste(selected_vars, collapse = ", "), "\n")
              cat("Error message: ", e$message, "\n")
              auc_vals[inner_index] <- NA  # Set to NA or some other indicator of failure
            })
          }
          
          mean_auc <- mean(auc_vals, na.rm = TRUE)
          cat("Mean AUC for params - kernel: ", grid$kernel[params], ", cost: ", grid$cost[params], ": ", mean_auc, "\n")
          if (!is.na(mean_auc) && mean_auc > best_auc_in_step) {
            best_auc_in_step <- mean_auc
            best_inner_model <- model
            best_inner_kernel <- grid$kernel[params]
            best_inner_cost <- grid$cost[params]
            best_inner_vars <- selected_vars
            improved <- TRUE
          }
        }
        
        if (improved) {
          vars <- best_inner_vars
          best_inner_auc <- best_auc_in_step
          cat("Vars improved: ", paste(vars, collapse=", "), "\n")
          cat("AUC: ", best_inner_auc, "\n")
        }
      }
      
      list(model = best_inner_model, auc = best_inner_auc, 
           kernel = best_inner_kernel, cost = best_inner_cost, 
           vars = best_inner_vars)
    }
    
    best_result <- stepAUC(best_vars, outer_train_data, target_name, grid, inner_folds)
    best_model <- best_result$model
    best_auc <- best_result$auc
    best_kernel <- best_result$kernel
    best_cost <- best_result$cost
    best_vars <- best_result$vars
    
    if (is.null(best_model)) {
      cat("Error: No valid model found for outer fold ", outer_index, "\n")
      next
    }

    # Evaluar el mejor modelo en el validation split del inner loop
    inner_auc_vals <- numeric()
    for (inner_index in seq_along(inner_folds)) {
      inner_train_data <- outer_train_data[-inner_folds[[inner_index]], ]
      inner_test_data <- outer_train_data[inner_folds[[inner_index]], ]
      
      predictions <- predict(best_model, inner_test_data, probability = TRUE)
      predictions_prob <- attr(predictions, "probabilities")[, 2]
      roc_curve <- roc(inner_test_data[[target_name]], predictions_prob)
      inner_auc_vals[inner_index] <- roc_curve$auc
      cat("Inner fold ", inner_index, " AUC: ", roc_curve$auc, "\n")

      pred_class <- ifelse(predictions_prob > threshold, "1", "0")
      confusion <- table(Actual = inner_test_data[[target_name]], Predicted = pred_class)
      TP <- ifelse(nrow(confusion) >= 2 && ncol(confusion) >= 2, confusion[2, 2], 0)
      TN <- ifelse(nrow(confusion) >= 1 && ncol(confusion) >= 1, confusion[1, 1], 0)
      FP <- ifelse(nrow(confusion) >= 1 && ncol(confusion) >= 2, confusion[1, 2], 0)
      FN <- ifelse(nrow(confusion) >= 2 && ncol(confusion) >= 1, confusion[2, 1], 0)
      auc <- roc_curve$auc
      
      inner_metrics <- setNames(c(outer_index, inner_index, TP, TN, FP, FN, auc), names(inner_fold_metrics))
      inner_fold_metrics <- rbind(inner_fold_metrics, as.data.frame(t(inner_metrics)))
    }
    
    best_inner_auc <- mean(inner_auc_vals, na.rm = TRUE)
    cat("Best inner AUC: ", best_inner_auc, "\n")
    best_model <- svm(as.formula(paste(target_name, "~", paste(best_vars, collapse = "+"))),
                      data = outer_train_data, kernel = best_kernel, cost = best_cost, 
                      probability = TRUE)
    predictions <- predict(best_model, outer_test_data, probability = TRUE)
    predictions_prob <- attr(predictions, "probabilities")[, 2]
    pred_class <- ifelse(predictions_prob > threshold, "1", "0")  # Asignar clases basadas en umbral
    confusion <- table(Actual = outer_test_data[[target_name]], Predicted = pred_class)
    
    if (nrow(confusion) < 2 || ncol(confusion) < 2) {
      TP <- ifelse(nrow(confusion) >= 2 && ncol(confusion) >= 2, confusion[2, 2], 0)
      TN <- ifelse(nrow(confusion) >= 1 && ncol(confusion) >= 1, confusion[1, 1], 0)
      FP <- ifelse(nrow(confusion) >= 1 && ncol(confusion) >= 2, confusion[1, 2], 0)
      FN <- ifelse(nrow(confusion) >= 2 && ncol(confusion) >= 1, confusion[2, 1], 0)
    } else {
      TP <- confusion[2, 2]
      TN <- confusion[1, 1]
      FP <- confusion[1, 2]
      FN <- confusion[2, 1]
    }
    auc <- roc(outer_test_data[[target_name]], predictions_prob)$auc
    
    metrics <- setNames(c(outer_index, TP, TN, FP, FN, auc, 
                          best_kernel, best_cost, paste(best_vars, collapse = ",")), 
                        names(performance_metrics))
    performance_metrics <- rbind(performance_metrics, as.data.frame(t(metrics)))
    
    cat("Performance:\n", paste(performance_metrics), "\n")
    cat(sprintf("Confusion Matrix for Fold %d:\n", outer_index))
    print(confusion)
    cat(sprintf("Metrics for Fold %d:\n", outer_index))
    cat(sprintf("TP: %d, TN: %d, FP: %d, FN: %d, AUC: %f, Best Kernel: %s, Best Cost: %f, Best Variables: %s\n",
                TP, TN, FP, FN, auc, best_kernel, best_cost, paste(best_vars, collapse = ",")))
  }
  
  list(performance_metrics = performance_metrics, inner_fold_metrics = inner_fold_metrics)
}

evaluate_with_seeds_SVM_lasso <- 
  function(evaluation_function, data, target_var, vars, threshold, seeds, outer, inner, kernels, costs) {
  results_list <- list()

  for (seed in seeds) {
    set.seed(seed)
    cat("Evaluating with seed:", seed, "\n")

    results <- evaluation_function(data = data, 
                                   target_name = target_var, 
                                   variables = vars, 
                                   threshold = threshold, 
                                   outer = outer, 
                                   inner = inner, 
                                   seed = seed,
                                   kernels = kernels,
                                   costs = costs)
    results_list[[paste0("seed_", seed)]] <- results
  }

  return(results_list)
  }

performance_SVM_lasso <- evaluate_with_seeds_SVM_lasso(
  evaluation_function = double_cross_validation_svm_lasso,
  data = data_factor,
  outer = 5,
  inner = 2,
  target_var = "PCR",
  kernels = c("linear", "poly"),
  costs = c(50, 100, 200),
  vars = todas_variables_factor,
  threshold = .35,
  seeds = seeds
)

load("RData_Files_Algorithms/SVM_Lasso.RData")

3.7 Artificial Neural Networks

3.7.1 Aparente

evaluate_aparent_performance_model_nnet <- function(data, target_var, model_func, vars = NULL, threshold = 0.5) {
  if (is.null(vars)) {
    vars <- setdiff(names(data), target_var)
  }
  
  data[[target_var]] <- factor(data[[target_var]], levels = c(0, 1))
  formula <- as.formula(paste(target_var, "~", paste(vars, collapse = "+")))
  
  model <- model_func(formula, data)
  predictions <- predict(model, newdata = data, type = "raw")
  
  actual_classes <- data[[target_var]]
  predicted_classes <- ifelse(predictions > threshold, 1, 0)
  confusion_matrix <- table(predicted_classes, actual_classes)
  
  tp <- ifelse("1" %in% rownames(confusion_matrix) && "1" %in% colnames(confusion_matrix), confusion_matrix["1", "1"], 0)
  tn <- ifelse("0" %in% rownames(confusion_matrix) && "0" %in% colnames(confusion_matrix), confusion_matrix["0", "0"], 0)
  fp <- ifelse("1" %in% rownames(confusion_matrix) && "0" %in% colnames(confusion_matrix), confusion_matrix["1", "0"], 0)
  fn <- ifelse("0" %in% rownames(confusion_matrix) && "1" %in% colnames(confusion_matrix), confusion_matrix["0", "1"], 0)
  
  accuracy <- (tp + tn) / (tp + tn + fp + fn)
  precision <- ifelse(tp + fp > 0, tp / (tp + fp), 0)
  recall <- ifelse(tp + fn > 0, tp / (tp + fn), 0)
  f1_score <- ifelse(precision + recall > 0, 2 * (precision * recall) / (precision + recall), 0)
  roc_obj <- pROC::roc(actual_classes, predictions)
  
  list(confusion_matrix = confusion_matrix, accuracy = accuracy, precision = precision, 
       recall = recall, f1_score = f1_score, roc_curve = roc_obj)
}


nn_model <- function(formula, data) {
  set.seed(90)
  nnet(formula = formula, data = data, size = 15, 
       decay = 0.1, maxit = 100, trace = FALSE, 
       linout = FALSE, random_state = 90)
}


resultadosAparentesNN <- evaluate_aparent_performance_model_nnet(data = data_factor, target_var = "PCR",
                                                             model_func = nn_model,
                                                             vars = ".",
                                                             threshold = .35)

load("RData_Files_Algorithms/ANN_Aparente.RData")

3.7.2 Todas las Variables

double_cross_validation_nnet <- function(data, target_name, outer, inner, sizes, decays, variables, threshold, seed) {
  set.seed(seed)
  performance_metrics <- data.frame(Fold = integer(), TP = integer(), TN = integer(),
                                    FP = integer(), FN = integer(), AUC = numeric(),
                                    Best_Size = integer(), Best_Decay = numeric(), stringsAsFactors = FALSE)
  data[[target_name]] <- as.numeric(data[[target_name]]) - 1
  outer_folds <- createFolds(data[[target_name]], k = outer, list = TRUE, returnTrain = FALSE)

  for (outer_index in seq_along(outer_folds)) {
    outer_train_data <- data[-outer_folds[[outer_index]], ]
    outer_test_data <- data[outer_folds[[outer_index]], ]

    cat(sprintf("Outer Fold %d\n", outer_index))
    cat("Training Data - Class 0: ", sum(outer_train_data[[target_name]] == 0), "\n")
    cat("Training Data - Class 1: ", sum(outer_train_data[[target_name]] == 1), "\n")
    cat("Training Data - Total: ", nrow(outer_train_data), "\n")
    cat("Testing Data - Class 0: ", sum(outer_test_data[[target_name]] == 0), "\n")
    cat("Testing Data - Class 1: ", sum(outer_test_data[[target_name]] == 1), "\n")
    cat("Testing Data - Total: ", nrow(outer_test_data), "\n")
    
    if (sum(outer_train_data[[target_name]] == 0) == 0 || sum(outer_train_data[[target_name]] == 1) == 0) {
      cat("Skipping outer fold due to lack of class diversity\n")
      next
    }

    inner_folds <- createFolds(outer_train_data[[target_name]], k = inner, list = TRUE, returnTrain = FALSE)
    best_model <- NULL
    best_auc <- 0
    best_size <- NULL
    best_decay <- NULL

    grid <- expand.grid(size = sizes, decay = decays)
    formula <- reformulate(variables, target_name)

    for (params in 1:nrow(grid)) {
      inner_auc <- numeric()

      for (inner_index in seq_along(inner_folds)) {
        inner_train_data <- outer_train_data[-inner_folds[[inner_index]], ]
        inner_test_data <- outer_train_data[inner_folds[[inner_index]], ]

        if (sum(inner_train_data[[target_name]] == 0) == 0 || sum(inner_train_data[[target_name]] == 1) == 0) {
          cat("Skipping inner fold due to lack of class diversity\n")
          next
        }

        result <- tryCatch({
          model <- nnet(formula, data = inner_train_data, size = grid$size[params], decay = grid$decay[params], linout = FALSE, maxit = 200)

          predictions <- predict(model, inner_test_data, type = "raw")
          roc_curve <- pROC::roc(inner_test_data[[target_name]], predictions)
          inner_auc[inner_index] <- roc_curve$auc
        }, error = function(e) {
          cat("Error with parameters: size =", grid$size[params], ", decay =", grid$decay[params], "\n")
          cat("Error message: ", e$message, "\n")
          inner_auc[inner_index] <- 0
        })
      }

      if (mean(inner_auc, na.rm = TRUE) > best_auc) {
        best_auc <- mean(inner_auc, na.rm = TRUE)
        best_model <- model
        best_size <- grid$size[params]
        best_decay <- grid$decay[params]
      }
    }
    model <- nnet(formula, data = outer_train_data, size = best_size, decay = best_decay, linout = FALSE, maxit = 200)

    predictions <- predict(best_model, outer_test_data, type = "raw")
    pred_class <- ifelse(predictions > threshold, 1, 0)  # Asignar clases basadas en umbral
    confusion <- table(Actual = outer_test_data[[target_name]], Predicted = pred_class)

    if (nrow(confusion) < 2 || ncol(confusion) < 2) {
      TP <- ifelse(nrow(confusion) >= 2 && ncol(confusion) >= 2, confusion[2, 2], 0)
      TN <- ifelse(nrow(confusion) >= 1 && ncol(confusion) >= 1, confusion[1, 1], 0)
      FP <- ifelse(nrow(confusion) >= 1 && ncol(confusion) >= 2, confusion[1, 2], 0)
      FN <- ifelse(nrow(confusion) >= 2 && ncol(confusion) >= 1, confusion[2, 1], 0)
    } else {
      TP <- confusion[2, 2]
      TN <- confusion[1, 1]
      FP <- confusion[1, 2]
      FN <- confusion[2, 1]
    }
    auc <- pROC::roc(outer_test_data[[target_name]], predictions)$auc

    metrics <- setNames(c(outer_index, TP, TN, FP, FN, auc, best_size, best_decay), names(performance_metrics))
    performance_metrics <- rbind(performance_metrics, as.data.frame(t(metrics)))
  }

  return(performance_metrics)
}
evaluate_with_seeds_ANN_TodasVariables_Asociacion <- 
  function(evaluation_function, data, target_var, vars, threshold, seeds, outer, inner, sizes, decays) {
  results_list <- list()

  for (seed in seeds) {
    set.seed(seed)
    cat("Evaluating with seed:", seed, "\n")

    results <- evaluation_function(data = data, 
                                   target_name = target_var, 
                                   variables = vars, 
                                   threshold = threshold, 
                                   outer = outer, 
                                   inner = inner, 
                                   seed = seed,
                                   sizes = sizes,
                                   decays = decays)
    results_list[[paste0("seed_", seed)]] <- results
  }

  return(results_list)
  }

performance_ANN_todasvariables <- evaluate_with_seeds_ANN_TodasVariables_Asociacion(
  evaluation_function = double_cross_validation_nnet,
  data = data_factor,
  outer = 5,
  inner = 2,
  target_var = "PCR",
  sizes = c(5,10,15,20),
  decays = c(0.1, 0.2, 0.3),
  vars = todas_variables_factor,
  threshold = .35,
  seeds = seeds
)

load("RData_Files_Algorithms/ANN_TodasVariables.RData")

3.7.3 Filtrado (Asociación)

performance_ANN_asociacion <- evaluate_with_seeds_ANN_TodasVariables_Asociacion(
  evaluation_function = double_cross_validation_nnet,
  data = data_factor,
  outer = 5,
  inner = 2,
  target_var = "PCR",
  sizes = c(5,10,15,20),
  decays = c(0.1, 0.2, 0.3),
  vars = asociacion,
  threshold = .35,
  seeds = seeds
)

load("RData_Files_Algorithms/ANN_Asociacion.RData")

3.7.4 Wrapped (StepAUC)

double_cross_validation_nnet_stepAUC <- function(data, target_name, outer, inner, sizes, decays, variables, threshold, seed) {
  set.seed(seed)
  performance_metrics <- data.frame(Fold = integer(), TP = integer(), TN = integer(),
                                    FP = integer(), FN = integer(), AUC = numeric(),
                                    Best_Size = integer(), Best_Decay = numeric(), Best_Variables = character(), stringsAsFactors = FALSE)
  
  inner_fold_metrics <- data.frame(Fold = integer(), Inner_Fold = integer(), 
                                   TP = integer(), TN = integer(), FP = integer(), 
                                   FN = integer(), AUC = numeric(), stringsAsFactors = FALSE)
  
  data[[target_name]] <- as.numeric(data[[target_name]]) - 1
  outer_folds <- createFolds(data[[target_name]], k = outer)
  
  for (outer_index in seq_along(outer_folds)) {
    outer_train_data <- data[-outer_folds[[outer_index]], ]
    outer_test_data <- data[outer_folds[[outer_index]], ]
    
    cat(sprintf("Outer Fold %d\n", outer_index))
    cat("Training Data - Class 0: ", sum(outer_train_data[[target_name]] == 0), "\n")
    cat("Training Data - Class 1: ", sum(outer_train_data[[target_name]] == 1), "\n")
    cat("Testing Data - Class 0: ", sum(outer_test_data[[target_name]] == 0), "\n")
    cat("Testing Data - Class 1: ", sum(outer_test_data[[target_name]] == 1), "\n")
    
    inner_folds <- createFolds(outer_train_data[[target_name]], k = inner)
    best_model <- NULL
    best_auc <- 0
    best_size <- NULL
    best_decay <- NULL
    best_vars <- variables
    
    grid <- expand.grid(size = sizes, decay = decays)
    
    stepAUC <- function(vars, data, target_name, grid, inner_folds) {
      current_vars <- vars
      best_inner_auc <- 0
      best_inner_model <- NULL
      best_inner_size <- NULL
      best_inner_decay <- NULL
      best_inner_vars <- current_vars
      improved <- TRUE
      
      while (improved && length(current_vars) > 1) {
        improved <- FALSE
        best_auc_in_step <- best_inner_auc
        
        for (var in current_vars) {
          temp_vars <- setdiff(current_vars, var)
          inner_auc <- numeric()
          
          for (params in 1:nrow(grid)) {
            auc_vals <- numeric()
            
            for (inner_index in seq_along(inner_folds)) {
              inner_train_data <- outer_train_data[-inner_folds[[inner_index]], ]
              inner_test_data <- outer_train_data[inner_folds[[inner_index]], ]
              
              result <- tryCatch({
                model <- nnet(as.formula(paste(target_name, "~", paste(temp_vars, collapse = "+"))),
                              data = inner_train_data, size = grid$size[params], decay = grid$decay[params], 
                              linout = FALSE, maxit = 200, trace = FALSE)
                
                predictions <- predict(model, inner_test_data, type = "raw")
                roc_curve <- roc(inner_test_data[[target_name]], predictions)
                auc_vals[inner_index] <- roc_curve$auc
                
                # Guardar los datos del inner fold
                pred_class <- ifelse(predictions > threshold, 1, 0)
                confusion <- table(Actual = inner_test_data[[target_name]], Predicted = pred_class)
                TP <- ifelse(nrow(confusion) >= 2 && ncol(confusion) >= 2, confusion[2, 2], 0)
                TN <- ifelse(nrow(confusion) >= 1 && ncol(confusion) >= 1, confusion[1, 1], 0)
                FP <- ifelse(nrow(confusion) >= 1 && ncol(confusion) >= 2, confusion[1, 2], 0)
                FN <- ifelse(nrow(confusion) >= 2 && ncol(confusion) >= 1, confusion[2, 1], 0)
                auc <- roc_curve$auc
                inner_metrics <- setNames(c(outer_index, inner_index, TP, TN, FP, FN, auc), names(inner_fold_metrics))
                inner_fold_metrics <<- rbind(inner_fold_metrics, as.data.frame(t(inner_metrics)))
              }, error = function(e) {
                cat("Error with variables:", paste(temp_vars, collapse = ", "), "\n")
                auc_vals[inner_index] <- 0  # Set to NA or some other indicator of failure
              })
            }
            
            if (mean(auc_vals, na.rm = TRUE) > best_auc_in_step) {
              best_auc_in_step <- mean(auc_vals, na.rm = TRUE)
              best_inner_model <- model
              best_inner_size <- grid$size[params]
              best_inner_decay <- grid$decay[params]
              best_inner_vars <- temp_vars
              improved <- TRUE
            }
          }
        }
        
        if (improved) {
          current_vars <- best_inner_vars
          best_inner_auc <- best_auc_in_step
          cat("Vars improved: ", paste(current_vars, collapse=", "), "\n")
          cat("AUC: ", paste(best_inner_auc), "\n")
        }
      }
      
      list(model = best_inner_model, auc = best_inner_auc, 
           size = best_inner_size, decay = best_inner_decay, 
           vars = best_inner_vars)
    }
    
    best_result <- stepAUC(variables, outer_train_data, target_name, grid, inner_folds)
    best_auc <- best_result$auc
    best_size <- best_result$size
    best_decay <- best_result$decay
    best_vars <- best_result$vars
    
    # Reentrenar el mejor modelo con todos los datos de outer_train_data
    best_model <- nnet(as.formula(paste(target_name, "~", paste(best_vars, collapse = "+"))),
                       data = outer_train_data, size = best_size, decay = best_decay, 
                       linout = FALSE, maxit = 200, trace = FALSE)
    
    predictions <- predict(best_model, outer_test_data, type = "raw")
    pred_class <- ifelse(predictions > threshold, 1, 0)  # Asignar clases basadas en umbral
    confusion <- table(Actual = outer_test_data[[target_name]], Predicted = pred_class)
    
    if (nrow(confusion) < 2 || ncol(confusion) < 2) {
      TP <- ifelse(nrow(confusion) >= 2 && ncol(confusion) >= 2, confusion[2, 2], 0)
      TN <- ifelse(nrow(confusion) >= 1 && ncol(confusion) >= 1, confusion[1, 1], 0)
      FP <- ifelse(nrow(confusion) >= 1 && ncol(confusion) >= 2, confusion[1, 2], 0)
      FN <- ifelse(nrow(confusion) >= 2 && ncol(confusion) >= 1, confusion[2, 1], 0)
    } else {
      TP <- confusion[2, 2]
      TN <- confusion[1, 1]
      FP <- confusion[1, 2]
      FN <- confusion[2, 1]
    }
    auc = roc(outer_test_data[[target_name]], predictions)$auc
    
    metrics <- setNames(c(outer_index, TP, TN, FP, FN, auc, 
                          best_size, best_decay, paste(best_vars, collapse = ",")), 
                        names(performance_metrics))
    performance_metrics <- rbind(performance_metrics, as.data.frame(t(metrics)))
    
    cat("Performance:", paste(performance_metrics))
    cat(sprintf("Confusion Matrix for Fold %d:\n", outer_index))
    print(confusion)
    cat(sprintf("Metrics for Fold %d:\n", outer_index))
    cat(sprintf("TP: %d, TN: %d, FP: %d, FN: %d, AUC: %f, Best Size: %d, Best Decay: %f, Best Variables: %s\n",
                TP, TN, FP, FN, auc, best_size, best_decay, paste(best_vars, collapse = ",")))
  }
  
  list(performance_metrics = performance_metrics, inner_fold_metrics = inner_fold_metrics)
}

evaluate_with_seeds_ANN_StepAUC <- 
  function(evaluation_function, data, target_var, vars, threshold, seeds, outer, inner, sizes, decays) {
  results_list <- list()

  for (seed in seeds) {
    set.seed(seed)
    cat("Evaluating with seed:", seed, "\n")

    results <- evaluation_function(data = data, 
                                   target_name = target_var, 
                                   variables = vars, 
                                   threshold = threshold, 
                                   outer = outer, 
                                   inner = inner, 
                                   seed = seed,
                                   sizes = sizes,
                                   decays = decays)
    results_list[[paste0("seed_", seed)]] <- results
  }

  return(results_list)
  }

performance_ANN_stepAUC <- evaluate_with_seeds_ANN_StepAUC(
  evaluation_function = double_cross_validation_nnet_stepAUC,
  data = data_factor,
  outer = 5,
  inner = 2,
  target_var = "PCR",
  sizes = c(5,10,15,20),
  decays = c(0.1, 0.2, 0.3),
  vars = todas_variables_factor,
  threshold = .35,
  seeds = seeds
)

load("RData_Files_Algorithms/ANN_StepAUC.RData")

3.7.5 Embedded (Lasso)

double_cross_validation_nnet_lasso <- function(data, target_name, outer, inner, sizes, decays, threshold, seed) {
  set.seed(seed)
  performance_metrics <- data.frame(Fold = integer(), TP = integer(), TN = integer(),
                                    FP = integer(), FN = integer(), AUC = numeric(),
                                    Best_Size = integer(), Best_Decay = numeric(), Best_Variables = character(), stringsAsFactors = FALSE)
  
  inner_fold_metrics <- data.frame(Fold = integer(), Inner_Fold = integer(), 
                                   TP = integer(), TN = integer(), FP = integer(), 
                                   FN = integer(), AUC = numeric(), stringsAsFactors = FALSE)
  
  # Convertir las variables del dataset con model.matrix
  target <- as.numeric(data[[target_name]]) - 1
  features <- data[, setdiff(names(data), target_name)]
  data_matrix <- as.data.frame(model.matrix(~ ., data=features)[,-1])  # Eliminar el intercepto
  data_matrix[[target_name]] <- target
  
  outer_folds <- createFolds(data_matrix[[target_name]], k = outer)
  
  for (outer_index in seq_along(outer_folds)) {
    outer_train_data <- data_matrix[-outer_folds[[outer_index]], ]
    outer_test_data <- data_matrix[outer_folds[[outer_index]], ]
    
    cat(sprintf("Outer Fold %d\n", outer_index))
    cat("Training Data - Class 0: ", sum(outer_train_data[[target_name]] == 0), "\n")
    cat("Training Data - Class 1: ", sum(outer_train_data[[target_name]] == 1), "\n")
    cat("Testing Data - Class 0: ", sum(outer_test_data[[target_name]] == 0), "\n")
    cat("Testing Data - Class 1: ", sum(outer_test_data[[target_name]] == 1), "\n")
    
    inner_folds <- createFolds(outer_train_data[[target_name]], k = inner)
    best_model <- NULL
    best_auc <- 0
    best_size <- NULL
    best_decay <- NULL
    best_vars <- colnames(outer_train_data)[!colnames(outer_train_data) %in% target_name]
    
    grid <- expand.grid(size = sizes, decay = decays)
    
    select_vars_lasso <- function(data, target_name) {
      target <- as.numeric(data[[target_name]])
      features <- data[, setdiff(names(data), target_name)]
      X <- model.matrix(~ ., data=features)[,-1]  # Eliminamos el intercepto
      y <- target
      
      # Ajustamos el modelo Lasso
      lasso_model <- cv.glmnet(X, y, alpha=1)
      lasso_coef <- coef(lasso_model, s = "lambda.min")
      
      # Seleccionamos las variables no nulas
      selected_vars <- rownames(lasso_coef)[lasso_coef[, 1] != 0]
      selected_vars <- setdiff(selected_vars, "(Intercept)")
      cat("Selected variables by Lasso: ", paste(selected_vars, collapse=", "), "\n")
      return(selected_vars)
    }
    
    stepAUC <- function(vars, data, target_name, grid, inner_folds) {
      best_inner_auc <- 0
      best_inner_model <- NULL
      best_inner_size <- NULL
      best_inner_decay <- NULL
      best_inner_vars <- vars
      improved <- TRUE
      
      while (improved && length(vars) > 1) {
        improved <- FALSE
        best_auc_in_step <- best_inner_auc
        
        for (params in 1:nrow(grid)) {
          auc_vals <- numeric()
          
          for (inner_index in seq_along(inner_folds)) {
            inner_train_data <- outer_train_data[-inner_folds[[inner_index]], ]
            inner_test_data <- outer_train_data[inner_folds[[inner_index]], ]
            
            selected_vars <- select_vars_lasso(inner_train_data, target_name)
            formula <- as.formula(paste(target_name, "~", paste(selected_vars, collapse = "+")))
            #cat("Testing with formula: ", deparse(formula), "\n")
            #cat("Grid parameters - size: ", grid$size[params], ", decay: ", grid$decay[params], "\n")
            
            result <- tryCatch({
              model <- nnet(formula, data = inner_train_data, size = grid$size[params], 
                            decay = grid$decay[params], linout = FALSE, maxit = 200, trace = FALSE)
              
              predictions <- predict(model, inner_test_data, type = "raw")
              roc_curve <- roc(inner_test_data[[target_name]], predictions)
              auc_vals[inner_index] <- roc_curve$auc
              #cat("AUC for inner fold ", inner_index, ": ", roc_curve$auc, "\n")
              
            }, error = function(e) {
              cat("Error with variables: ", paste(selected_vars, collapse = ", "), "\n")
              cat("Error message: ", e$message, "\n")
              auc_vals[inner_index] <- NA  # Set to NA or some other indicator of failure
            })
          }
          
          mean_auc <- mean(auc_vals, na.rm = TRUE)
          #cat("Mean AUC for params - size: ", grid$size[params], ", decay: ", grid$decay[params], ": ", mean_auc, "\n")
          if (!is.na(mean_auc) && mean_auc > best_auc_in_step) {
            best_auc_in_step <- mean_auc
            best_inner_model <- model
            best_inner_size <- grid$size[params]
            best_inner_decay <- grid$decay[params]
            best_inner_vars <- selected_vars
            improved <- TRUE
          }
        }
        
        if (improved) {
          vars <- best_inner_vars
          best_inner_auc <- best_auc_in_step
          cat("Vars improved: ", paste(vars, collapse=", "), "\n")
          cat("AUC: ", best_inner_auc, "\n")
        }
      }
      
      list(model = best_inner_model, auc = best_inner_auc, 
           size = best_inner_size, decay = best_inner_decay, 
           vars = best_inner_vars)
    }
    
    best_result <- stepAUC(best_vars, outer_train_data, target_name, grid, inner_folds)
    best_model <- best_result$model
    best_auc <- best_result$auc
    best_size <- best_result$size
    best_decay <- best_result$decay
    best_vars <- best_result$vars
    
    if (is.null(best_model)) {
      cat("Error: No valid model found for outer fold ", outer_index, "\n")
      next
    }

    # Evaluar el mejor modelo en el validation split del inner loop
    inner_auc_vals <- numeric()
    for (inner_index in seq_along(inner_folds)) {
      inner_train_data <- outer_train_data[-inner_folds[[inner_index]], ]
      inner_test_data <- outer_train_data[inner_folds[[inner_index]], ]
      
      predictions <- predict(best_model, inner_test_data, type = "raw")
      roc_curve <- roc(inner_test_data[[target_name]], predictions)
      inner_auc_vals[inner_index] <- roc_curve$auc
      #cat("Inner fold ", inner_index, " AUC: ", roc_curve$auc, "\n")

      pred_class <- ifelse(predictions > threshold, 1, 0)
      confusion <- table(Actual = inner_test_data[[target_name]], Predicted = pred_class)
      TP <- ifelse(nrow(confusion) >= 2 && ncol(confusion) >= 2, confusion[2, 2], 0)
      TN <- ifelse(nrow(confusion) >= 1 && ncol(confusion) >= 1, confusion[1, 1], 0)
      FP <- ifelse(nrow(confusion) >= 1 && ncol(confusion) >= 2, confusion[1, 2], 0)
      FN <- ifelse(nrow(confusion) >= 2 && ncol(confusion) >= 1, confusion[2, 1], 0)
      auc <- roc_curve$auc
      
      inner_metrics <- setNames(c(outer_index, inner_index, TP, TN, FP, FN, auc), names(inner_fold_metrics))
      inner_fold_metrics <- rbind(inner_fold_metrics, as.data.frame(t(inner_metrics)))
    }
    
    best_inner_auc <- mean(inner_auc_vals, na.rm = TRUE)
    #cat("Best inner AUC: ", best_inner_auc, "\n")
    # Reentrenar el mejor modelo con todos los datos de outer_train_data
    best_model <- nnet(as.formula(paste(target_name, "~", paste(best_vars, collapse = "+"))),
                       data = outer_train_data, size = best_size, decay = best_decay, 
                       linout = FALSE, maxit = 200, trace = FALSE)
    predictions <- predict(best_model, outer_test_data, type = "raw")
    pred_class <- ifelse(predictions > threshold, 1, 0)  # Asignar clases basadas en umbral
    confusion <- table(Actual = outer_test_data[[target_name]], Predicted = pred_class)
    
    if (nrow(confusion) < 2 || ncol(confusion) < 2) {
      TP <- ifelse(nrow(confusion) >= 2 && ncol(confusion) >= 2, confusion[2, 2], 0)
      TN <- ifelse(nrow(confusion) >= 1 && ncol(confusion) >= 1, confusion[1, 1], 0)
      FP <- ifelse(nrow(confusion) >= 1 && ncol(confusion) >= 2, confusion[1, 2], 0)
      FN <- ifelse(nrow(confusion) >= 2 && ncol(confusion) >= 1, confusion[2, 1], 0)
    } else {
      TP <- confusion[2, 2]
      TN <- confusion[1, 1]
      FP <- confusion[1, 2]
      FN <- confusion[2, 1]
    }
    auc <- roc(outer_test_data[[target_name]], predictions)$auc
    
    metrics <- setNames(c(outer_index, TP, TN, FP, FN, auc, 
                          best_size, best_decay, paste(best_vars, collapse = ",")), 
                        names(performance_metrics))
    performance_metrics <- rbind(performance_metrics, as.data.frame(t(metrics)))
    
    cat("Performance:\n", paste(performance_metrics), "\n")
    cat(sprintf("Confusion Matrix for Fold %d:\n", outer_index))
    print(confusion)
    cat(sprintf("Metrics for Fold %d:\n", outer_index))
    cat(sprintf("TP: %d, TN: %d, FP: %d, FN: %d, AUC: %f, Best Size: %d, Best Decay: %f, Best Variables: %s\n",
                TP, TN, FP, FN, auc, best_size, best_decay, paste(best_vars, collapse = ",")))
  }
  
  list(performance_metrics = performance_metrics, inner_fold_metrics = inner_fold_metrics)
}

evaluate_with_seeds_ANN_lasso <- 
  function(evaluation_function, data, target_var, vars, threshold, seeds, outer, inner, sizes, decays) {
  results_list <- list()

  for (seed in seeds) {
    set.seed(seed)
    cat("Evaluating with seed:", seed, "\n")

    results <- evaluation_function(data = data, 
                                   target_name = target_var, 
                                   threshold = threshold, 
                                   outer = outer, 
                                   inner = inner, 
                                   seed = seed,
                                   sizes = sizes,
                                   decays = decays)
    results_list[[paste0("seed_", seed)]] <- results
  }

  return(results_list)
  }

performance_ANN_lasso <- evaluate_with_seeds_ANN_lasso(
  evaluation_function = double_cross_validation_nnet_lasso,
  data = data_factor,
  outer = 5,
  inner = 2,
  target_var = "PCR",
  sizes = c(5,10,15,20),
  decays = c(0.1, 0.2, 0.3),
  threshold = .35,
  seeds = seeds
)

load("RData_Files_Algorithms/ANN_Lasso.RData")

4 Discusión

#  Función que combina los dataframes de métricas de rendimiento (performance_metrics) 
# e inner_fold_metrics de cada semilla en un solo dataframe, añadiendo una columna 
# "semilla" para diferenciar cada uno

process_experiment_data <- function(raw_data, hyperparameters) {
  # Combinar performance_metrics en un solo dataframe
  performance_metrics <- bind_rows(
    lapply(names(raw_data), function(seed) {
      df <- raw_data[[seed]]$performance_metrics
      df$seed <- seed
      return(df)
    })
  )
  
  # Convertir las columnas relevantes a tipo numérico
  performance_metrics <- performance_metrics %>%
    mutate(AUC = as.numeric(AUC),
           Fold = as.numeric(Fold),
           TP = as.numeric(TP),
           TN = as.numeric(TN),
           FP = as.numeric(FP),
           FN = as.numeric(FN),
           Seed = as.factor(seed))
  
  # Convertir los hiperparámetros relevantes a tipo numérico
  for (param in hyperparameters) {
    performance_metrics[[param]] <- as.numeric(performance_metrics[[param]])
  }
  
  # Combinar inner_fold_metrics en un solo dataframe
  inner_fold_metrics <- bind_rows(
    lapply(names(raw_data), function(seed) {
      df <- raw_data[[seed]]$inner_fold_metrics
      df$seed <- seed
      return(df)
    })
  )
  
  # Convertir las columnas relevantes a tipo numérico
  inner_fold_metrics <- inner_fold_metrics %>%
    mutate(AUC = as.numeric(AUC),
           Fold = as.numeric(Fold),
           Inner_Fold = as.numeric(Inner_Fold),
           TP = as.numeric(TP),
           TN = as.numeric(TN),
           FP = as.numeric(FP),
           FN = as.numeric(FN),
           Seed = as.factor(seed))
  
  inner_fold_metrics$seed <- NULL
  performance_metrics$seed <- NULL
  
  return(list(performance_metrics = performance_metrics, inner_fold_metrics = inner_fold_metrics))
}

# ANN
hyperparameters_ann <- c("Best_Size", "Best_Decay")

ANN_TodasVariables_Metrics <- process_experiment_data(performance_ANN_todasvariables, hyperparameters_ann)
ANN_Asociacion_Metrics <- process_experiment_data(performance_ANN_asociacion, hyperparameters_ann)
ANN_StepAUC_Metrics <- process_experiment_data(performance_ANN_stepAUC, hyperparameters_ann)
ANN_Lasso_Metrics <- process_experiment_data(performance_ANN_lasso, hyperparameters_ann)

# SVM
hyperparameters_svm <- c("Best_Cost", "Best_Kernel")

SVM_TodasVariables_Metrics <- process_experiment_data(performance_SVM_todasvariables, hyperparameters_svm)
SVM_Asociacion_Metrics <- process_experiment_data(performance_svm_asociacion,hyperparameters_svm)
SVM_StepAUC_Metrics <- process_experiment_data(performance_SVM_stepauc, hyperparameters_svm)
SVM_Lasso_Metrics <- process_experiment_data(performance_SVM_lasso, hyperparameters_svm)


# KNN
hyperparameters_knn <- c("Best_k", "Best_l")

KNN_TodasVariables_Metrics <- process_experiment_data(performance_knn_todasvariables, hyperparameters_knn)
KNN_Asociacion_Metrics <- process_experiment_data(performance_knn_Asociacion, hyperparameters_knn)
KNN_StepAUC_Metrics <- process_experiment_data(performance_knn_stepauc, hyperparameters_knn)
KNN_Lasso_Metrics <- process_experiment_data(performance_knn_lasso, hyperparameters_knn)


# RandomForest
hyperparameters_RandomForest <- c("Best_ntree", "Best_mtry")

RandomForest_TodasVariables_Metrics <- process_experiment_data(performance_RandomForest_todasvariables, hyperparameters_RandomForest)
RandomForest_Asociacion_Metrics <- process_experiment_data(performance_RandomForest_Asociacion, hyperparameters_RandomForest)
RandomForest_StepAUC_Metrics <- process_experiment_data(performance_RandomForest_StepAUC, c("Best_Ntree", "Best_Mtry"))

RandomForest_StepAUC_Metrics$performance_metrics <- RandomForest_StepAUC_Metrics$performance_metrics %>% rename(Best_ntree = Best_Ntree)
RandomForest_StepAUC_Metrics$performance_metrics <- RandomForest_StepAUC_Metrics$performance_metrics %>% rename(Best_mtry = Best_Mtry)


# NaiveBayes
hyperparameters_NaiveBayes <- c()

NaiveBayes_TodasVariables_Metrics <- process_experiment_data(performance_naiveBayes_todasVariables, hyperparameters_NaiveBayes)
NaiveBayes_Asociacion_Metrics <- process_experiment_data(performance_naiveBayes_Asociacion, hyperparameters_NaiveBayes)
NaiveBayes_StepAUC_Metrics <- process_experiment_data(performance_naiveBayes_stepAUC, hyperparameters_NaiveBayes)
NaiveBayes_Lasso_Metrics <- process_experiment_data(performance_naiveBayes_lasso, hyperparameters_NaiveBayes)


# AdaBoost
hyperparameters_AdaBoost <- c("Best_Iteration")

AdaBoost_TodasVariables_Metrics <- process_experiment_data(performance_adaboost_todasvariables, hyperparameters_AdaBoost)
AdaBoost_Asociacion_Metrics <- process_experiment_data(performance_adaboost_Asociacion, hyperparameters_AdaBoost)
AdaBoost_StepAUC_Metrics <- process_experiment_data(performance_adaboost_stepAUC, c("Best_Iterations"))
AdaBoost_Lasso_Metrics <- process_experiment_data(performance_adaboost_lasso, hyperparameters_AdaBoost)

AdaBoost_StepAUC_Metrics$performance_metrics <- AdaBoost_StepAUC_Metrics$performance_metrics %>% rename(Best_Iteration = Best_Iterations)

# RPart
hyperparameters_RPart <- c("Best_cp", "Best_minsplit")

RPart_TodasVariables_Metrics <- process_experiment_data(performance_rpart_todasvariables, hyperparameters_RPart)
RPart_Asociacion_Metrics <- process_experiment_data(performance_rpart_Asociacion, hyperparameters_RPart)
RPart_StepAUC_Metrics <- process_experiment_data(performance_rpart_StepAUC, c("Best_CP", "Best_Minsplit"))
RPart_Lasso_Metrics <- process_experiment_data(performance_rpart_lasso, hyperparameters_RPart)

RPart_StepAUC_Metrics$performance_metrics <- RPart_StepAUC_Metrics$performance_metrics %>% rename(Best_cp = Best_CP)
RPart_StepAUC_Metrics$performance_metrics <- RPart_StepAUC_Metrics$performance_metrics %>% rename(Best_minsplit = Best_Minsplit)

El primer paso que se seguirá consiste en estandarizar los datos de los experimentos realizados en el apartado Resultados de tal forma que tengamos un único objeto data.frame por cada método de selección de variables. Este objeto data.frame contentrá métricas de rendimiento del modelo (TP, TN, AUC, etc.), las variables seleccionadas por el método aplicado y los hiperparámetros que han dado lugar al mejor rendimiento del modelo, estos datos se repiten por cada outer fold, por cada semilla aplicada, de tal forma que, como hemos aplicado 5 semillas (27, 3, 123, 50, 90), y hemos realizado una doble validación cruzada 5x2, tendremos un data.frame de 25 filas por cada método de selección de variables aplicado. Por otra parte, también se han estandarizado de la misma forma las métricas de rendmiento de la evaluación del mejor modelo encontrado en los folds internos con los conjuntos de validación de los folds internos. Estos datos servirán para validar que no haya habido sobre entrenamientos.

A continuación se muestra una porción de cada uno de estos data.frame:

Conjunto performance_metrics de los Folds Externos

Se muestra el data.frame generado con los datos de rendimiento para los folds externos por semillas y folds para el método wrapped StepAUC para las redes de neuronas artificiales.

knitr::kable(head(ANN_StepAUC_Metrics$performance_metrics, 10))

Fold	TP	TN	FP	FN	AUC	Best_Size	Best_Decay	Best_Variables	Seed
1	10	66	14	11	0.7160714	5	0.3	Edad,REst,Estadio,NodAfec,Grado,Fenotipo	seed_27
2	3	79	5	15	0.8121693	10	0.3	Edad,REst,Grado	seed_27
3	6	74	11	10	0.7257353	5	0.3	Edad,NodAfec,Grado,Fenotipo	seed_27
4	6	70	9	17	0.7660980	20	0.3	Edad,REst,Estadio,NodAfec,Grado,Fenotipo	seed_27
5	10	62	19	11	0.6740153	5	0.3	Edad,REst,NodAfec,Grado,Fenotipo	seed_27
1	10	66	20	5	0.7953488	5	0.3	REst,Estadio,Grado,Fenotipo	seed_3
2	8	68	11	15	0.7438085	5	0.3	Edad,REst,Grado,Fenotipo	seed_3
3	11	59	19	13	0.6907051	10	0.1	REst,RPro,Estadio,Grado,Fenotipo	seed_3
4	5	73	12	11	0.7481618	20	0.3	Edad,REst,RPro,Estadio,Fenotipo	seed_3
5	7	74	7	14	0.7589653	10	0.3	Edad,REst,RPro,Estadio,Grado,Fenotipo	seed_3

Conjunto inner_performance de los Folds Internos

Se muestra el data.frame generado con los datos de rendimiento para los folds iternos por semillas y folds para el método wrapped StepAUC para las redes de neuronas artificiales.

knitr::kable(head(ANN_StepAUC_Metrics$inner_fold_metrics, 15))

Fold	Inner_Fold	TP	TN	FP	FN	AUC	Seed
1	1	16	149	14	24	0.7951687	seed_27
1	2	9	153	13	29	0.7770292	seed_27
2	1	0	170	0	33	0.7796791	seed_27
2	2	0	155	0	48	0.7002016	seed_27
3	1	13	154	10	27	0.7945122	seed_27
3	2	6	153	7	37	0.7715116	seed_27
4	1	0	169	0	34	0.7836756	seed_27
4	2	0	161	0	42	0.7758060	seed_27
5	1	8	157	7	31	0.8074578	seed_27
5	2	3	156	8	36	0.7578174	seed_27
1	1	0	164	0	40	0.7891006	seed_3
1	2	0	159	0	44	0.7339909	seed_3
2	1	24	133	25	21	0.7802391	seed_3
2	2	13	145	27	18	0.7815079	seed_3
3	1	19	138	29	17	0.8119594	seed_3

4.1 Estabilidad de AUC entre Folds

Esta sección se dedica al estudio de la estabilidad de las métricas de rendimiento para todos los modelos. Una de las principales características deseadas en un modelo de Aprendizaje Automático es la robustez en sus predicciones. Un modelo puede alcanzar un AUC máximo alto; sin embargo, si la varianza del AUC correspondiente a sus predicciones es muy alta, esto indica que el modelo es demasiado sensible a los datos de entrada, lo que podría resultar en predicciones inesperadas.

Para abordar esta cuestión, se llevará a cabo un estudio comparativo para cada método de selección de variables utilizado por cada modelo, considerando las diferentes semillas y folds externos aplicados en la validación interna.

Un modelo será considerado óptimo en esta sección si logra un alto AUC medio y si la variabilidad de sus AUC en función de la semilla y el fold es baja. Además, se evaluará si todos los métodos de selección de variables muestran un patrón consistente de AUC a lo largo de los diferentes folds, independientemente de la semilla. Esto indicaría que las fluctuaciones observadas se deben a la divisiones de los datos y no a la estructura del modelo, señalando una limitación impuesta por las muestras utilizadas.

A continuación, se comentan las gráficas que serán generadas para el estudio comparativo

AUC Medio entre Inner Folds y Outer Folds

El objetivo del gráfico comparativo del AUC medio de los inner folds y los outer folds es visualizar y comparar el rendimiento del modelo durante el proceso de validación cruzada. Específicamente, el gráfico de barras muestra los valores promedio del AUC tanto para los inner folds como para los outer folds, agrupados por fold. Al hacerlo, permite evaluar si el rendimiento del modelo es consistente entre las dos fases de validación y detectar posibles discrepancias o sobreajustes. La visualización ayuda a entender cómo el modelo generaliza a datos no vistos y proporciona una medida de estabilidad y confianza en el rendimiento del modelo.

# Gráfico de barras comparativo del AUC medio de los inner folds y el AUC medio de los outer folds
generate_inner_outer_auc_comparison <- function(inner_fold_metrics, outer_fold_metrics) {
  # Calcular AUC medio por fold para inner folds
  mean_inner_auc <- inner_fold_metrics %>%
    group_by(Fold) %>%
    summarize(mean_AUC_inner = mean(AUC, na.rm = TRUE))
  
  # Calcular AUC medio por fold para outer folds
  mean_outer_auc <- outer_fold_metrics %>%
    group_by(Fold) %>%
    summarize(mean_AUC_outer = mean(AUC, na.rm = TRUE))
  
  # Unir los dataframes por Fold
  comparison_data <- merge(mean_inner_auc, mean_outer_auc, by = "Fold")
  
  # Crear la gráfica de barras comparativa
  fig <- plot_ly(data = comparison_data, x = ~Fold) %>%
    add_trace(
      y = ~mean_AUC_inner, 
      type = 'bar', 
      name = 'AUC Medio Inner Folds',
      marker = list(color = 'blue')
    ) %>%
    add_trace(
      y = ~mean_AUC_outer, 
      type = 'bar', 
      name = 'AUC Medio Outer Folds',
      marker = list(color = 'red')
    ) %>%
    layout(
      title = 'Comparación de AUC Medio: Inner Folds vs Outer Folds',
      xaxis = list(title = 'Fold'),
      yaxis = list(title = 'AUC Medio', range = c(0, 1)),
      barmode = 'group'
    )
  
  return(fig)
}

AUC Medio por Outer Fold

El objetivo de la gráfica de violín del AUC medio por outer fold es proporcionar una visualización clara y detallada de la distribución del AUC para cada fold, distinguiendo los resultados por diferentes semillas de validación. Esta visualización es útil para detectar si hay algún fold específico que muestra una mayor variabilidad o si el rendimiento del modelo es consistentemente estable a través de todos los folds y semillas. Además, los elementos de la gráfica de violín, como la línea de la media y los puntos de datos individuales, facilitan la comprensión de la distribución del AUC y de los posibles outliers en los resultados.

#  Gráfico de violín del AUC medio por outer fold con colores distintos para cada fold
generate_outer_auc_violin_plot <- function(outer_fold_metrics) {
  # Colores para cada fold
  colors <- RColorBrewer::brewer.pal(n = length(unique(outer_fold_metrics$Fold)), name = "Set1")
  
  # Crear la gráfica de violín
  fig <- plot_ly(data = outer_fold_metrics, y = ~AUC, x = ~factor(Fold), type = 'violin', 
                 box = list(visible = TRUE), meanline = list(visible = TRUE), points = 'all', spanmode = 'hard',
                 color = ~factor(Fold), colors = colors) %>%
    layout(
      title = 'Distribución del AUC por Outer Fold',
      xaxis = list(title = 'Fold'),
      yaxis = list(title = 'AUC')
    )
  
  return(fig)
}

Comparativa entre los métodos de Selección de Variables

El objetivo de esta gráfica 3D es comparar la eficacia de diferentes métodos de selección de características (Todas Variables, Asociación, StepAUC y Lasso) en términos de su rendimiento medido por el AUC (Area Under the Curve). La visualización representa los resultados de múltiples ejecuciones (diferentes semillas) y los pliegues (folds) de validación cruzada para cada método. Cada punto en la gráfica indica el AUC obtenido en un pliegue específico para una semilla dada, y estos puntos se conectan con líneas para mostrar la continuidad dentro del mismo método y semilla. Los diferentes métodos se distinguen por colores, y en la leyenda se muestra el AUC medio global para cada método, proporcionando una comparación clara y visual de cómo cada método performa en promedio y cómo varía su rendimiento a través de diferentes pliegues y ejecuciones. Esto facilita la identificación de métodos consistentes y estables, así como aquellos que pueden necesitar ajustes adicionales.

generate_combined_auc_3d_plot <- function(dataframes, names) {
  # Definir los colores para cada dataframe
  colors <- brewer.pal(n = length(dataframes), name = "Set1")
  
  # Calcular el AUC promedio global para cada dataframe
  mean_aucs <- sapply(dataframes, function(df) mean(df$AUC, na.rm = TRUE))
  
  # Crear las trazas 3D para cada dataframe y semilla
  fig <- plot_ly()
  
  for (i in seq_along(dataframes)) {
    df <- dataframes[[i]]
    unique_seeds <- unique(df$Seed)
    for (seed in unique_seeds) {
      seed_data <- df %>% filter(Seed == seed)
      fig <- fig %>%
        add_trace(data = seed_data, x = ~Fold, y = ~Seed, z = ~AUC, type = 'scatter3d', mode = 'lines+markers',
                  line = list(color = colors[i]), marker = list(size = 5, color = colors[i]),
                  showlegend = ifelse(seed == unique_seeds[1], TRUE, FALSE),
                  name = ifelse(seed == unique_seeds[1], paste(names[i], "- AUC medio global:", round(mean_aucs[i], 3)), ""))
    }
  }
  
  fig <- fig %>%
    layout(
      title = 'Comparación de AUC por Semilla y Fold para Diferentes Métodos',
      scene = list(
        xaxis = list(title = 'Fold'),
        yaxis = list(title = 'Semilla'),
        zaxis = list(title = 'AUC')
      )
    )
  
  return(fig)
}

Adicionalmente, se mostrará la distribución de AUC a través de Folds y Semillas. El objetivo de esta gráfica es proporcionar una visualización tridimensional detallada de la distribución del AUC a través de diferentes folds y semillas para un modelo específico. La gráfica permite observar cómo varía el rendimiento del modelo (medido por el AUC) en cada fold y cómo esta variación se distribuye a lo largo de varias ejecuciones independientes (diferentes semillas). Cada color representa una semilla distinta, lo que facilita la identificación de patrones de estabilidad y variabilidad en el rendimiento del modelo. Es gráfica se dejará de manera complementaria a la comparativa general, para que el lector pueda observar las diferencias específicas para cada modelo, sin embargo, no se comentará.

# Función para generar la distribución de AUC a través de folds para todas las semillas usadas de un modelo
generate_auc_per_seed_and_fold <- function(performance_metrics) {
  # Calcular el AUC medio por semilla
  mean_auc_by_seed <- performance_metrics %>%
    group_by(Seed) %>%
    summarize(mean_AUC = mean(AUC, na.rm = TRUE))
  
  # Crear la visualización en 3D con plotly
  fig <- plot_ly()
  
  # Colores para cada semilla
  colors <- RColorBrewer::brewer.pal(n = length(unique(performance_metrics$Seed)), name = "Set1")
  
  # Añadir los datos de cada semilla como histogramas
  for (i in seq_along(unique(performance_metrics$Seed))) {
    seed <- unique(performance_metrics$Seed)[i]
    subset_data <- performance_metrics %>% filter(Seed == seed)
    
    fig <- fig %>%
      add_trace(
        x = subset_data$Fold, 
        y = rep(i, nrow(subset_data)), 
        z = subset_data$AUC, 
        type = 'scatter3d',
        mode = 'lines+markers',
        line = list(width = 2, color = colors[i]),
        marker = list(size = 5, color = colors[i]),
        name = paste('Semilla', gsub("seed_", "", seed), 'AUC medio:', round(mean_auc_by_seed$mean_AUC[i], 3))
      )
  }
  
  # Añadir detalles al gráfico
  fig <- fig %>% layout(scene = list(
    xaxis = list(title = 'Fold'),
    yaxis = list(title = 'Semilla'),
    zaxis = list(title = 'AUC'),
    title = "Visualización 3D de Folds, Semillas y AUC"
  ))
  
  # Mostrar el gráfico
  return(fig)
}

4.1.1 Naive Bayes

4.1.1.1 AUC Fold Interno vs Externo

TodasVariables

generate_inner_outer_auc_comparison(NaiveBayes_TodasVariables_Metrics$inner_fold_metrics, NaiveBayes_TodasVariables_Metrics$performance_metrics)

Asociacion

generate_inner_outer_auc_comparison(NaiveBayes_Asociacion_Metrics$inner_fold_metrics, NaiveBayes_Asociacion_Metrics$performance_metrics)

StepAUC

generate_inner_outer_auc_comparison(NaiveBayes_StepAUC_Metrics$inner_fold_metrics, NaiveBayes_StepAUC_Metrics$performance_metrics)

Lasso

generate_inner_outer_auc_comparison(NaiveBayes_Lasso_Metrics$inner_fold_metrics, NaiveBayes_Lasso_Metrics$performance_metrics)

El valor de AUC medio para todas las semillas por cada Fold es excepcional con el método del algoritmo de Naive Bayes. Como se puede observar, independientemente del método de selección de características aplicado, el rendimiento en los folds internos y externos es prácticamente el mismo, sin que el del externo supere al interno, que podría ser un signo de que el conjunto de datos ha sido mal construído, o de overfitting. En concreto se deben destacar los métodos de TodasVariables y StepAUC, por mantener el AUC promedio más alto, y el método de Asociación, por mantener el mismo rendimiento en los folds internos y externos, exhibiendo una gran capacidad de generalización y predicción.

4.1.1.2 AUC Medio Folds Externos

TodasVariables

generate_outer_auc_violin_plot(NaiveBayes_TodasVariables_Metrics$performance_metrics)

Asociacion

generate_outer_auc_violin_plot(NaiveBayes_Asociacion_Metrics$performance_metrics)

StepAUC

generate_outer_auc_violin_plot(NaiveBayes_StepAUC_Metrics$performance_metrics)

Lasso

generate_outer_auc_violin_plot(NaiveBayes_Lasso_Metrics$performance_metrics)

El método más destacable de este apartado podría ser el método de selección de características Lasso, ya que la dispersión del AUC se mantiene muy baja en los folds 1 y 5, y bastante estable en los folds 2 y 3, con la excepción del cuarto fold, el cual parece mostrar un valor de AUC bajo. La dispersión del resto de métodos es, en general, alta, sin embargo, hay que destacar que los métodos TodasVariables y StepAuc tienen la mediana en el valor de AUC alto, indicando que el, o los datos que causan la dispersión, pueden ser excepciones.

4.1.1.3 Comparativa Naive Bayes

fig <- generate_combined_auc_3d_plot(
  list(
    NaiveBayes_TodasVariables_Metrics$performance_metrics,
    NaiveBayes_Asociacion_Metrics$performance_metrics,
    NaiveBayes_StepAUC_Metrics$performance_metrics,
    NaiveBayes_Lasso_Metrics$performance_metrics
  ),
  c("Todas Variables", "Asociacion", "StepAUC", "Lasso")
)

fig

Los valores de AUC, para todos los métodos de selección de características aplicados, muestran una notable estabilidad, pareciendo que, independientemente de la semilla, las fluctuaciones entre folds se mantienen, indicando que el modelo aprende bien de los datos, siendo robusto a la incorporación de nuevos datos, y realizando buenas predicciones a futuro. La técnica que mejores métricas muestra es la de TodasVariables, seguida muy cerca por StepAUC, la cual muestra la ventaja de usar menos cantidad de características.

TodasVariables

generate_auc_per_seed_and_fold(NaiveBayes_TodasVariables_Metrics$performance_metrics)

Asociacion

generate_auc_per_seed_and_fold(NaiveBayes_Asociacion_Metrics$performance_metrics)

StepAUC

generate_auc_per_seed_and_fold(NaiveBayes_StepAUC_Metrics$performance_metrics)

Lasso

generate_auc_per_seed_and_fold(NaiveBayes_Lasso_Metrics$performance_metrics)

4.1.2 AdaBoost

4.1.2.1 AUC Fold Interno vs Externo

TodasVariables

generate_inner_outer_auc_comparison(AdaBoost_TodasVariables_Metrics$inner_fold_metrics, AdaBoost_TodasVariables_Metrics$performance_metrics)

Asociacion

generate_inner_outer_auc_comparison(AdaBoost_Asociacion_Metrics$inner_fold_metrics, AdaBoost_Asociacion_Metrics$performance_metrics)

StepAUC

generate_inner_outer_auc_comparison(AdaBoost_StepAUC_Metrics$inner_fold_metrics, AdaBoost_StepAUC_Metrics$performance_metrics)

Lasso

generate_inner_outer_auc_comparison(AdaBoost_Lasso_Metrics$inner_fold_metrics, AdaBoost_Lasso_Metrics$performance_metrics)

Se puede observar que todos los métodos de selección de características mantienen un rendimiento mayor en los folds internos que externos, asegurando así que el conjunto de datos está bien construído. Sin embargo, el método donde más pareados parecen estar estos dos valores, y por ello donde el modelo final podría tener una mayor capacidad de generalización, es en el método de selección de variables Lasso.

4.1.2.2 AUC Medio Folds Externos

TodasVariables

generate_outer_auc_violin_plot(AdaBoost_TodasVariables_Metrics$performance_metrics)

Asociacion

generate_outer_auc_violin_plot(AdaBoost_Asociacion_Metrics$performance_metrics)

StepAUC

generate_outer_auc_violin_plot(AdaBoost_StepAUC_Metrics$performance_metrics)

Lasso

generate_outer_auc_violin_plot(AdaBoost_Lasso_Metrics$performance_metrics)

En el método TodasVariables se puede observar una alta dispersión en todos los folds, con el cuarto siendo el que menos muestra, pero para el cual la mediana es de menor valor. Para el método Asociación, se puede observar que el segundo fold muestra una muy baja dispersión y un valor de AUC medio alto, sin embargo, los folds 1 y 5 muestran características totalmente opuestas. El método stepAUC podría ser una mezcla entre los dos anteriores, mostrando una dispersión hacia valores altos de AUC, pero, debido a que la mayoría de muestras caen en valores algo más mediocres de AUC, la media está más centrada en la escala. Por último, El método Lasso muestra una estabilidad en la dispersión de los AUC, ya que todos los folds están igual de dispersos, manteniendo las medianas en los mismos valores aproximadamente, sin embargo, también parecer ser el método con más outliers. Podemos concluir que parece que los métodos Lasso y Asociación son algo más apropiados en cuanto a la robustez de los datos de entrada al modelo.

4.1.2.3 Comparativa AdaBoost

fig <- generate_combined_auc_3d_plot(
  list(
    AdaBoost_TodasVariables_Metrics$performance_metrics,
    AdaBoost_Asociacion_Metrics$performance_metrics,
    AdaBoost_StepAUC_Metrics$performance_metrics,
    AdaBoost_Lasso_Metrics$performance_metrics
  ),
  c("Todas Variables", "Asociacion", "StepAUC", "Lasso")
)

fig

Se puede observar que la dispersión de las medidas de AUC a través de los Folds no es demasiado alta, manteniendo estos valores entre .65 y .8, adicionalmente, parece que todos los puntos siguen un patrón a lo largo de las semillas para los folds generados, transmitiendo que el modelo es robusto a la aleatoriedad de la división de los datos, y concluyendo que es robusto ante nuevos datos, con una buena capacidad de predicción. Estos patrones que se repiten pueden incluir una tendencia a la disminución del AUC en el fold 1, una inestabilidad en el fold 2, una estabilidad en el 3, y una tendencia a la mejora en el 4 y 5. Por otra aprte, el modelo que parece haber obtenido mejores métricas es el que ha empleado la selección de características Lasso, concordando con lo concluído en las anteriores secciones.

TodasVariables

generate_auc_per_seed_and_fold(AdaBoost_TodasVariables_Metrics$performance_metrics)

Asociacion

generate_auc_per_seed_and_fold(AdaBoost_Asociacion_Metrics$performance_metrics)

StepAUC

generate_auc_per_seed_and_fold(AdaBoost_StepAUC_Metrics$performance_metrics)

Lasso

generate_auc_per_seed_and_fold(AdaBoost_Lasso_Metrics$performance_metrics)

4.1.3 RandomForest

4.1.3.1 AUC Fold Interno vs Externo

TodasVariables

generate_inner_outer_auc_comparison(RandomForest_TodasVariables_Metrics$inner_fold_metrics, RandomForest_TodasVariables_Metrics$performance_metrics)

Asociacion

generate_inner_outer_auc_comparison(RandomForest_Asociacion_Metrics$inner_fold_metrics, RandomForest_Asociacion_Metrics$performance_metrics)

StepAUC

generate_inner_outer_auc_comparison(RandomForest_StepAUC_Metrics$inner_fold_metrics, RandomForest_StepAUC_Metrics$performance_metrics)

Se puede observar que el rendimiento en para los Folds externos es notablemente más bajo que el rendimiento para los internos a lo largo de todos los métodos de selección de características. Este hecho denota que el modelo es estable y no se han cometido errores a la hora de evaluarlo, sin embargo, también resalta la pobre capacidad de generalización, incluso habiendo ajustado los hiperparámetros a sus mejores valores.

4.1.3.2 AUC Medio Folds Externos

TodasVariables

generate_outer_auc_violin_plot(RandomForest_TodasVariables_Metrics$performance_metrics)

Asociacion

generate_outer_auc_violin_plot(RandomForest_Asociacion_Metrics$performance_metrics)

StepAUC

generate_outer_auc_violin_plot(RandomForest_StepAUC_Metrics$performance_metrics)

En el método TodasVariables se puede observar que la dispersión del Fold 1 es muy baja y centrada en un valor alto del AUC, sin embargo, el resto de folds tienen una varianza alta, demostrando una baja capacidad de generalización y predicción. Por otra parte, el método Asociación muestra una baja dispersión en los Folds 1, 2 y 5, con alguna excepción en los Folds 3 y 4, adicionalmente, todos los Folds están centrados alrededor de un valor de AUC común, que se sitúa entre el 0.62 y el 0.68, mostrando una estabilidad. Finalmente, el método StepAUC muestra una alta dispersión en el primer y último Fold, pero un alto valor de AUC promedio, al igual que el primer método. Se puede concluir que, para este modelo, el método de selección de características más adecuado podría ser el de TodasVariables.

4.1.3.3 Comparativa RandomForest

fig <- generate_combined_auc_3d_plot(
  list(
    RandomForest_TodasVariables_Metrics$performance_metrics,
    RandomForest_Asociacion_Metrics$performance_metrics,
    RandomForest_StepAUC_Metrics$performance_metrics
  ),
  c("Todas Variables", "Asociacion", "StepAUC")
)

fig

Los métodos de selección más competitivos son StepAUC y TodasVariables, con AUC’s de 0.677 y 0.682, respectivamente. Parece que, aunque el StepAUC muestre un mayor AUC global, la estabilidad del método hace que la dispersión del AUC sea más alta, teniendo un mínimo de 0.5 en la semila 27 para el Fold 1. Se podría concluir que el método más adecuado para el rendimiento con respecto al AUC del RandomForest es el de TodasVariables.

TodasVariables

generate_auc_per_seed_and_fold(RandomForest_TodasVariables_Metrics$performance_metrics)

Asociacion

generate_auc_per_seed_and_fold(RandomForest_Asociacion_Metrics$performance_metrics)

StepAUC

generate_auc_per_seed_and_fold(RandomForest_StepAUC_Metrics$performance_metrics)

4.1.4 SVM

4.1.4.1 AUC Fold Interno vs Externo

TodasVariables

generate_inner_outer_auc_comparison(SVM_TodasVariables_Metrics$inner_fold_metrics, SVM_TodasVariables_Metrics$performance_metrics)

Asociacion

generate_inner_outer_auc_comparison(SVM_Asociacion_Metrics$inner_fold_metrics, SVM_Asociacion_Metrics$performance_metrics)

StepAUC

generate_inner_outer_auc_comparison(SVM_StepAUC_Metrics$inner_fold_metrics, SVM_StepAUC_Metrics$performance_metrics)

Lasso

generate_inner_outer_auc_comparison(SVM_Lasso_Metrics$inner_fold_metrics, SVM_Lasso_Metrics$performance_metrics)

En el método de todas las variables, el AUC medio en los inner folds es consistentemente mayor que en los outer folds. El método de todas las variables muestra valores notablemente bajos de AUC en el fold externo que en el interno, indicando estabilidad interna pero una pobre generalización. En el método de asociación, las diferencias entre los AUC medios de inner y outer folds son menores, sugiriendo una mejor generalización y consistencia en el rendimiento del modelo, llegando, de manera puntual y seguramente por outliers, a sobrepasar el rendimiento con test al rendimiento con la validación. El método stepAUC muestra una variabilidad similar a la de todas las variables, con inner folds generalmente superiores, pero con menos diferencia en comparación. Finalmente, el método Lasso presenta una paridad notable entre los AUC de inner y outer folds, reflejando una mayor estabilidad y generalización, con los valores más cercanos entre sí. En resumen, los métodos de asociación y Lasso parecen ofrecer una mejor generalización y estabilidad en el rendimiento del modelo en comparación con las técnicas que utilizan todas las variables y stepAUC.

4.1.4.2 AUC Medio Folds Externos

TodasVariables

generate_outer_auc_violin_plot(SVM_TodasVariables_Metrics$performance_metrics)

Asociacion

generate_outer_auc_violin_plot(SVM_Asociacion_Metrics$performance_metrics)

StepAUC

generate_outer_auc_violin_plot(SVM_StepAUC_Metrics$performance_metrics)

Lasso

generate_outer_auc_violin_plot(SVM_Lasso_Metrics$performance_metrics)

En el método de todas las variables, se observa una mayor variabilidad en los valores de AUC, especialmente en el Fold 1, indicando inestabilidad en el rendimiento del modelo. El método de asociación presenta una distribución más compacta en general, con menos variabilidad, especialmente en los folds 1, 3 y 5, lo que sugiere una mayor consistencia en el rendimiento. El método stepAUC muestra una mayor variabilidad en el Fold 2, con un rango amplio de AUC, y menor consistencia en los folds 4 y 5. Finalmente, el método Lasso exhibe una distribución más equilibrada a través de todos los folds, con un AUC medio más alto, pero con una dispersión muy grande a pesar de todo. En conjunto, podríamos concluir que, los métodos de asociación y, de manera secundaria, Lasso demuestran un rendimiento más confiable y menos variable en comparación con todas las variables y stepAUC.

4.1.4.3 Comparativa SVM

fig <- generate_combined_auc_3d_plot(
  list(
    SVM_TodasVariables_Metrics$performance_metrics,
    SVM_Asociacion_Metrics$performance_metrics,
    SVM_StepAUC_Metrics$performance_metrics,
    SVM_Lasso_Metrics$performance_metrics
  ),
  c("Todas Variables", "Asociacion", "StepAUC", "Lasso")
)

fig

Las gráficas comparativas del AUC para todos los métodos muestran una notable inestabilidad en el rendimiento del SVM independientemente del método de selección de características utilizado. Aunque el método Lasso presenta el AUC medio global más alto (0.64), seguido por la asociación (0.623), stepAUC (0.618) y todas las variables (0.617), en general, se observa una considerable variabilidad a través de diferentes semillas y folds. En todas las gráficas, los valores de AUC fluctúan significativamente, especialmente en los métodos de todas las variables y stepAUC, donde los puntos individuales muestran amplios rangos de dispersión. Incluso los métodos de selección de características que generalmente mejoran la estabilidad, como Lasso y asociación, no eliminan completamente la variabilidad, sugiriendo que el SVM es inherentemente inestable en este contexto. Esta inestabilidad destaca la dificultad del SVM para mantener un rendimiento consistente, lo que puede estar influenciado por la naturaleza de los datos y la sensibilidad del modelo a las particiones específicas de entrenamiento y prueba.

TodasVariables

generate_auc_per_seed_and_fold(SVM_TodasVariables_Metrics$performance_metrics)

Asociacion

generate_auc_per_seed_and_fold(SVM_Asociacion_Metrics$performance_metrics)

StepAUC

generate_auc_per_seed_and_fold(SVM_StepAUC_Metrics$performance_metrics)

Lasso

generate_auc_per_seed_and_fold(SVM_Lasso_Metrics$performance_metrics)

4.1.5 ANN

4.1.5.1 AUC Fold Interno vs Externo

TodasVariables

generate_inner_outer_auc_comparison(ANN_TodasVariables_Metrics$inner_fold_metrics, ANN_TodasVariables_Metrics$performance_metrics)

Asociacion

generate_inner_outer_auc_comparison(ANN_Asociacion_Metrics$inner_fold_metrics, ANN_Asociacion_Metrics$performance_metrics)

StepAUC

generate_inner_outer_auc_comparison(ANN_StepAUC_Metrics$inner_fold_metrics, ANN_StepAUC_Metrics$performance_metrics)

Lasso

generate_inner_outer_auc_comparison(ANN_Lasso_Metrics$inner_fold_metrics, ANN_Lasso_Metrics$performance_metrics)

Todas las gráficas muestran que en ningún momento el AUC del outer fold supera al del inner fold, lo que indica que no ha habido errores en la división de datos, confirmando la correcta implementación de la validación cruzada. En general, el AUC medio de los inner folds es consistentemente mayor o igual al de los outer folds a través de todos los folds, sugiriendo que el modelo generaliza adecuadamente sin sobreajustarse a los datos de entrenamiento. La similitud en los AUC medios entre los inner y outer folds en la mayoría de los casos refleja una buena estabilidad y consistencia en el rendimiento de las redes neuronales, con algunas diferencias observadas que son esperadas y muestran la variabilidad en el rendimiento del modelo dependiendo de la partición de los datos.

4.1.5.2 AUC Medio Folds Externos

TodasVariables

generate_outer_auc_violin_plot(ANN_TodasVariables_Metrics$performance_metrics)

Asociacion

generate_outer_auc_violin_plot(ANN_Asociacion_Metrics$performance_metrics)

StepAUC

generate_outer_auc_violin_plot(ANN_StepAUC_Metrics$performance_metrics)

Lasso

generate_outer_auc_violin_plot(ANN_Lasso_Metrics$performance_metrics)

Los gráficos de violín comparan la distribución del AUC por fold para los métodos de todas las variables, filtrado por análisis de asociación, stepAUC y embedded Lasso. En general, los métodos que usan todas las variables y el filtrado por análisis de asociación muestran distribuciones de AUC relativamente estrechas y consistentes en la mayoría de los folds, con valores centrados alrededor de 0.75 a 0.8, lo que indica un rendimiento estable. El método stepAUC muestra una variabilidad más amplia, especialmente en los folds 4 y 5, con AUC que varían desde aproximadamente 0.55 hasta 0.75, indicando una inestabilidad significativa en estos folds. Similarmente, el método embedded Lasso presenta una mayor dispersión en los folds 4 y 5, con rangos de AUC que van desde 0.6 a 0.75, sugiriendo una mayor inestabilidad y sensibilidad a las particiones de datos. En resumen, mientras que los métodos de todas las variables y el filtrado por análisis de asociación mantienen una consistencia aceptable, los métodos stepAUC y Lasso muestran mayor variabilidad y menor estabilidad en ciertos folds.

4.1.5.3 Comparativa ANN

fig <- generate_combined_auc_3d_plot(
  list(
    ANN_TodasVariables_Metrics$performance_metrics,
    ANN_Asociacion_Metrics$performance_metrics,
    ANN_StepAUC_Metrics$performance_metrics,
    ANN_Lasso_Metrics$performance_metrics
  ),
  c("Todas Variables", "Asociacion", "StepAUC", "Lasso")
)

fig

Los gráficos 3D muestran un análisis comparativo del AUC a través de los métodos de todas las variables, filtrado por asociación, stepAUC y Lasso, destacando sus desempeños a lo largo de diferentes semillas y folds. En general, el método Lasso presenta el AUC medio global más alto (0.746), seguido por el filtrado por asociación (0.74), todas las variables (0.733) y stepAUC (0.731). El método Lasso muestra una mayor consistencia en los valores de AUC, particularmente en las semillas 27 y 90, con menos variabilidad en comparación con los otros métodos. El filtrado por asociación también demuestra una consistencia razonable, aunque con ligeras fluctuaciones. Por otro lado, los métodos de todas las variables y stepAUC muestran una mayor variabilidad, especialmente en ciertos folds y semillas, indicando inestabilidades en su rendimiento. La visualización revela que, aunque todos los métodos tienen rendimientos comparables, Lasso y filtrado por asociación parecen ofrecer una mayor estabilidad y un rendimiento ligeramente superior.

TodasVariables

generate_auc_per_seed_and_fold(ANN_TodasVariables_Metrics$performance_metrics)

Asociacion

generate_auc_per_seed_and_fold(ANN_Asociacion_Metrics$performance_metrics)

StepAUC

generate_auc_per_seed_and_fold(ANN_StepAUC_Metrics$performance_metrics)

Lasso

generate_auc_per_seed_and_fold(ANN_Lasso_Metrics$performance_metrics)

4.2 Estabilidad para Hiperparámetros

generate_combined_hyperparameter_stability_analysis <- function(stepauc_metrics, lasso_metrics, asociacion_metrics, todasvariables_metrics, hyperparameters) {
  # Combinar los dataframes
  performance_metrics <- bind_rows(
    stepauc_metrics %>% mutate(method = "stepauc"),
    lasso_metrics %>% mutate(method = "lasso"),
    asociacion_metrics %>% mutate(method = "asociacion"),
    todasvariables_metrics %>% mutate(method = "todasvariables")
  )
  
  # Verificar el tamaño de la lista de hiperparámetros
  if (length(hyperparameters) == 2) {
    # Extraer los nombres de los hiperparámetros y sus valores posibles
    hyper_x <- names(hyperparameters)[1]
    hyper_y <- names(hyperparameters)[2]
    possible_values_x <- hyperparameters[[hyper_x]]
    possible_values_y <- hyperparameters[[hyper_y]]
    
    # Calcular la frecuencia de selección de cada combinación de hiperparámetros
    hyperparameter_grid <- expand.grid(possible_values_x, possible_values_y)
    colnames(hyperparameter_grid) <- c(hyper_x, hyper_y)
    
    frequency_data <- performance_metrics %>%
      count(!!sym(hyper_x), !!sym(hyper_y)) %>%
      right_join(hyperparameter_grid, by = c(hyper_x, hyper_y)) %>%
      replace_na(list(n = 0))
    
    # Crear el mapa de calor
    heatmap_plot <- plot_ly(
      data = frequency_data,
      x = ~get(hyper_x),
      y = ~get(hyper_y),
      z = ~n,
      type = 'heatmap',
      colorscale = 'Reds',
      colorbar = list(title = 'Frecuencia')
    ) %>%
      layout(
        title = paste('Mapa de Calor de Frecuencia de Hiperparámetros:', hyper_x, 'vs', hyper_y),
        xaxis = list(title = hyper_x),
        yaxis = list(title = hyper_y)
      )
    
    # Crear los gráficos de violín para mostrar la varianza de cada hiperparámetro
    violin_plot_x <- plot_ly(data = performance_metrics, y = ~get(hyper_x), type = 'violin', box = list(visible = TRUE), meanline = list(visible = TRUE), split = ~method) %>%
      layout(
        title = paste('Distribución del Hiperparámetro:', hyper_x),
        yaxis = list(title = hyper_x),
        xaxis = list(title = "Método")
      )
    
    violin_plot_y <- plot_ly(data = performance_metrics, y = ~get(hyper_y), type = 'violin', box = list(visible = TRUE), meanline = list(visible = TRUE), split = ~method) %>%
      layout(
        title = paste('Distribución del Hiperparámetro:', hyper_y),
        yaxis = list(title = hyper_y),
        xaxis = list(title = "Método")
      )
    
    return(list(heatmap_plot = heatmap_plot, violin_plot_x = violin_plot_x, violin_plot_y = violin_plot_y))
    
  } else if (length(hyperparameters) == 1) {
    # Extraer el nombre del único hiperparámetro y sus valores posibles
    hyper_x <- names(hyperparameters)[1]
    possible_values_x <- hyperparameters[[hyper_x]]
    
    # Crear los gráficos de violín para mostrar la varianza del único hiperparámetro
    violin_plot_x <- plot_ly(data = performance_metrics, y = ~get(hyper_x), type = 'violin', box = list(visible = TRUE), meanline = list(visible = TRUE), split = ~method) %>%
      layout(
        title = paste('Distribución del Hiperparámetro:', hyper_x),
        yaxis = list(title = hyper_x),
        xaxis = list(title = "Método")
      )
    
    return(list(violin_plot_x = violin_plot_x))
  } else {
    stop("La lista de hiperparámetros debe contener uno o dos elementos.")
  }
}

4.2.1 Decision Tree

hyperparameters <- list(
  Best_cp = c(0.001, 0.01, 0.1),
  Best_minsplit = c(2, 5, 10, 15, 20)
)

#RPart_StepAUC_Metrics$performance_metrics <- RPart_StepAUC_Metrics$performance_metrics %>% rename(Best_cp = Best_CP)
#RPart_StepAUC_Metrics$performance_metrics <- RPart_StepAUC_Metrics$performance_metrics %>% rename(Best_minsplit = Best_Minsplit)



result <- generate_combined_hyperparameter_stability_analysis(RPart_StepAUC_Metrics$performance_metrics, 
                                                              RPart_Lasso_Metrics$performance_metrics, 
                                                              RPart_Asociacion_Metrics$performance_metrics,
                                                              RPart_TodasVariables_Metrics$performance_metrics,
                                                              hyperparameters)

result$heatmap_plot

result$violin_plot_x

result$violin_plot_y

4.2.2 SVM

hyperparameters <- list(
  Best_Kernel = c(1, 2), # 1: Linear, 2: Poly
  Best_Cost = c(50, 100, 200)
)

result <- generate_combined_hyperparameter_stability_analysis(SVM_StepAUC_Metrics$performance_metrics, 
                                                              SVM_Lasso_Metrics$performance_metrics, 
                                                              SVM_Asociacion_Metrics$performance_metrics,
                                                              SVM_TodasVariables_Metrics$performance_metrics,
                                                              hyperparameters)

result$heatmap_plot

result$violin_plot_x

result$violin_plot_y

4.2.3 ANN

hyperparameters <- list(
  Best_Size = c(5, 10, 15, 20),
  Best_Decay = c(0.1, 0.2, 0.3)
)

result <- generate_combined_hyperparameter_stability_analysis(ANN_StepAUC_Metrics$performance_metrics, 
                                                              ANN_Lasso_Metrics$performance_metrics, 
                                                              ANN_Asociacion_Metrics$performance_metrics,
                                                              ANN_TodasVariables_Metrics$performance_metrics,
                                                              hyperparameters)

result$heatmap_plot

result$violin_plot_x

result$violin_plot_y

4.3 Estabilidad de Variables

# Función para contar la frecuencia de selección de variables agrupadas
contar_seleccion_variables_agrupadas <- function(df, variables_principales) {
  variable_counts <- setNames(rep(0, length(variables_principales)), variables_principales)
  
  for (variables in df$Best_Variables) {
    var_list <- unlist(strsplit(variables, ","))
    for (principal in variables_principales) {
      if (any(grepl(principal, var_list))) {
        variable_counts[principal] <- variable_counts[principal] + 1
      }
    }
  }
  
  result <- as.data.frame(variable_counts)
  result$Variable <- rownames(result)
  rownames(result) <- NULL
  colnames(result) <- c("Frecuencia", "Variable")
  result <- result[order(-result$Frecuencia), ]
  return(result)
}

# Función combinada
contar_seleccion_combined <- function(df1, df2, variables_principales) {
  # Contar frecuencias en el primer dataframe
  result_step_auc <- contar_seleccion_variables_agrupadas(df1, variables_principales)
  result_step_auc$Metodo <- "StepAUC"
  
  # Contar frecuencias en el segundo dataframe con variables agrupadas
  result_lasso <- contar_seleccion_variables_agrupadas(df2, variables_principales)
  result_lasso$Metodo <- "Lasso"
  
  # Combinar ambos resultados
  combined_result <- rbind(result_step_auc, result_lasso)
  return(combined_result)
}

# Función para combinar todos los dataframes y añadir la columna de Modelo
combinar_dataframes <- function(dataframes, modelos) {
  combined_df <- do.call(rbind, dataframes)
  
  # Asegurar que la longitud de los modelos coincida con la de los dataframes combinados
  combined_df$Modelo <- rep(modelos, each = nrow(dataframes[[1]]) / length(modelos))
  
  return(combined_df)
}

crear_visualizaciones <- function(combined_df) {
  # Métricas por Método de selección
  plot_metodo <- plot_ly(
    combined_df, 
    x = ~Variable, 
    y = ~Frecuencia, 
    color = ~Metodo, 
    type = 'bar'
  ) %>% layout(
    title = 'Métricas por Método de Selección',
    xaxis = list(title = 'Variable'),
    yaxis = list(title = 'Frecuencia')
  )
  
  # Métricas por Modelo
  plot_modelo <- plot_ly(
    combined_df, 
    x = ~Variable, 
    y = ~Frecuencia, 
    color = ~Modelo, 
    type = 'bar'
  ) %>% layout(
    title = 'Métricas por Modelo de Machine Learning',
    xaxis = list(title = 'Variable'),
    yaxis = list(title = 'Frecuencia')
  )
  
  return(list(plot_metodo = plot_metodo, plot_modelo = plot_modelo))
}

variablesFactor <- c("Edad", "REst", "RPro", "Estadio", "Grado", "Fenotipo", "NodAfec")


# ANN
ANN_Variables <- contar_seleccion_combined(ANN_StepAUC_Metrics$performance_metrics, 
                          ANN_Lasso_Metrics$performance_metrics, 
                          variablesFactor)
ANN_Variables$Modelo <- "ANN"

# SVM
SVM_Variables <- contar_seleccion_combined(SVM_StepAUC_Metrics$performance_metrics, 
                          SVM_Lasso_Metrics$performance_metrics, 
                          variablesFactor)
SVM_Variables$Modelo <- "SVM"

# DT
DT_Variables <- contar_seleccion_combined(RPart_StepAUC_Metrics$performance_metrics, 
                          RPart_Lasso_Metrics$performance_metrics, 
                          variablesFactor)
DT_Variables$Modelo <- "DT"

# ADA
ADA_Variables <- contar_seleccion_combined(AdaBoost_StepAUC_Metrics$performance_metrics, 
                          AdaBoost_Lasso_Metrics$performance_metrics, 
                          variablesFactor)
ADA_Variables$Modelo <- "ADA"

# NaiveBayes
Naive_Variables <- contar_seleccion_combined(NaiveBayes_StepAUC_Metrics$performance_metrics, 
                          NaiveBayes_Lasso_Metrics$performance_metrics, 
                          variablesFactor)
Naive_Variables$Modelo <- "Naive"


dataframes <- list(Naive_Variables, ADA_Variables, DT_Variables, ANN_Variables, SVM_Variables)
modelos <- c("NAIVE", "ADA", "DT", "ANN", "SVM")

combined_df <- combinar_dataframes(dataframes, modelos)
plots <- crear_visualizaciones(combined_df)
plots$plot_metodo

plots$plot_modelo

4.4 “Mejor” Modelo

5 Conclusión

González, Mario Pascual. 2024. «Modelado Predictivo: Cáncer de Mama». https://rpubs.com/mariopascuma/informecancermamamineria.