La clasificación supervisada es una técnica en la que el usuario proporciona información previa sobre las clases en forma de una serie de lugares previamente clasificados. Estos lugares o “sitios de entrenamiento” consisten en áreas homogéneas que pertenecen a una clase particular. El patrón de respuesta espectral de estos sitios se utiliza para clasificar toda una imagen, ya que los píxeles de la imagen con un patrón de respuesta espectral similar al de los sitios de entrenamiento de una clase particular son agrupados en esa misma clase. El uso de varios sitios de entrenamiento por cada clase se explica por la variabilidad en el patrón de respuesta espectral de una clase en particular.
La clasificación supervisada de imágenes se puede realizar mediante varios métodos. Como métodos tradicionales se podrían citar el clasificador de máxima verosimilitud (“maximum likehood”), el clasificador de distancia mínima a la media (“minimum distance to mean”)y el clasificador de paralelepípedos. Actualmente, cada vez hay un uso más generalizado de algoritmos de aprendizaje automático, que son entrenados inicialmente con la información proporcionada por los sitios de entrenamiento. Este modelo entrenado se valida luego con conjuntos de datos independientes, para evaluar su precisión; si esta es aceptable, dicho modelo se usará para la clasificación del resto de la imagen.
En el mundo actual los algoritmos de aprendizaje automático tienen numerosas aplicaciones en una gran cantidad de campos que van desde los negocios, el entretenimiento, la ingeniería de software, la medicina, las redes de comunicación, la seguridad informática, el diagnóstico de fallas, el reconocimiento de voz, la minería web, etc. Por lo tanto, las aplicaciones de estos algoritmos en el análisis de datos de detección de remate son obvias. Según Svensson y Soderberg (2008), ‘el aprendizaje automático se ocupa del diseño y desarrollo de algoritmos y técnicas que permiten a las computadoras ’aprender’.
El enfoque principal de la investigación mediante procedimientos de aprendizaje automático es extraer información de los datos automáticamente, mediante métodos computacionales y estadísticos. Por lo tanto, está estrechamente relacionado con la minería de datos y la estadística. El aprendizaje automático es una rama de la inteligencia artificial. Por lo tanto, los algoritmos de aprendizaje automático no requieren un conjunto de ecuaciones predeterminaas, sino que construye algoritmos automáticamente a partir de los datos de entrada, proporcionando una salida. La salida se actualiza cuando el usuario agrega un nuevo conjunto de datos de entrada.
Se pueden usar muchos algoritmos de aprendizaje automático para la clasificación supervisada de datos ráster. En este capítulo se analizan algunos de los algoritmos enumerados a continuación y su implementación en el software R.
Para la implementación de los algoritmos mencionados anteriormente en el software R, se utilizan los datos ráster de entrada importados anteriormente. Los datos de muestra para entrenamiento y validación están presentes como archivos de texto en formato “*.csv”.
En este capítulo se mostrarán algunas técnicas de clasificación supervisada. Por lo tanto, partimos de un conocimiento previo sobre el tipo de usos del suelo, obtenido a través de trabajos de campo, de bases de datos de referencia (por ejemplo, Corine Land Cover) o a través de la interpretación de imágenes de alta resolución (como las disponibles en los mapas de Google). En todos los casos, el procedimiento de clasificación requiere los siguientes pasos:
Es fundamental disponer de un fichero con los sitios de entrenamiento, que incluyen una muestra de píxeles clasificados en seis clases o categorías:
Antes de revisar los distintos procedimientos de análisis, son necesarias varias acciones.
setwd("D:/G174_Teledeteccion_2022/OCW/Lab_10_Clasificacion_imagenes/Clasificacion_supervisada")
library(rgdal)
library(raster) # Importación y análisis de imágenes
library(RStoolbox) # Visualización de imágenes
library(ggplot2) # Visualización de imágenes
library(sp)
# library(sf) # Importación de datos formato vectorial
RasterBrick
que incluye 9 bandas.<- brick("./datos/input_img.tif")
img img
Para un mejor desarrollo de la actividad, las diferentes bandas del
RasterBrick
se identificarán con su nombre.
names(img) <- c("Coastal_aerosol", "Blue", "Green", "Red", "NIR", "SWIR1", "SWIR2", "cirrus", "TIRS1", "TIRS2")
img
A continuación se representa gráficamente la escena en Falso color
ggRGB(img, r =4, g=3, b=2, stretch ="lin") +
ggtitle("Imagen original") + # 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 =15)) # Tamaño de las etiquetas de los ejes
Este tipo de método, conocido como CART por sus siglas en inglés (“Classification and Regression Tree”) 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, tanto para la clasificación como para la regresión.
CART es un modelo predictivo que los usuarios pueden interpretar fácilmente y puede manejar tanto datos no paramétricos y no lineales de manera eficiente.
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)
Esta técnica requiere utilizar un fichero con datos de entrenamiento:
<- read.csv("./datos/muestra_datos.csv",
muestra_datos header =TRUE)
head(muestra_datos)
La primera columna o “Clase” es la variable dependiente; las restantes, que contienen los valores de reflectancia en las bandas respectivas, son las variables independientes,
Posteriormente, es necesario subdividir los datos de la muestra en un
conjunto de “entrenamiento” y otro de
“verificación”. Para ello se utiliza la función
sample()
que crea otro objeto, denominado denominado
muestra
, del mismo tamaño que muestra_datos
(nrow(muestra_datos) que contiene unos y doses aleatoriamente
<-sample(2,
muestra nrow(muestra_datos),
replace =TRUE,
prob = c(0.7,0.3))
<- muestra_datos[muestra==1, ]
entrenamiento <- muestra_datos[muestra==2, ]
verificacion summary(entrenamiento)
Posteriormente, se aplicará el algoritmo CART al conjunto de datos de
entrenamiento usando la función rpart()
. Esta función
presenta los siguientes argumentos:
Como se ha señalado, la variable “Clase” es la variable
dependiente y 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 los estamos usando 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 argumento minsplit
en la función
rpart.control()
indica el número mínimo de observaciones
que deben estar presentes en un nodo para intentar una
división.
El argumento cp
representa el valor del parámetro de
complejidad que determina la división: a mayor valor un árbol más
pequeño y viceversa.
<- rpart(Clase ~.,
cart data = entrenamiento,
method = "class",
control = rpart.control(cp = 0.0, minsplit = 2))
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(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 9
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.
<- rpart(Clase ~.,
cart_mod data = entrenamiento,
method = "class",
control = rpart.control(cp = 0.0, minsplit = 9))
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(cart_mod,
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.
Une vez entrenado el modelo, el siguiente paso es comparar la clase predicha con la clase real. La bondad del modelo se puede determinar mediante el cálculo de la matriz de confusión y de la precisión. Para ello, se aplica CART para predecir las clases en el conjunto de entrenamiento.
<- predict(cart_mod,
cart_entrenamiento
entrenamiento, type = "class")
A continuación se calcula la matriz de confusión
<- table(actual = entrenamiento[ ,1], predictions = cart_entrenamiento)
tr tr
La precisión es el resultado de
sum(diag(tr))/sum(tr)
Dado que los resultados son muy buenos, el siguiente paso consiste en repetir el mismo procedimiento, pero ahora sobre el conjunto de datos de verificación.
<- predict(cart_mod,
cart_verificacion type = "class") verificacion,
A continuación se vuelve a calcular la matriz de confusión, pero sobre el
<- table(actual = verificacion[ ,1], predictions = cart_verificacion)
ts ts
La precisión es el resultado de
sum(diag(ts))/sum(ts)
De acuerdo con estos resultados, la precisión del modelo CART desarrollado es buena, por lo que podemos usar el modelo obtenido para la clasificación de toda la imagen original.
<- predict(img,
raster_cart
cart_mod, type= "class")
raster_cart
A continuación se dibujará la imagen con los píxeles clasificados según los resultados de CART.
ggR(raster_cart, geom_raster = TRUE) + # PLot with ggplot2
ggtitle("Clasificación con CART") + # Título
labs(x="Longitud (m)", y="Latitud (m)") + # Etiquetas de los ejes X e Y
theme(plot.title=element_text(hjust =0.5, # Título: alineación en el centro
size =30), # Título: tamaño
axis.title=element_text(size=20), # Etiquetas de los ejes: tamaño
legend.key.size=unit(1, "cm"), # Leyenda: tamaño.
legend.title=element_text(size =20), # Leyenda: tamaño del título.
legend.text =element_text(size =15)) + # Leyenda: tamaño del texto.
scale_fill_manual(values= c("grey", "burlywood", "wheat", "black", "darkgreen", "skyblue"),
name="Class")
Para no favorecer la acumulación de objetos en el Global environment, se eliminan a continuación.
rm(list=ls())
Random forest es otra técnica de aprendizaje automático que construye múltiples árboles de decisión aleatorios no correlacionados a partir del subconjunto de datos de entrenamiento y agrupa sus resultados para proporcionar una decisión final o resultado, evitando el problema del sobreajuste. Este algorítmo también se utiliza para la clasificación y la regresión; el resultado final en caso de clasificación utiliza la predicción de la moda, mientras que utiliza la media en el caso de la regresión.
En R se puede implementar este procedimiento de dos maneras:
El primer procedimiento es similar al visto en el caso de CART, y consiste en la conversión de la imagen en una matriz y la posterior aplicación de Random Forest usando el paquete randomForest.
La segunda es mediante una función integrada en el paquete RStoolbox directamente a la imagen original.
Ambos métodos serán analizados en esta sección.
En primer lugar, deben activarse los paquetes necesarios para realizar este ejercicio.
library(randomForest) # Para ejecutar el primer procedimiento
library(caret) # Para la creación de la matriz de confusión
Usaremos el mismo ráster multibanda utilizado en el caso de la clasificación con CART
<- brick("./datos/input_img.tif")
img names(img) <- c("Coastal_aerosol", "Blue", "Green", "Red", "NIR", "SWIR1", "SWIR2", "cirrus", "TIRS1", "TIRS2")
img
Para aplicar Random Forest usando el paquete RStoolbox, se debe importar el fichero que contiene los sitios de entrenamiento. Este fichero será importado en formado vectorial, por lo que deberá estar proyectado en el mismo Sistema de referencia de coordenadas (CRS) de la imagen ráster.
<- readOGR(dsn = "./datos", "training_points")
datos_entrenamiento datos_entrenamiento
El objeto datos_entrenamiento
tiene varios slots
slotNames(datos_entrenamiento)
El fichero contiene 400 puntos de entrenamiento, correspondientes a seis clases. Para determinar si ambos ficheros comparten el mismo Sisdtema de coordenadas de referencia,
compareCRS(img,
datos_entrenamiento)
Dado que poseen sistemas de proyección diferentes, se reproyectará el fichero en formato vectorial al mismo CRS del fichero raster.
<- spTransform(datos_entrenamiento,
datos_entrenamiento CRS(projection(img)))
compareCRS(img, datos_entrenamiento)
Una vez ambos tipos de ficheros comparten el mismo sistema de
proyección, se puede iniciar el procedimiento de clasificación. En la
función superClass()
debe asignarse el argumento
model = "rf"
. Además, debe utilizarse un 70 % de los datos
para el entrenamiento y el 30 % de los datos restantes para la
verificación asignando el argumento trainPartition = 0.7
.
La descripción de los restantes parámetros se proporciona en la sección
dedicada al procedimiento de máxima verosimilitud.
<- superClass(img,
img_rf trainData = datos_entrenamiento,
responseCol = "Class",
nSamples = 50,
model = "rf",
mode = "classification",
trainPartition = 0.7,
tuneLength = 3,
kfold =10)
img_rf
La salida del código R muestra que el modelo randomForest desarrollado tiene muy buena precisión. Podemos ver los atributos de la imagen raster ya clasificada y las etiquetas correspondintes.
$map
img_rf
$classMapping img_rf
Como se ha señalado en líneas anteriores, para aplicar el paquete randomForest a la clasificación de una imagen, se deben convertir los valores de cada banda de esta última en una matriz.
<- getValues(img) # Transformación de los valores de los píxeles en una matriz
img_matrix <- na.omit(img_matrix) # Eliminación en dicha matriz de todos los valores NA. img_matrix
En este caso, importaremos los datos de las muestras mediante un fichero de texto “.csv”, y los dividiremos en dos subconjuntos, de entrenamiento (70 %) y de verificación (30%).
<- read.csv("./datos/muestra_datos.csv", header =TRUE)
datos_muestra
<-sample(2,
muestra nrow(datos_muestra),
replace =TRUE,
prob = c(0.7,0.3))
<- datos_muestra[muestra==1, ]
entrenamiento <- datos_muestra[muestra==2, ] verificacion
Ahora aplicaremos el modelo randomForest al conjunto de datos de
entrenamiento con, digamos, 2000 árboles de decisión (argumento
ntree
) de la función randomForest()
. Como el
propósito es la clasificación, el método de "Class"
se
asigna en el argumento strata
. La columna “Class” es la
variable dependiente (que debermos transformar en factor) mientras que
las columnas restantes son las independientes.
<- randomForest(as.factor(Clase)~.,
rf data = entrenamiento,
ntree = 2000,
strata = "Class",
importance = T)
rf
La salida de la función randomForest proporciona una estimación OOB de la tasa de error y la matriz de confusión. La estimación de error OOB o Out-Of-Bag proporciona una tasa de clasificación errónea. La matriz de confusión y estimación de error OOB muestra que la precisión del modelo RandomForest construido es buena. Se puede mejorar aún más el rendimiento del modelo RandomForest seleccionando el número óptimo de árboles. Esto se puede obtener de la gráfica de error del modelo de bosque aleatorio como se muestra a continuación.
plot(rf, lwd=rep(2, 3))
legend("right",
legend = c("OOB Error", "FPR", "FNR"),
lwd=rep(2, 3),
lty = c(1,2,3),
col = c("black", "red", "green"))
La gráfica muestra que la tasa de error es constante a partir de unos 700 árboles. Por lo tanto, es posible mejorar el modelo alterando el número de árboles.
<- randomForest(as.factor(Clase)~.,
rf data = entrenamiento,
ntree = 700,
strata = "Class",
importance = T)
rf
Al haber proporcionado un número óptimo de árboles, la estimación OOB de la tasa de error ha disminuido del 1,6 % al 1,29 %, por lo que se puede aplicar este último modelo entrenado sobre el conjunto de datos de verificación y evaluar su precisión.
<- predict(rf, verificacion) pred
Creación de la matriz de confusión del conjunto de datos de verificacion
confusionMatrix(pred, as.factor(verificacion$Clase))
La salida de la función confusionMatrix muestra claramente que el la precisión de clasificación del modelo evaluado sobre el conjunto de datos de verificación es buena.
La importancia de cada variable independiente en el proceso de clasificación se obtiene creando el siguiente gráfico.
varImpPlot(rf)
En este gráfico, MeanDecrease Accuracy representa cuánto disminuye la precisión del modelo si se omite esa variable en particular. De manera similar, MeanDecreaseGini mide la importancia de la variable según el índice de impureza de Gini, que se utiliza para calcular las divisiones en los árboles. Un valor alto indica una mayor importancia de esa variable en el modelo. Ambos gráficos muestran que la reflectancia en las bandas NIR y SWIR-1 tiene la mayor importancia en comparación con la reflectancia en las restantes bandas en el modelo desarrollado.
La matriz que refleja la importancia de cada variable en cada clase
se puede obtener de la función importance()
.
importance(rf)
El número de veces que cada variable es usada por el algorítmo
randomForest puede ser obtenido con la funcion varUsed
.
varUsed(rf)
Si desea conocer el efecto de los valores de reflectancia de una banda particular al predecir una clase específica, se puede crear un diagrama de dependencia parcial. Este gráfico representa el efecto marginal de una variable en la probabilidad de una clase. Supongamos que desea verificar el efecto de la variación en los valores de reflectancia NIR al predecir la clase “vegetación” en el modelo.
partialPlot(rf, entrenamiento, NIR, "Vegetacion")
partialPlot(rf, entrenamiento, Red, "Agua")
En este gráfico, podemos ver que cuando el valor de reflectancia NIR es superior a 0,35, existe la mayor probabilidad de predicción de la clase ‘Vegetación’. De manera similar, se puede crear un gráfico de dependencia parcial de la clase “Water” en la banda roja
Si desea obtener un árbol de decisión único del modelo randomForest,
debe utilizarse la función getTree()
. En el siguiente
ejemplo se extrae el primer árbol de decisión del modelo
desarrollado.
getTree(rf, 1, labelVar = TRUE)
De esta forma, se pueden derivar varios parámetros del modelo desarrollado. Ahora usaremos este modelo para clasificar la matriz.
<- predict(rf, img_matrix) matrix_rf
A continuación esta matriz clasficada se convertirá a formato raster.
Primero se crea un fichero raster vacío.
<- raster(img)
raster_rf # DispLay attri.butes of the cr·eateci raster
raster_rf
summary(raster_rf) #
A continuación, el objeto raster vacío se rellenará con los píxiles
clasificados según el modelo randomForest, eliminando inicialmente todos
los valores perdidos de la matriz img_matrix
.
<- which(!is.na(img_matrix))
i
<- matrix_rf
raster_rf[i]
#Display attributes of classi.fied raster raster_rf
A continuación se dibujará la imagen con los píxeles clasificados según los resultados de randomForest.
ggR(raster_rf, geom_raster = TRUE) + # PLot with ggplot2
ggtitle("Clasificación con randomForest") + # Título
labs(x="Longitud (m)", y="Latitud (m)") + # Etiquetas de los ejes X e Y
theme(plot.title=element_text(hjust =0.5, # Título: alineación en el centro
size =30), # Título: tamaño
axis.title=element_text(size=20), # Etiquetas de los ejes: tamaño
legend.key.size=unit(1, "cm"), # Leyenda: tamaño.
legend.title=element_text(size =20), # Leyenda: tamaño del título.
legend.text =element_text(size =15)) + # Leyenda: tamaño del texto.
scale_fill_gradientn(colors = c("grey", "burlywood", "wheat", "white", "darkgreen", "skyblue"),
name="Clase", labels = levels(datos_muestra$Clase))
Para no favorecer la acumulación de objetos en el Global environment, se eliminan a continuación.
rm(list=ls())
Muchos algoritmos convencionales se utilizan para la clasificación supervisada de datos ráster, como el clasificador de distancia mínima respecto a la media (Minimum-distance-to-means classifier), el clasificador de paralelepípedo (parallelepiped classifier) y el clasificador de máxima verosimilitud (Maximum likelihood classifier). En este ejercicio sólo se trabajará con el último de ellos.
La técnica de clasificación de máxima verosimilitud clasifica los píxeles desconocidos en función de su probabilidad. Este método asume que los datos de cada clase en cada banda poseen una distribución gaussiana o normal. Este método calcula la función de densidad de probabilidad para cada clase en cada banda, además del vector medio y la matriz de covarianza de cada clase. Para la clasificación de un píxel, se calcula la probabilidad del valor del píxel para cada clase, asignándose el píxel se asignará a la clase que tenga la probabilidad más alta. De esta forma, este método considera tanto la varianza como la covarianza de las clases al clasificar el píxel desconocido. Debido a ello, es un método computacionalmente intensivo.
Usaremos el mismo ráster multibanda utilizado en el caso de la clasificación con CART
<- brick("./datos/input_img.tif")
img names(img) <- c("Coastal_aerosol", "Blue", "Green", "Red", "NIR", "SWIR1", "SWIR2", "cirrus", "TIRS1", "TIRS2")
img
<- readOGR(dsn = "./datos", "training_points")
datos_entrenamiento datos_entrenamiento
Dado que poseen sistemas de proyección diferentes, se reproyectará el fichero en formato vectorial al mismo CRS del fichero raster.
<- spTransform(datos_entrenamiento,
datos_entrenamiento CRS(projection(img)))
compareCRS(img, datos_entrenamiento)
Para la implementación de este método se utiliza la función
superClass()
del paquete RStoolbox y se
asigna el método en su argumento model = "mlc"
.
trainData
. Si no se proporcionan datos de validación por
separado, los datos se dividen en datos de entrenamiento y validación.
En el ejemplo de código R dado, se asigna un valor de 0,8 al argumento
trainPartition
, lo que indica que el 80 % de los datos se
usarán para el entrenamiento y el 20 % de los datos para la
validación.responseCol
.nSamples
.mode
.kfold
.verbose
se establece como VERDADERO, el
progreso y las estadísticas se imprimen durante la ejecución del
código.<- superClass(img,
img_mlc trainData = datos_entrenamiento,
responseCol = "Class",
nsamples = 30,
model = "mlc",
mode = "classification",
trainPartition =0.8,
kfold=5,
verbose =TRUE)
El resultado proporciona la matriz de confusión, la precisión y el coeficiente Kappa de los datos de entrenamiento y validación, y las estadísticas por clase de los datos de validación. En el ejemplo presente, la matriz de confusión, la precisión y el coeficiente Kappa indican que la precisión de la clasificación es buena.
Las estadísticas por clase incluyen otros parámetros, como sensibilidad, especificidad, valor predictivo positivo, valor predictivo negativo, prevalencia, tasa de detección, prevalencia de detección y precisión equilibrada.
Los valores de todos los parámetros discutidos indican que el método utilizado ha clasificado la imagen ráster con un nivel aceptable de precisión.
$map # Atributos del fichero resultante del proceso de clasificación
img_mlc
$classMapping # Correspondencia entre los identificadores de clase img_mlc
Cuando la imagen ráster clasificada y las etiquetas respectivas están disponibles, el usuario no necesita especificar la clase en función de la interpretación visual, a diferencia de la clasificación no supervisada. A continuación se dibuja la imagen clasificada y se asignan los colores respectivos a las etiquetas de clase.
ggR(img_mlc$map, geom_raster = TRUE) +
ggtitle("Clasificación supervisada según el método de máxima verosimilitud") + # 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"), # size of the legend
legend.title = element_text(size = 20), # size of Legend tittle
legend.text = element_text (size = 10)) + # size of legend text
# Dibuja los colores según un gradiente
scale_fill_manual(values = c("grey", "burlywood", "wheat", "white","darkgreen","skyblue"),
name = "Class")