LABORATORIO 7. CLASIFICACIÓN DE IMÁGENES DE SATÉLITE (MÉTODOS SUPERVISADOS): CART
G174: FOTOINTERPRETACIÓN Y TELEDETECCIÓN
Domingo F. Rasilla
2026-04-26
INTRODUCCIÓN A CART
El término CART es el acrónimo de Classification and Regression Trees (Árboles de Clasificación y Regresión). Es uno de los algoritmos de aprendizaje automático más clásicos y potentes para la clasificación de imágenes de satélite.
Es una técnica de aprendizaje automático que construye una estructura gráfica similar a un árbol a partir de un conjunto de datos, mediante el desarrollo de múltiples reglas de decisión en su nodo raíz, sus ramas y nodos. Trabaja mediante un diagrama de flujo de decisiones, tomando como referencia los datos de entrenamiento y haciendo preguntas binarias (sí/no) sobre los valores de los píxeles. Por ejemplo, ¿es la reflectancia de un pixel en el infrarrojo cercano (NIR) es menor a 0.2? Si la respuesta es Sí, probablemente ese pixel corresponda a un “Suelo quemado” o a “Agua”; si la respuesta es No probablemente corresponda a “Vegetación”. Este procedimiento se repite banda por banda hasta llegar a una respuesta final (a qué clase corresponde ese pixel.
Este algorítmo puede aplicarse tanto para la clasificación de un dataset como para un procedimiento de regresión. Las diferencias son las siguientes:
Clasificación: se usa cuando la variable de salida es una categoría (ej: Bosque, Urbano, Quemado).
Regresión: se emplea cuando se quiere predecir un valor numérico continuo (ej: la concentración de partículas PM2.5 o la severidad del incendio en una escala de 0 a 100).
Para la clasificación de usos de suelo tiene varias ventajas:
Interpretabilidad: a diferencia de otros procedimientos, no es una “caja negra”; en CART se obtienen las reglas utilizadas por el algoritmo para separar unas categorías respecto a otras.
No requiere normalidad: no importa si los datos no siguen una campana de Gauss (a diferencia de otros métodos, como el de “Máxima Verosimilitud”).
Maneja **datos complejos*: funciona bien aunque existan bandas con escalas muy diferentes.
PROCEDIMIENTOS INICIALES
Se comienza la actividad activando, si fuera necesario, los siguientes paquetes
library(terra) # Importación y análisis de imágenes## terra 1.8.93library(RStoolbox) # Visualización de imágenes## This is version 1.0.2.2 of RStoolboxlibrary(ggplot2) # Visualización de imágenes
library(sf) # Importación de datos formato vectorial## Linking to GEOS 3.11.2, GDAL 3.8.2, PROJ 9.3.1; sf_use_s2() is TRUElibrary(tidyverse)## Warning: package 'tidyverse' was built under R version 4.3.2## Warning: package 'readr' was built under R version 4.3.2## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr 1.2.0 ✔ readr 2.1.5
## ✔ forcats 1.0.1 ✔ stringr 1.6.0
## ✔ lubridate 1.9.5 ✔ tibble 3.3.1
## ✔ purrr 1.2.1 ✔ tidyr 1.3.2## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ tidyr::extract() masks terra::extract()
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag() masks stats::lag()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errorsA continuación se importa el fichero que contiene los valores de
reflectividad según coberturas de suelo con la función
read.csv2, a la que se añade el argumento
rownames = 1 para convertir la primera columna en nombres
de fila.
df_muestra <- read.csv2("D:/G174_2026/LABORATORIO_7_Clasificacion_imagenes_PIXEL/datos/supervisada/df_muestra.csv", row.names = 1)Puede verficarse a continuación
df_muestraTambién se importará la escena Landsat.
imagen <- rast("D:/G174_2026/LABORATORIO_7_Clasificacion_imagenes_PIXEL/datos/supervisada/imagen.tif")
imagen ## class : SpatRaster
## size : 1245, 1497, 7 (nrow, ncol, nlyr)
## resolution : 30, 30 (x, y)
## extent : 594090, 639000, 4190190, 4227540 (xmin, xmax, ymin, ymax)
## coord. ref. : WGS 84 / UTM zone 10N (EPSG:32610)
## source : imagen.tif
## names : coast~rosol, blue, green, red, NIR, SWIR1, ...
## min values : 0.09641791, 0.0748399, 0.04259216, 0.02084067, 0.0008457669, -0.007872183, ...
## max values : 0.73462820, 0.7177562, 0.69246972, 0.78617686, 1.0124315023, 1.043204546, ...ENTRENAMIENTO DEL MODELO USANDO EL FICHERO DE ENTRENAMIENTO
Para implementar la clasificación supervisada de datos ráster aplicando CART, se requiere el paquete rpart, mientras que para mostrar gráficamente los modelos se requiere rpart.plot
library(rpart)
library(rpart.plot) Posteriormente, se aplicará el algoritmo CART al conjunto de datos de
entrenamiento usando la función rpart(). Esta función
presenta los siguientes argumentos:
La variable
categoríaes la variable dependiente. Las restantes, que contienen los valores de reflectancia en las bandas respectivas, conforman el conjunto de variables independientes (o predictoras), y se mencionan en el código como ~. para resumir.Como CART se usa con fines de clasificación, se asigna el argumento
method = "class".Los parámetros de control para el algoritmo CART se mencionan en el argumento
control. El argumentominspliten la funciónrpart.control()indica el número mínimo de observaciones que deben estar presentes en un nodo para intentar una división. El argumentocprepresenta el valor del parámetro de complejidad que determina la división: a mayor valor un árbol más pequeño y viceversa.
modelo_cart <- rpart(as.factor(categoria)~.,
data = df_muestra,
method = 'class',
control = rpart.control(cp = 0.0, minsplit = 2))Una de las principales razones para elegir el modelo de CART es
visualizar el modelo de clasificación. A diferencia de otros algoritmos,
CART ofrece una forma muy sencilla de inspeccionar y trazar la
estructura del modelo. Por ejemplo, la función plotcp()
dibuja un diagrama del que se obtiene información sobre el número óptimo
de divisiones y el parámetro de complejidad, gracias a que muestra
gráficamente un diagrama de error
plotcp(modelo_cart)
En este diagrama se observa que el error relativo mínimo se alcanza
con un valor de cp igual a 0 y con un árbol con 5
subdivisiones. Esto permitirá “podar” (es decir, parar el proceso de
clasificación) cuando se alcancen estos parámetros óptimos. Ambos
parámetros se incluyen ahora y se repite el procedimiento.
cart_mod <- rpart(as.factor(categoria)~.,
data = df_muestra,
method = 'class',
control = rpart.control(cp = 0.0, minsplit = 5))A continuación se muestra los resultados (numéricos) del modelo.
print(modelo_cart)## n= 200
##
## node), split, n, loss, yval, (yprob)
## * denotes terminal node
##
## 1) root 200 122 agua (0.28 0.39 0.1 0.095 0.12)
## 2) NIR< 0.08322171 78 0 agua (0 1 0 0 0) *
## 3) NIR>=0.08322171 122 65 abierto (0.47 0 0.17 0.16 0.2)
## 6) SWIR1>=0.2775647 63 6 abierto (0.9 0 0 0 0.095)
## 12) SWIR2< 0.2706468 59 2 abierto (0.97 0 0 0 0.034)
## 24) blue< 0.1694362 58 1 abierto (0.98 0 0 0 0.017)
## 48) NIR>=0.2590879 57 0 abierto (1 0 0 0 0) *
## 49) NIR< 0.2590879 1 0 urbano (0 0 0 0 1) *
## 25) blue>=0.1694362 1 0 urbano (0 0 0 0 1) *
## 13) SWIR2>=0.2706468 4 0 urbano (0 0 0 0 1) *
## 7) SWIR1< 0.2775647 59 38 barbecho (0 0 0.36 0.32 0.32)
## 14) blue< 0.1334584 40 19 barbecho (0 0 0.52 0.47 0)
## 28) blue>=0.09720947 21 0 barbecho (0 0 1 0 0) *
## 29) blue< 0.09720947 19 0 cultivos (0 0 0 1 0) *
## 15) blue>=0.1334584 19 0 urbano (0 0 0 0 1) *La estructura del árbol de decisión puede visualizarse gráficamente
mediante la función prp() del paquete
rpart.plot. El argumento varlen = 0
incluye los nombres completos de las variables en el gráfico, mientras
que el argumento cex determina el tamaño del texto.
prp(modelo_cart,
varlen =0,
yesno =TRUE,
cex =0.8)
El árbol de decisión muestra claramente cómo los valores de reflectancia de las diferentes bandas en la imagen ráster determinan cada una de las clases.
Se puede consultar ?rpart.control para establecer diferentes
parámetros para construir el modelo. También se puede utilizar
plotcp) para conocer el conocido como
coste-complejidad (argumento cp en rpart).
CLASIFICACIÓN DE LA IMAGEN A PARTIR DEL MODELO ENTRENADO
Una vez entrenado el algoritmo y obtenido el modelo de clasificación
(modelo_cart), éste se utiliza para clasificar todas las
celdas de la imagen original. Importante : Los nombres de las capas en
el SpatRaster deberían coincidir exactamente con los que se usaron para
entrenar el modelo. Esto ocurrirá si se usó el mismo objeto SpatRaster
(mediante extract) para obtener los valores que ajustan al modelo. Si
no, tienes que especificar los nombres que coinciden.
raster_cart <- predict(imagen,
modelo_cart,
na.rm = TRUE)## |---------|---------|---------|---------|========================================= raster_cart## class : SpatRaster
## size : 1245, 1497, 5 (nrow, ncol, nlyr)
## resolution : 30, 30 (x, y)
## extent : 594090, 639000, 4190190, 4227540 (xmin, xmax, ymin, ymax)
## coord. ref. : WGS 84 / UTM zone 10N (EPSG:32610)
## source : spat_67187d3825fb_26392_o8up7yuDdJ6XHBX.tif
## names : abierto, agua, barbecho, cultivos, urbano
## min values : 0, 0, 0, 0, 0
## max values : 1, 1, 1, 1, 1plot(raster_cart)
El resultado es un objeto ráster con 5 capas, cada una de las cuales representa la probabilidad de que ese pixel pertenezca a una clase de uso de suelo concreta. Se puede simplificar este resultado creando un nuevo ráster que muestra, para cada pixel, el uso de suelo más probable.
raster_cart_final <- which.max(raster_cart)
raster_cart_final## class : SpatRaster
## size : 1245, 1497, 1 (nrow, ncol, nlyr)
## resolution : 30, 30 (x, y)
## extent : 594090, 639000, 4190190, 4227540 (xmin, xmax, ymin, ymax)
## coord. ref. : WGS 84 / UTM zone 10N (EPSG:32610)
## source(s) : memory
## name : which.max
## min value : 1
## max value : 5EVALUACIÓN DE LA CLASIFICACIÓN
El último aspecto a revisar es la evaluación de la precisión del modelo para deducir si el mapa temático refleja la distribución real de las categorías o usos de suelo. Téngase en cuenta que el proceso de clasificación se ha realizado a partir de una pequeña muestra, que luego se extiende a toda la imagen.
La evaluación de la clasificación supervisada es un paso crítico para garantizar que el mapa generado sea una representación fiel de la realidad. Este proceso se divide principalmente en dos enfoques: el análisis cualitativo (inspección visual) y el análisis cuantitativo (estadístico).
La Matriz de Confusión (o de Errores)
Es la herramienta fundamental para la evaluación cuantitativa. Compara la clasificación del modelo frente a datos de referencia (verdad terreno) en una tabla cruzada.
Precisión Global (Overall Accuracy): Es el porcentaje total de píxeles correctamente clasificados (suma de la diagonal dividida por el total).
Precisión del Productor (Omisión): Indica qué tan bien se clasificó un área real específica. Un valor bajo significa que se omitieron áreas que pertenecían a esa clase.
Precisión del Usuario (Comisión): Indica la fiabilidad del mapa para el usuario final. Un valor bajo significa que el mapa incluye erróneamente áreas que no pertenecen a esa clase.
Estadísticas de Concordancia y Ajuste
Más allá del conteo directo, se utilizan índices que ajustan la precisión por el azar o el desbalance de clases:
Coeficiente Kappa: Mide la concordancia entre la clasificación y los valores reales, descontando la probabilidad de aciertos por puro azar. Un valor de 1 es concordancia perfecta; cercano a 0 indica que la clasificación no es mejor que una asignación aleatoria.
F1-Score: Muy utilizado cuando hay clases desbalanceadas (ej. pocas zonas urbanas frente a mucho bosque). Es la media armónica entre la precisión y el recall (sensibilidad), buscando un equilibrio entre ambas.
Técnicas de Validación y Muestreo
La fiabilidad de la evaluación depende de cómo se eligen los datos para probar el modelo:
Conjunto de Validación Independiente: Es vital evaluar el error con datos que el clasificador no vio durante el entrenamiento para evitar el sobreajuste (overfitting).
Validación Cruzada (Cross-Validation): Divide los datos en varios grupos para entrenar y validar repetidamente, asegurando que los resultados no dependan de una partición de datos específica.
Muestreo Estratégico: Se recomienda un muestreo controlado con suficientes puntos por categoría para evitar confusiones entre clases con respuestas espectrales similares.
Evaluación Visual y Espacial
Análisis Cualitativo: Consiste en examinar visualmente el mapa para detectar patrones de error obvios, como el “efecto sal y pimienta” (píxeles aislados mal clasificados).
Mapas de Error: Generar un mapa que muestre específicamente dónde falló la clasificación permite identificar si los errores son geográficos (ej. fallos sistemáticos en laderas de sombra)
Para evaluar cualquier modelo, puedes usar la validación cruzada k-fold. En esta técnica, los datos utilizados para ajustar el modelo se dividen en k grupos (normalmente 5 grupos). A su vez, uno de los grupos se utilizará para pruebas de modelos, mientras que el resto de los datos se destinará al entrenamiento (ajuste) del modelo.
set.seed(99) # number of folds
k <- 5
j <- sample(rep(1:k, each = round(nrow(df_muestra))/k))
table(j)## j
## 1 2 3 4 5
## 40 40 40 40 40Ahora entrenamos y probamos el modelo cinco veces, cada vez calculando las predicciones y almacenando eso con los valores reales en una lista. Más adelante usamos la lista para calcular la obtención final.
x <- list()
for (k in 1:5) {
train <- df_muestra[j!= k, ]
test <- df_muestra[j == k, ]
cart <- rpart(as.factor(categoria)~., data=train, method = 'class',
minsplit = 5)
pclass <- predict(cart, test, na.rm = TRUE)
# assign class to maximum probablity
pc <- apply(pclass, 1, which.max)
# use labels instead of numbers
pc <- colnames(pclass)[pc]
# create a data.frame using the reference and prediction
x[[k]] <- cbind(test$categoria, pc)
}Ahora combina los cinco elementos de la lista en un solo data.frame, usando do.call y calcula una matriz de confusión.
y <- do.call(rbind, x)
y <- data.frame(y)
colnames(y) <- c('observado', 'predicho')
cart_matriz_confusion <- table(y) # confusion matrix
print(cart_matriz_confusion)## predicho
## observado abierto agua barbecho cultivos urbano
## abierto 52 0 1 0 4
## agua 0 78 0 0 0
## barbecho 0 0 21 0 0
## cultivos 0 0 0 19 0
## urbano 4 0 1 0 20- Precisión general:
Primero se calcula el número total de casos/muestras
n <- sum(cart_matriz_confusion)
n## [1] 200Número de casos correctamente clasificados por clase
diag <- diag(cart_matriz_confusion)Precisión general
OA <- sum(diag) / n
OA## [1] 0.95- Kappa: Casos observados (verdaderos) por clase
rowsums <- apply(cart_matriz_confusion, 1, sum)
p <- rowsums / nCasos predichos por clase
colsums <- apply(cart_matriz_confusion, 2, sum)
q <- colsums / n
expAccuracy <- sum(p*q)
kappa <- (OA - expAccuracy) / (1 - expAccuracy)
kappa## [1] 0.931694Precisión tanto de productores como de **usuarios*
Producer Accuracy (precisión del productor): Mide el error de omisión. Es la probabilidad de que un valor real (un píxel en el terreno) haya sido clasificado correctamente. Si de 100 píxeles que son realmente “Bosque”, el mapa identifica 90, la exactitud del productor es del 90%. El 10% restante son píxeles de bosque que el mapa “omitió”.
PA <- diag / colsums- Users Accuracy (precisión del usuario): mide el error de comisión. Es la probabilidad de que un píxel clasificado en el mapa represente realmente esa categoría en el terreno. Si el mapa dice que 100 píxeles son “Urbano”, en el terreno solo son 80, la exactitud del usuario es del 80%. El 20% restante son píxeles que el mapa “cometió” el error de incluir en esa clase sin serlo.
UA <- diag / rowsums
outAcc <- data.frame(producerAccuracy = PA, userAccuracy = UA)
outAcc## producerAccuracy userAccuracy
## abierto 0.9285714 0.9122807
## agua 1.0000000 1.0000000
## barbecho 0.9130435 1.0000000
## cultivos 1.0000000 1.0000000
## urbano 0.8333333 0.8000000Los resultados son excelentes para un modelo basado en CART. Interpretación:
Agua y Cultivos (1.00): separación espectral perfecta. Esto es lógico porque el agua absorbe casi todo el infrarrojo (NIR) y los cultivos tienen un pico de reflectancia muy distintivo. No hay errores ni de omisión ni de comisión.
Urbano (0.83 / 0.80): es la clase más “débil” (aunque sigue siendo un valor muy bueno). Exactitud del Productor (0.83): El modelo “perdió” un 17% de las zonas urbanas reales (las clasificó como otra cosa, probablemente “abierto” o “barbecho”). Exactitud del Usuario (0.80): Cuando el mapa dice “Urbano”, hay un 20% de probabilidad de que sea otra clase (error de comisión). Esto suele pasar porque el asfalto o el hormigón pueden confundirse espectralmente con suelos desnudos o rocosos.
Abierto y Barbecho: valores por encima del 90% indican una clasificación “robusta” para tipos de cobertura que suelen ser muy parecidos entre sí.
Por último, el fichero ráster conteniendo la clasificación se salvaría como fichero tiff.
writeRaster(raster_cart_final, "D:/G174_2026/LABORATORIO_7_Clasificacion_imagenes_PIXEL/datos/supervisada/raster_CART_final.tif)VISUALIZACIÓN DE LA CLASIFICACIÓN.
Una vez obtenida y validada la clasificación, se interpretar desde el punto de vista de las coberturas del suelo (clases de información). Para ello, los píxeles clasificados (con una transparencia del 0.4, my alta) se superponen a la imagen en color verdadero
plotRGB(imagen, r=3, g=2, b=1, axes=TRUE, stretch="lin")
plot(raster_cart_final, col = rainbow(5),
main = 'Visualización preliminar clasificación con CART', add=TRUE, alpha = 0.4, legend = TRUE)
Una vez realizada una interpretación preliminar, hay que identificar
los grupos o clusters (categorías espectrales) en forma de coberturas
reales del territorio (clases de información). Dado que el objeto ráster
raster_kmeans es un ráster cuantitativo, se debe
transformar los números en variables categóricas, convirtiéndolos en
factor (variable categórica), Una vez conozcamos la relación con las
clases de información, se dará nombre a estas categorías
raster_cart_final <- as.factor(raster_cart_final)Si ya supiéramos a qué clases de información corresponden las clases espectrales, se deberían asignar etiquetas, por ejemplo
levels(raster_cart_final) <- data.frame(ID = 1:5,
clase = c("abierto","agua","barbecho", "cultivos", "urbano"))Una vez creado este nuevo ráster con los resultados de la clasificación, se puede representar gráficamente, asignando etiquetas de color y clase según la interpretación visual (modo ‘classes’)
colores_clasificacion_cart <- c("abierto" = "lightgreen", "agua" = "blue", "barbecho" = "burlywood", "cultivos" = "darkgreen", "urbano" = "red")
levels(raster_cart_final)## [[1]]
## ID clase
## 1 1 abierto
## 2 2 agua
## 3 3 barbecho
## 4 4 cultivos
## 5 5 urbanoFinalmente, representaríamos gráficamente el mapa con ggplot2
ggR(raster_cart_final, geom_raster = TRUE, stretch = "lin", forceCat = TRUE) +
ggtitle("Clasificación no supervisada. Método CART") + # Título
labs(x ="Longitud (m)", y ="Latitud (m)") + # Etiquetas de los ejes
theme(plot.title = element_text(hjust =0.5, # Alineación del título
size =20), # Tamaño del título
axis.title = element_text(size =10), # Tamaño de las etiquetas de los ejes
legend.key.size =unit(1,"cm"), # Tamaño de la leyenda
legend.title = element_text(size = 20), # Tamaño del título de la leyenda
legend.text = element_text (size = 10)) + # Tamaño del texto de la leyenda
# Dibuja los colores según un gradiente
scale_fill_manual(values = colores_clasificacion_cart,
name = "Clase") +
theme_minimal()