Una industria arrocera clasifica granos de arroz en clase 1 (evidentemente más grande en tamaño) y clase 2 (más pequeño). Esto lo hace por medios mecánicos y para controlar este proceso mide de forma digital (por imágen) distintos parámetros listados a continuación: -Area → Área -MajorAxisLength → Longitud del Eje Mayor -MinorAxisLength → Longitud del Eje Menor -Eccentricity → Excentricidad -ConvexArea → Área Convexa -Perimeter → Perímetro
Luego en función del tipo de arroz ya clasificado previamente requiere encontrar cual es el mejor método de clasificación para los lotes futuros. Para esto se realiza a continuación un análisis de datos.
Aquí cargaremos el excel dispuesto por la empresa “X” con los datos mencionados.
arroz<- read_excel("C:/Users/marco/Documents/Rscripts/arroz.xlsx")
## New names:
## • `` -> `...1`
Como primera medida se oberva la relación entre las variables por medio de gráficos de dispersión, en donde se mostrarán ambas clases en colores distintos ya que esta solo puede tomar valores binarios.
# Convertir la variable de respuesta en factor
arroz$Class <- as.factor(arroz$Class)
# Seleccionar las variables numéricas para el gráfico
variables_grafico <- arroz[, c("Area", "MajorAxisLength", "MinorAxisLength", "Eccentricity", "ConvexArea", "Perimeter")]
# Crear una paleta de colores para cada clase
colores <- c("blue", "red") # Usar 'blue' para la Clase 1 y 'red' para la Clase 2
# Crear el gráfico de pares con colores diferenciados para cada clase
pairs(variables_grafico,
col = colores[arroz$Class], # Asignar colores según la clase
pch = 19,
main = "Gráfico de dispersión de variables de arroz por clase")
Obs: podemos ver que los datos se encuentran uniformemente distribuidos por clases, sin embargo puntos de dos clases distintas comparten las mismas áreas, es decir no son heterogéneos en su totalidad. Una segunda observación es que la clase 1 (azul) es de mayor tamaño que la clase 2(roja), se evidencia viendo que el área es mayor en el grupo azul.
Como primer método de clasificación pondremos a prueba el KNN (K-Nearest Neightbors) en español Método de los Vecinos mas Cercanos, el cual en función de un punto ya definido como clase 1 o 2 respectivamente, busca los “K” puntos mas cercanos y los define como observaciones de la misma clase. En este caso se usan valores de K desde 2 (minimo posible) hasta 10, se calcular la correspondiente matriz de confusión, ecm, sensibilidad, especificidad, y precisión.
# Fijar la semilla para reproducibilidad
set.seed(123)
# Separar los datos en conjuntos de entrenamiento y prueba (70% entrenamiento, 30% prueba)
variables_grafico <- arroz[, c("Area", "MajorAxisLength", "MinorAxisLength", "Eccentricity", "ConvexArea", "Perimeter", "Class")]
trainIndex <- createDataPartition(variables_grafico$Class, p = 0.7, list = FALSE)
entrenamiento <- variables_grafico[trainIndex, ]
test <- variables_grafico[-trainIndex, ]
# Dividir en variables predictoras y de respuesta
train_x <- entrenamiento[, -which(names(entrenamiento) == "Class")]
train_y <- entrenamiento$Class
test_x <- test[, -which(names(test) == "Class")]
test_y <- test$Class
# Inicializar los valores de k y vectores para almacenar los resultados
k_values <- seq(2, 10, by = 1) # Valores de k desde 2 hasta 10, de 1 en 1
ecm_values <- numeric(length(k_values))
accuracy_values <- numeric(length(k_values))
sensitivity_values <- numeric(length(k_values))
specificity_values <- numeric(length(k_values))
cost_values <- numeric(length(k_values)) # Vector para almacenar los costos
# Inicializar una lista para almacenar las matrices de confusión
confusion_matrices <- list()
# Matriz de costo proporcionada
cost_matrix <- matrix(c(1000, -30, -50, 2000), nrow = 2, byrow = TRUE)
# Iterar sobre cada valor de k
for (i in seq_along(k_values)) {
k <- k_values[i]
# Aplicar el modelo KNN para el valor actual de k
predicciones <- knn(train = train_x, test = test_x, cl = train_y, k = k)
# Evaluar la clasificación mediante una matriz de confusión
t <- table(Predicted = predicciones, Actual = test_y)
# Almacenar la matriz de confusión en la lista
confusion_matrices[[i]] <- t
# Calcular el Error Cuadrático Medio (ECM) para el valor actual de k
ecm_values[i] <- mean((as.numeric(predicciones) - as.numeric(test_y))^2)
# Calcular la precisión (accuracy) para el valor actual de k
accuracy_values[i] <- sum(diag(t)) / sum(t)
# Calcular la sensibilidad (recall) para el valor actual de k
TP <- t[2, 2] # Verdaderos positivos
FN <- t[2, 1] # Falsos negativos
sensitivity_values[i] <- TP / (TP + FN)
# Calcular la especificidad para el valor actual de k
TN <- t[1, 1] # Verdaderos negativos
FP <- t[1, 2] # Falsos positivos
specificity_values[i] <- TN / (TN + FP)
# Calcular el costo total para el valor actual de k
total <- sum(t) # Total de muestras
P_VP <- TP / total # Probabilidad de Verdaderos Positivos
P_FP <- FP / total # Probabilidad de Falsos Positivos
P_FN <- FN / total # Probabilidad de Falsos Negativos
P_VN <- TN / total # Probabilidad de Verdaderos Negativos
total_cost <- (P_VP * cost_matrix[2, 2]) + (P_FP * cost_matrix[1, 2]) +
(P_FN * cost_matrix[2, 1]) + (P_VN * cost_matrix[1, 1])
cost_values[i] <- total_cost # Almacenar el costo esperado en el vector
}
# Graficar ECM, Precisión, Sensibilidad, Especificidad y Costo en función de k
par(mfrow = c(3, 2)) # Dividir el área de gráficos en tres filas y dos columnas
# Graficar ECM
plot(k_values, ecm_values, type = "o", col = "blue", pch = 16,
xlab = "Valor de k", ylab = "Error Cuadrático Medio (ECM)",
main = "ECM en función de k")
# Graficar Precisión
plot(k_values, accuracy_values, type = "o", col = "green", pch = 16,
xlab = "Valor de k", ylab = "Precisión (Accuracy)",
main = "Precisión en función de k")
# Graficar Sensibilidad
plot(k_values, sensitivity_values, type = "o", col = "red", pch = 16,
xlab = "Valor de k", ylab = "Sensibilidad",
main = "Sensibilidad en función de k")
# Graficar Especificidad
plot(k_values, specificity_values, type = "o", col = "purple", pch = 16,
xlab = "Valor de k", ylab = "Especificidad",
main = "Especificidad en función de k")
# Mostrar las matrices de confusión para cada k
for (i in seq_along(k_values)) {
cat("Matriz de Confusión para k =", k_values[i], "\n")
print(confusion_matrices[[i]])
cat("\n")
}
## Matriz de Confusión para k = 2
## Actual
## Predicted 0 1
## 0 2347 45
## 1 113 2950
##
## Matriz de Confusión para k = 3
## Actual
## Predicted 0 1
## 0 2345 29
## 1 115 2966
##
## Matriz de Confusión para k = 4
## Actual
## Predicted 0 1
## 0 2331 33
## 1 129 2962
##
## Matriz de Confusión para k = 5
## Actual
## Predicted 0 1
## 0 2330 27
## 1 130 2968
##
## Matriz de Confusión para k = 6
## Actual
## Predicted 0 1
## 0 2313 25
## 1 147 2970
##
## Matriz de Confusión para k = 7
## Actual
## Predicted 0 1
## 0 2305 27
## 1 155 2968
##
## Matriz de Confusión para k = 8
## Actual
## Predicted 0 1
## 0 2308 28
## 1 152 2967
##
## Matriz de Confusión para k = 9
## Actual
## Predicted 0 1
## 0 2307 27
## 1 153 2968
##
## Matriz de Confusión para k = 10
## Actual
## Predicted 0 1
## 0 2298 28
## 1 162 2967
Podemos observar en los gráficos Que si utilizamos k=3 optimicamos la
precisión, es decir la inversa del ECM La sensibilidad maxima se da en
K=2, La especificidad maxima se da en k=6. Por ende a simple vista no
podemos determinar el valor óptimo de K.
| ara elegir que método es mejor se empleó la matriz de costo: |
|---|
| | | Predicción 0 | Predicción 1 | |
| Real 0| 1000 | -30 | <- Costo de predicción correcta (1000) y falsa positiva (-30) | Real 1| -50 | 2000 | <- Costo de falsa negativa (-50) y predicción correcta (2000) ————————————– Explicación: - Cuando el modelo predice correctamente la clase 0 (verdaderos negativos), el costo es 1000, es decir ganamos $1000 - Cuando el modelo predice incorrectamente la clase 1 (falsos positivos), el costo es -30, es decir perdemos $30 - Cuando el modelo predice incorrectamente la clase 0 (falsos negativos), el costo es -50, es decir perdemos $50 - Cuando el modelo predice correctamente la clase 1 (verdaderos positivos), el costo es 2000. es decir ganamos $2000
A partir de la misma se calculo el costo promedio en función de k y se proyectó en el siguiente gráfico.
# Graficar Costo Total
plot(k_values, cost_values, type = "o", col = "orange", pch = 16,
xlab = "Valor de k", ylab = "Costo Total",
main = "Costo Total en función de k")
# Imprimir los costos en función del valor de k
cat("Costos en función del valor de k:\n")
## Costos en función del valor de k:
for (i in seq_along(k_values)) {
cat("k =", k_values[i], ": Costo Total =", cost_values[i], "\n")
}
## k = 2 : Costo Total = 1510.541
## k = 3 : Costo Total = 1516.11
## k = 4 : Costo Total = 1511.927
## k = 5 : Costo Total = 1513.967
## k = 6 : Costo Total = 1511.439
## k = 7 : Costo Total = 1509.155
## k = 8 : Costo Total = 1509.36
## k = 9 : Costo Total = 1509.54
## k = 10 : Costo Total = 1507.435
En función del gráfico de costos, se determina que si utilizamos un K=3 Maximizamos la ganancia obtenida (ya que el costo es positivo). Analiticamente si k = 3 el Costo Total Promedio = 1516.11
Como segundo método de clasificación se pondrá a prueba el método logístico, el cual modela la probabilidad de que una observación pertenezca a una clase específica, usualmente codificada como 0 o 1, por lo que no puede aplicarse el método de regresión.
# Fijar la semilla para reproducibilidad
set.seed(123)
# Separar los datos en conjuntos de entrenamiento y prueba (70% entrenamiento, 30% prueba)
variables_grafico <- arroz[, c("Area", "MajorAxisLength", "MinorAxisLength", "Eccentricity", "ConvexArea", "Perimeter", "Class")]
trainIndex <- createDataPartition(variables_grafico$Class, p = 0.7, list = FALSE)
entrenamiento <- variables_grafico[trainIndex, ]
test <- variables_grafico[-trainIndex, ]
# Dividir en variables predictoras y de respuesta
train_x <- entrenamiento[, -which(names(entrenamiento) == "Class")]
train_y <- entrenamiento$Class
test_x <- test[, -which(names(test) == "Class")]
test_y <- test$Class
# Modelo de Regresión Logística (glm) con el conjunto de entrenamiento
modelo.rl <- glm(Class ~ Area + MajorAxisLength + MinorAxisLength + Eccentricity + ConvexArea + Perimeter,
data = entrenamiento, family = binomial)
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
# Calcular las probabilidades estimadas en el conjunto de prueba
prob.est <- predict(modelo.rl, newdata = test_x, type = "response")
# Inicializar vectores para almacenar los resultados
umbral_values <- seq(0, 1, by = 0.1) # Valores del umbral de 0 a 1 con paso de 0.1
accuracy_values <- numeric(length(umbral_values))
sensitivity_values <- numeric(length(umbral_values))
specificity_values <- numeric(length(umbral_values))
ecm_values <- numeric(length(umbral_values)) # Vector para almacenar el ECM
cost_values <- numeric(length(umbral_values)) # Vector para almacenar los costos promedio por unidad
# Definir la matriz de costos
cost_matrix <- matrix(c(1000, -30, -50, 2000), nrow = 2, byrow = TRUE)
# Iterar sobre cada valor de umbral
for (i in seq_along(umbral_values)) {
umbral <- umbral_values[i]
# Clasificar las observaciones según el umbral
clas.pred <- ifelse(prob.est > umbral, 1, 0)
# Crear la matriz de confusión para evaluación (asegurando que tenga las 2 clases)
t <- table(Predicted = factor(clas.pred, levels = c(0, 1)), Actual = factor(test_y, levels = c(0, 1)))
# Imprimir la matriz de confusión
cat("\nMatriz de confusión para umbral =", umbral, ":\n")
print(t)
# Calcular las métricas de evaluación
TP <- ifelse(!is.na(t[2, 2]), t[2, 2], 0) # Verdaderos positivos
TN <- ifelse(!is.na(t[1, 1]), t[1, 1], 0) # Verdaderos negativos
FP <- ifelse(!is.na(t[1, 2]), t[1, 2], 0) # Falsos positivos
FN <- ifelse(!is.na(t[2, 1]), t[2, 1], 0) # Falsos negativos
# Calcular precisión, sensibilidad, especificidad y ECM
accuracy_values[i] <- sum(diag(t)) / sum(t)
sensitivity_values[i] <- TP / (TP + FN)
specificity_values[i] <- TN / (TN + FP)
# Calcular el Error Cuadrático Medio (ECM)
ecm_values[i] <- mean((clas.pred - as.numeric(test_y))^2) # ECM para clasificaciones binarias
# Calcular el costo promedio por unidad (por observación)
cost_values[i] <- (FP * cost_matrix[1, 2] + FN * cost_matrix[2, 1] +
TN * cost_matrix[1, 1] + TP * cost_matrix[2, 2]) / length(test_y)
}
##
## Matriz de confusión para umbral = 0 :
## Actual
## Predicted 0 1
## 0 0 0
## 1 2460 2995
##
## Matriz de confusión para umbral = 0.1 :
## Actual
## Predicted 0 1
## 0 2333 4
## 1 127 2991
##
## Matriz de confusión para umbral = 0.2 :
## Actual
## Predicted 0 1
## 0 2370 5
## 1 90 2990
##
## Matriz de confusión para umbral = 0.3 :
## Actual
## Predicted 0 1
## 0 2388 11
## 1 72 2984
##
## Matriz de confusión para umbral = 0.4 :
## Actual
## Predicted 0 1
## 0 2401 15
## 1 59 2980
##
## Matriz de confusión para umbral = 0.5 :
## Actual
## Predicted 0 1
## 0 2417 19
## 1 43 2976
##
## Matriz de confusión para umbral = 0.6 :
## Actual
## Predicted 0 1
## 0 2426 34
## 1 34 2961
##
## Matriz de confusión para umbral = 0.7 :
## Actual
## Predicted 0 1
## 0 2432 42
## 1 28 2953
##
## Matriz de confusión para umbral = 0.8 :
## Actual
## Predicted 0 1
## 0 2435 71
## 1 25 2924
##
## Matriz de confusión para umbral = 0.9 :
## Actual
## Predicted 0 1
## 0 2443 129
## 1 17 2866
##
## Matriz de confusión para umbral = 1 :
## Actual
## Predicted 0 1
## 0 2460 2995
## 1 0 0
# Graficar las métricas en función del umbral
par(mfrow = c(2, 2)) # Dividir el área de gráficos en cuatro paneles
# Graficar precisión
plot(umbral_values, accuracy_values, type = "o", col = "green", pch = 16,
xlab = "Umbral", ylab = "Precisión (Accuracy)",
main = "Precisión en función del umbral")
# Graficar sensibilidad
plot(umbral_values, sensitivity_values, type = "o", col = "red", pch = 16,
xlab = "Umbral", ylab = "Sensibilidad",
main = "Sensibilidad en función del umbral")
# Graficar especificidad
plot(umbral_values, specificity_values, type = "o", col = "purple", pch = 16,
xlab = "Umbral", ylab = "Especificidad",
main = "Especificidad en función del umbral")
# Graficar ECM
plot(umbral_values, ecm_values, type = "o", col = "blue", pch = 16,
xlab = "Umbral", ylab = "Error Cuadrático Medio (ECM)",
main = "ECM en función del umbral")
A simple vista y como es de esperarse, los mejores umbrales se ubican desde 0.1 hasta 0.9, donde los parámetros de ECM, especificidad y sensibilidad se mantienen casi contínuos, por lo que a continuación determinaremos que umbral es mejor en función de la matriz de costros.
Nuevamente al igual que en el método de KNN, determinaremos el costo unitario en función de cada umbral del método logístico.
# Graficar el costo promedio por unidad en función del umbral
plot(umbral_values, cost_values, type = "o", col = "orange", pch = 16,
xlab = "Umbral", ylab = "Costo Promedio por Unidad",
main = "Costo Promedio en función del umbral")
# Imprimir el costo promedio por unidad en función del umbral
cat("\nCosto Promedio por Unidad en función del Umbral:\n")
##
## Costo Promedio por Unidad en función del Umbral:
for (i in seq_along(umbral_values)) {
cat("Umbral:", umbral_values[i], "- Costo Promedio por Unidad:", cost_values[i], "\n")
}
## Umbral: 0 - Costo Promedio por Unidad: 1075.527
## Umbral: 0.1 - Costo Promedio por Unidad: 1523.104
## Umbral: 0.2 - Costo Promedio por Unidad: 1529.853
## Umbral: 0.3 - Costo Promedio por Unidad: 1531.085
## Umbral: 0.4 - Costo Promedio por Unidad: 1532.099
## Umbral: 0.5 - Costo Promedio por Unidad: 1533.69
## Umbral: 0.6 - Costo Promedio por Unidad: 1529.841
## Umbral: 0.7 - Costo Promedio por Unidad: 1528.018
## Umbral: 0.8 - Costo Promedio por Unidad: 1517.804
## Umbral: 0.9 - Costo Promedio por Unidad: 1497.76
## Umbral: 1 - Costo Promedio por Unidad: 434.4913
Obs: El Umbral: 0.5 cuyo Costo Promedio por Unidad: 1533.69 es el mayor, por ende es el que mayor “ganancias” genera.
El método de KNN con K=3 es el de su tipo que más “ganancias” o “ahorros” genera con un Costo Promedio = 1516.11 El método logístico con un umbral 0.5 arroja un Costo Promedio por Unidad: 1533.69, siendo el mejor en su tipo
Comparando estos dos resultados se concluye que el método logístico es mejor en este caso en particular