LABORATORIO 7. CLASIFICACIÓN DE IMÁGENES DE SATÉLITE: INTRODUCCIÓN
G174: FOTOINTERPRETACIÓN Y TELEDETECCIÓN
Domingo F. Rasilla
2025-05-12
💡 OBJETIVOS:
Comprender los fundamentos de la clasificación no supervisada y entender en qué se diferencia de la supervisada.
Aplicar diferentes algoritmos de clasificación no supervisada (k-means, ISODATA, Random Forest etc.) con R para clasificar los píxeles en grupos o clases espectrales.
Definir adecuadamente el número de clases (k) y justificarlo.
Evaluar y visualizar los resultados, representándolos espacialmente.
Comparar diferentes resultados variando parámetros como el número de clases.
Identificar visualmente qué zonas se agrupan en cada clase espectral.
Generar productos y exportar resultados en forma de ficheros ráster, cartografía temática o tablas con estadísticas básicas de los grupos clasificados (áreas por clase, histogramas de frecuencia, etc.).
Interpretar las clases resultantes en función de la imagen original y su contexto.
💡 MATERIALES PARA LA ACTIVIDAD:
Los materiales para el desarrollo de esta actividad son los siguientes.
► Imágenes satélite. Es la misma imagen multibanda utilizada en la actividad previa.
INTRODUCCIÓN
La clasificación no supervisada es una forma de clasificación basada en píxeles y es esencialmente una clasificación automatizada, si bien el usuario tiene control sobre ciertas entradas. Esto incluye el número de clases, el número máximo de iteraciones (que es la cantidad de veces que se ejecuta el algoritmo de clasificación) y el porcentaje de umbral de cambio, que especifica cuándo finalizar el procedimiento de clasificación
La clasificación no supervisada es la técnica de clasificación en la que el usuario no proporciona información previa. Las imágenes se clasifican en función de sus características espectrales, de modo que los píxeles que tienen propiedades espectrales similares se agrupan en la misma clase o categoría.
Las clases espectrales se crean únicamente en función de la información numérica de los datos (es decir, los valores de píxeles para cada una de las bandas o índices). Los algoritmos de agrupamiento se utilizan para determinar la agrupación estadística natural de los datos. Los píxeles se agrupan según su similitud espectral. La computadora utiliza el espacio de características para analizar y agrupar los datos en clases. Pase el cursor sobre la imagen siguiente para ver cómo la computadora podría usar el espacio de funciones para agrupar los datos en diez clases.
Una vez clasificados los datos, el usuario debe interpretar, etiquetar y codificar con colores las clases en consecuencia.
Ventajas
Bastante rápida y fácil de ejecutar.
No requiere trabajo de campo ni conocimiento previo extenso del área.
Las clases se crean basándose únicamente en información espectral, por lo que no son tan subjetivas como la interpretación visual manual.
Útil en áreas nuevas, poco conocidas o con pocos datos.
Desventajas
Las clases espectrales no siempre corresponden a clases reales.
El usuario también tiene que dedicar tiempo a denominar, interpretar y asociar las categorías obtenidas a coberturas reales.
Las propiedades espectrales de las clases también pueden cambiar con el tiempo, por lo que no siempre se puede utilizar la misma información de clase al pasar de una imagen a otra.
Gran incertidumbre si hay muchas cubiertas con firmas similares.
Casos ideales de uso
Exploración preliminar del paisaje.
Estudios en zonas desconocidas.
Cuando no se cuenta con datos de referencia.
Los algorítmos de agrupamiento
La clasificación no supervisada utiliza algoritmos de agrupamiento. El propósito de estos algoritmos es minimizar la variabilidad dentro de cada categoría y maximizar la separación entre categorías. Es labor del usuario interpretar y etiquetar las cateogorías, convirtiendo las clases espectrales en clases de información. Este método de clasificación requiere pocos datos, por lo que es un procedimiento rápido y sencillo, que elimina la subjetividad y que puede ser especialmente útil cuando el usuario no tiene ningún conocimiento previo sobre el área de estudio, pero requiere un esfuerzo suplementario en la identificación de las clases de información, ya que podría no reconocerse a simple vista.
K-means es el método más utilizado para realizar una clasificación no supervisada, al que se añadirán en este apartado otros métodos de uso cada vez más frecuentes, como CLARA o Random forest.
PASOS PRELIMINARES
Para realizar una clasificación no supervisada de datos ráster en R, se requieren los siguientes paquetes.
library(terra) # Manejo de datos en formato ráster.## Warning: package 'terra' was built under R version 4.3.3## terra 1.8.21library(ggplot2) # Visualización de datos ráster## Warning: package 'ggplot2' was built under R version 4.3.3library(RStoolbox) # Clasificación y visualización de imágenes siguiendo el método K-means## Warning: package 'RStoolbox' was built under R version 4.3.3## This is version 1.0.0 of RStoolboxlibrary(cluster) # Clasificación según el según el método CLARA.
library(randomForest) # Clasificación siguiendo el método RandomForest.## Warning: package 'randomForest' was built under R version 4.3.3## randomForest 4.7-1.2## Type rfNews() to see new features/changes/bug fixes.##
## Attaching package: 'randomForest'## The following object is masked from 'package:ggplot2':
##
## marginLos datos que usaremos en este tema se encuentran en la carpeta
/datos/no_supervisada/ que se accede mediante la siguiente
ruta, que se comprueba
A continuación se importará la imagen en formato ráster sobre la que se realizará una clasificación no supervisada
imagen <- rast("D:/G174_2025/LABORATORIO_7_Clasificacion_imagenes/datos/no_supervisada/ejemplo_clas_no_supervisada.tif")
imagen## class : SpatRaster
## dimensions : 678, 772, 6 (nrow, ncol, nlyr)
## resolution : 30, 30 (x, y)
## extent : 432345, 455505, 2905185, 2925525 (xmin, xmax, ymin, ymax)
## coord. ref. : WGS 84 / UTM zone 44N (EPSG:32644)
## source : ejemplo_clas_no_supervisada.tif
## names : Band 1, Band 2, Band 3, Band 4, Band 5, Band 6Esa imagen contiene una escena de Landsat 8 con 6 bandas, cuyos nombres se cambian a continuación.
names(imagen) <- c("blue", "green", "red", "nir", "swir1", "swir2")
imagen## class : SpatRaster
## dimensions : 678, 772, 6 (nrow, ncol, nlyr)
## resolution : 30, 30 (x, y)
## extent : 432345, 455505, 2905185, 2925525 (xmin, xmax, ymin, ymax)
## coord. ref. : WGS 84 / UTM zone 44N (EPSG:32644)
## source : ejemplo_clas_no_supervisada.tif
## names : blue, green, red, nir, swir1, swir2A continuación será representada en forma de imagen en falso color
ggRGB(imagen, 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
Las imágenes de satélite deben convertirse en una matriz de datos
antes de aplicar cualquiera de los algoritmos disponibles. Este
dataframe se crea con la función terra:as.data.frame(),
teniendo la precaución de eliminar los valores ausentes (No data) con el
argumento .
df <- as.data.frame(imagen, # Nombre del ráster multibanda.
xy = TRUE, # Añade las coordendas x e y
na.rm = TRUE) # Elimina valores ausentesPara grandes imágenes, se puede extraer una muestra para una clasificación exploratoria:
set.seed(123)
df_muestra <- df[sample(1:nrow(df), 10000), ] # por ejemplo, 10.000 píxelesCLASIFICACIÓN USANDO EL ALGORITMO K-MEANS
Es el método más utilizado para la clasificación no supervisada de ráster. En el algoritmo K-Means, el número de grupos (\(k\)) debe ser predefinido por el usuario. El algorítmo agrupa los píxeles dentro de cada cateogría de acuerdo a su proximidad al centroide más cercano, que es la posición media de todos los casos de un grupo, en el espacio de características multiespectrales. Una vez asignados todos los píxeles a cada uno de los grupos, se recalculan de nuevo los centroides del grupo. Este procedimiento se itera hasta que no se encuentre una diferencia significativa en los centroides o se alcance el número predefinido de iteraciones. De esta forma, este algoritmo optimiza la posición de los centroides. El inconveniente de este algorítmo es su sensibilidad a los valores atípicos.
Con el paquete cluster
Recordar de nuevo que el procedimiento de clasificación se realiza sobre las bandas del dataframe (no debe incluirse las coordeandas):
bandas_kmeans <- df[, c("blue", "green", "red", "nir", "swir1", "swir2")] # Seleccionar solo columnas espectralesLa agrupación de k-means en R se puede realizar utilizando la función
cluster::kmeans(), que se acompaña de los siguientes
argumentos:
El número de grupos se define en el argumento
centers.El número máximo de iteraciones permitidas se menciona en el argumento
iter.max.El algoritmo a utilizar se menciona en el argumento
algorithm. El algoritmo de Lloyd es el más utilizado.
set.seed(999) # Semillas aleatorias
kmeans <- kmeans(bandas_kmeans, # Bandas a clasificar
centers = 5, # Número de grupos o clusters
nstart = 10, # Número de puntos de partida.
iter.max = 500, # número de iteraciones
algorithm = 'Lloyd') # Algoritmo
str(kmeans)## List of 9
## $ cluster : Named int [1:523416] 4 4 4 4 4 4 4 4 4 3 ...
## ..- attr(*, "names")= chr [1:523416] "1" "2" "3" "4" ...
## $ centers : num [1:5, 1:6] 0.158 0.165 0.149 0.153 0.137 ...
## ..- attr(*, "dimnames")=List of 2
## .. ..$ : chr [1:5] "1" "2" "3" "4" ...
## .. ..$ : chr [1:6] "blue" "green" "red" "nir" ...
## $ totss : num 224993
## $ withinss : num [1:5] 4918 6294 6465 4638 2925
## $ tot.withinss: num 25240
## $ betweenss : num 2e+05
## $ size : int [1:5] 196425 60447 81913 175528 9103
## $ iter : int 102
## $ ifault : NULL
## - attr(*, "class")= chr "kmeans"La salida de la función kmeans proporciona una gran
cantidad de información:
kmeans$cluster # Número de casos
kmeans$centers # Aparecen los valores de los centroides
kmeans$totss # Total sum of squares
kmeans$withinss # Within-cluster sum of squares
kmeans$totalwithins # Total within-cluster sum of squares
kmeans$betweenss # Between-cluster sum of squares
kmeans$size # Número de píxeles en cada cluster
kmeans$iter # Número de iteraciones
kmeans$ifault # Indicador de un posible fallo del algorítmo Lo más relevante es lo siguiente:
Tamaño de cada grupo.
Valor medio de las variables (bandas) que conforman cada grupo.
El vector de agrupación, es decir la identidad de grupo de cada celda o píxel clasificado, dentro del grupo suma de cuadrados por grupo.
La “suma de cuadrados por conglomerado dentro del conglomerado” proporciona la medida de la varianza total. Los resultados indican que este algoritmo captura el 91,8% de la varianza de los datos.
Para transformar la matriz de grupos en un objeto SpatRaster debemos seguir el siguiente procedimiento
knr <- imagen[[1]]
values(knr) <- kmeans$cluster
knr## class : SpatRaster
## dimensions : 678, 772, 1 (nrow, ncol, nlyr)
## resolution : 30, 30 (x, y)
## extent : 432345, 455505, 2905185, 2925525 (xmin, xmax, ymin, ymax)
## coord. ref. : WGS 84 / UTM zone 44N (EPSG:32644)
## source(s) : memory
## varname : ejemplo_clas_no_supervisada
## name : blue
## min value : 1
## max value : 5Este objeto ráster es un ráster cuantitativo. Para transformar los números en variables categóricas debemos realizar lo siguiente:
- Convertir a factor (variable categórica)
knr <- as.factor(knr)- Asignar etiquetas opcionales a cada clase
levels(knr) <- data.frame(ID = 1:5,
clase = c("Arena","Suelos salinos","Edificaciones","Agua", "Cultivos"))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’)
mis_colores <- c("Arena" = "#1f78b4", # azul
"Suelos salinos" = "#33a02c", # verde
"Edificaciones" = "#ff7f00", # naranja
"Agua" = "#6a3d9a", # púrpura
"Cultivos" = "#e31a1c" # rojo
)
# Verifica los niveles del raster como factor
levels(knr)## [[1]]
## ID clase
## 1 1 Arena
## 2 2 Suelos salinos
## 3 3 Edificaciones
## 4 4 Agua
## 5 5 CultivosggR(knr, geom_raster = TRUE, stretch = "lin", forceCat = TRUE) +
ggtitle("Clasificación no supervisada. Método K-means") + # 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 = mis_colores,
name = "Clase") +
theme_minimal()
Con el paquete RStoolbox
Una posible alternativa es realizar directamente la clasificación con
K-means sobre el objeto ráster mediante el argumento
algorithm = "Lloyd"en la función
RStoolbox::unsuperClass() del paquete
RStoolbox.
kmeans_img <- unsuperClass(imagen,
nSamples = 10000,
nClasses = 5,
nStarts = 200,
algorithm ="Lloyd")## Warning: did not converge in 100 iterations
## Warning: did not converge in 100 iterations
## Warning: did not converge in 100 iterations
## Warning: did not converge in 100 iterations
## Warning: did not converge in 100 iterations
## Warning: did not converge in 100 iterations
## Warning: did not converge in 100 iterations
## Warning: did not converge in 100 iterations
## Warning: did not converge in 100 iterations
## Warning: did not converge in 100 iterations
## Warning: did not converge in 100 iterations
## Warning: did not converge in 100 iterations
## Warning: did not converge in 100 iterations
## Warning: did not converge in 100 iterations
## Warning: did not converge in 100 iterations
## Warning: did not converge in 100 iterations
## Warning: did not converge in 100 iterations
## Warning: did not converge in 100 iterations
## Warning: did not converge in 100 iterations
## Warning: did not converge in 100 iterations
## Warning: did not converge in 100 iterations
## Warning: did not converge in 100 iterations
## Warning: did not converge in 100 iterations
## Warning: did not converge in 100 iterations
## Warning: did not converge in 100 iterations
## Warning: did not converge in 100 iterationskmeans_img$model## K-means clustering with 5 clusters of sizes 1312, 3723, 3251, 186, 1528
##
## Cluster means:
## blue green red nir swir1 swir2
## 1 0.1639963 0.1702606 0.1856635 0.2782294 0.31017683 114.3327
## 2 0.1577765 0.1618060 0.1744639 0.2761161 0.29728440 113.7515
## 3 0.1530809 0.1533615 0.1598658 0.2664835 0.26868863 113.2652
## 4 0.1373154 0.1290701 0.1147543 0.1269839 0.09502989 110.6453
## 5 0.1489939 0.1447215 0.1428879 0.2564182 0.22703302 112.6845
##
## Clustering vector:
## NULL
##
## Within cluster sum of squares by cluster:
## [1] 121.60054 88.48678 83.50666 63.79018 120.87192
## (between_SS / total_SS = 89.0 %)
##
## Available components:
##
## [1] "centers" "totss" "withinss" "tot.withinss" "betweenss"
## [6] "size" "iter" "ifault"kmeans_img$map## class : SpatRaster
## dimensions : 678, 772, 1 (nrow, ncol, nlyr)
## resolution : 30, 30 (x, y)
## extent : 432345, 455505, 2905185, 2925525 (xmin, xmax, ymin, ymax)
## coord. ref. : WGS 84 / UTM zone 44N (EPSG:32644)
## source(s) : memory
## name : class_unsupervised
## min value : 1
## max value : 5CLASIFICACIÓN APLICANDO EL ALGORÍTMO “CLARA”.
CLARA es la abreviatura de Clustering LARge Application,
método de agrupamiento propuesto por Kaufmann y Rousseeuw (1990). A
diferencia de K-means, es más robusto con los valores atípicos porque
toma medoides en lugar de centroides. A diferencia de los
centroides, el meoide de una clase consiste en la posición
real de un punto que tiene una diferencia mínima con el resto de los
puntos de esa clase, es decir, es el punto central del clúster.
El procedimiento para calcular los meoides consiste en la extracción de numerosas muestras procedentes del conjunto de datos inicial, sobre los que se aplica el algoritmo PAM (Paititioning Around Medoids) para encontrar el número de medoides especificado por el usuario (equivalente al número de grupos o clústers). Este algoritmo calcula los medoides óptimos tomando un conjunto aleatorio inicial de medoides, reemplazadando iterativamente cada uno de ellos hasta obtener el costo de intercambio más bajo. Después de obtener los medoides óptimos, cada observación individual se asigna a su medoide más cercano y se calcula la disimilitud de las observaciones con su meoide más cercano. Este proceso de muestreo y el agrupamiento se itera sucesivas veces y se selecciona como resultado final el agrupamiento con una diferencia mínima.
CLARA se aplica usando la función cluster:clara
disponible en el paquete cluster. Esta función se
acompaña de los siguientes argumentos:
Una matrix creada a partir de los datos raster.
kes el número de conglomerados (clases).sampleses el número de muestras extraidas del conjunto inicial de datos.metrices la unidad de medida que sirve para la separación entre puntos.
clara_fun <- clara(bandas_kmeans,
k = 5,
samples = 1000,
metric = "euclidean")
clara_fun## Call: clara(x = bandas_kmeans, k = 5, metric = "euclidean", samples = 1000)
## Medoids:
## blue green red nir swir1 swir2
## 489409 0.1521807 0.1565304 0.1715786 0.2826200 0.29381654 113.7009
## 504392 0.1503380 0.1524266 0.1640332 0.2671577 0.28779092 113.3135
## 146157 0.1562033 0.1536578 0.1521729 0.2575472 0.21940006 112.8104
## 189700 0.1356411 0.1276448 0.1117124 0.1079210 0.07844345 110.7207
## 326313 0.1614168 0.1699130 0.1927464 0.2832502 0.31742409 114.1121
## Objective function: 0.1637284
## Clustering vector: Named int [1:523416] 1 2 2 2 2 2 2 2 2 3 3 3 3 3 3 3 3 3 ...
## - attr(*, "names")= chr [1:523416] "1" "2" "3" "4" "5" "6" "7" ...
## Cluster sizes: 158404 148794 100480 9843 105895
## Best sample:
## [1] 2852 47234 51674 56658 69796 78326 100641 117645 146157 151428
## [11] 154423 157266 178607 189700 198510 208134 226982 306450 318614 323350
## [21] 324300 326313 330107 337742 341192 353571 372045 372668 375439 385566
## [31] 397634 401100 414382 426801 427001 441018 442112 447271 454803 466176
## [41] 479442 481830 483699 489409 504392 504576 507715 509208 514975 520158
##
## Available components:
## [1] "sample" "medoids" "i.med" "clustering" "objective"
## [6] "clusinfo" "diss" "call" "silinfo" "data"La salida proporcionada por la función incluye también una gran cantidad de información, entre los que cabe destacar:
la función objetivo para el agrupamiento final de todo el conjunto de datos,
el vector de agrupamiento que representa la identidad del grupo de cada píxel clasificado,
el tamaño de cada grupo,
la mejor muestra, es decir, el número de píxel o celda de la mejor muestra utilizada por CLARA algoritmo para la partición final y componentes disponibles.
clara_fun$sample # Número de los píxeles que corresponden a la mejor muestra usada para la partición final## [1] "2852" "47234" "51674" "56658" "69796" "78326" "100641" "117645"
## [9] "146157" "151428" "154423" "157266" "178607" "189700" "198510" "208134"
## [17] "226982" "306450" "318614" "323350" "324300" "326313" "330107" "337742"
## [25] "341192" "353571" "372045" "372668" "375439" "385566" "397634" "401100"
## [33] "414382" "426801" "427001" "441018" "442112" "447271" "454803" "466176"
## [41] "479442" "481830" "483699" "489409" "504392" "504576" "507715" "509208"
## [49] "514975" "520158"clara_fun$medoids # Medoids de los grupos## blue green red nir swir1 swir2
## 489409 0.1521807 0.1565304 0.1715786 0.2826200 0.29381654 113.7009
## 504392 0.1503380 0.1524266 0.1640332 0.2671577 0.28779092 113.3135
## 146157 0.1562033 0.1536578 0.1521729 0.2575472 0.21940006 112.8104
## 189700 0.1356411 0.1276448 0.1117124 0.1079210 0.07844345 110.7207
## 326313 0.1614168 0.1699130 0.1927464 0.2832502 0.31742409 114.1121clara_fun$i.med # Índices de los medoids## [1] 489409 504392 146157 189700 326313clara_fun$objective # Funcion objetiva para la clasificación final de todo el conjunto de datos## [1] 0.1637284clara_fun$clusinfo # Información de cada cluster## size max_diss av_diss isolation
## [1,] 158404 0.3018528 0.1096418 0.7782516
## [2,] 148794 0.3756003 0.1251780 0.9683909
## [3,] 100480 1.0602446 0.2095762 2.0871385
## [4,] 9843 1.1133797 0.5441528 0.5300553
## [5,] 105895 1.7885972 0.2199378 4.3337588El componente clusinfo proporciona información
sobre:
El número de observaciones en cada grupo (su tamaño).
La disimilitud máxima y su promedio entre las observaciones en cada grupo y su respectivo medoide, “max_diss” y “av_diss”’ respectivamente.
La ratio entre la máxima disimilitud entre las observaciones del grupo y el medoide del grupo con respecto a la mínima disimilitud entre el medoide del grupo y el medoide de cualquier otro grupo (“isolation”). Un valor pequeño de este último indicador indica que el clúster está bien separado de los otros clústeres. La matriz de disimilitud se puede obtener del componente ‘diss’ como se muestra a continuación.
clara_fun$dissLa lista con información acerca del “silouette width” para la mejor
muestra es proporcionada por el argumento silinfo.
clara_fun$silinfo## $widths
## cluster neighbor sil_width
## 489409 1 2 0.80953556
## 341192 1 2 0.80630632
## 483699 1 2 0.80605881
## 442112 1 2 0.79988408
## 426801 1 2 0.78364428
## 507715 1 2 0.78329714
## 454803 1 2 0.76383695
## 481830 1 5 0.75895883
## 520158 1 2 0.75701657
## 323350 1 2 0.70563637
## 414382 1 2 0.69003656
## 479442 1 2 0.63375375
## 372045 1 5 0.59695546
## 178607 1 2 0.46322058
## 375439 1 2 0.44309481
## 504392 2 1 0.65211713
## 226982 2 3 0.64760179
## 117645 2 1 0.63794639
## 447271 2 3 0.62713734
## 56658 2 3 0.59191172
## 397634 2 1 0.56378160
## 324300 2 1 0.16181161
## 51674 3 2 0.60484607
## 100641 3 2 0.60166094
## 514975 3 2 0.59466430
## 146157 3 2 0.57122256
## 69796 3 2 0.54860573
## 330107 3 2 0.46104964
## 151428 3 2 0.43638858
## 198510 3 2 0.41146637
## 47234 3 4 0.33138285
## 154423 3 2 0.11050774
## 441018 3 2 -0.09128124
## 157266 4 3 0.93806770
## 189700 4 3 0.93399380
## 372668 5 1 0.58458717
## 353571 5 1 0.58377245
## 509208 5 1 0.58205623
## 504576 5 1 0.57528556
## 337742 5 1 0.54245089
## 326313 5 1 0.46743240
## 306450 5 1 0.45560434
## 208134 5 1 0.44718061
## 466176 5 1 0.44633789
## 78326 5 1 0.43431627
## 401100 5 1 0.35476396
## 2852 5 1 0.29219550
## 318614 5 1 0.08905023
## 385566 5 1 -0.15832724
## 427001 5 1 -0.28717212
##
## $clus.avg.widths
## [1] 0.7067491 0.5546154 0.4164103 0.9360307 0.3606356
##
## $avg.width
## [1] 0.5269131El resultado del procedimiento CLARA también es una matriz, por lo que deberás ser convertido en un raster, siguiendo el mismo procedimiento mencionado en K-means.
clara_raster <- imagen[[1]]
values(clara_raster) <- clara_fun$clustering
clara_raster## class : SpatRaster
## dimensions : 678, 772, 1 (nrow, ncol, nlyr)
## resolution : 30, 30 (x, y)
## extent : 432345, 455505, 2905185, 2925525 (xmin, xmax, ymin, ymax)
## coord. ref. : WGS 84 / UTM zone 44N (EPSG:32644)
## source(s) : memory
## varname : ejemplo_clas_no_supervisada
## name : blue
## min value : 1
## max value : 5plot(clara_raster)
Después de crear el archivo ráster clasificado, se sepresenta cartográficamente, se identifican las clases y se asignan las etiquetas a cada clase a través de la interpretación visual.
ggR(clara_raster, geom_raster = TRUE, stretch = "lin") +
ggtitle("Clasificación no supervisada según el método CLARA") + # 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 gradiente
scale_fill_gradientn(name = "Class",
colours = c("white","forestgreen", "gray", "black", "wheat", "burlywood"),
labels = c("Suelos salinos","Vegetacion","Edificaciones","Suelos arenosos", "Cultivos"))
EVALUACIÓN VISUAL DE LAS CATEGORÍAS
El problema fundamental al que se enfrenta el investigador, sobre todo si no tiene acceso a la zona de estudio, es determinar la correspondencia entre las clases espectrales obtenidas a partir del proceso de clasificación, y los diferentes usos del suelo.
Se muestra a continuación un procedimiento para clarificar este proceso.
umbrales <- c(1, NA,
2, NA,
3, NA,
4, 1,
5, NA,
6, NA)
matriz <- matrix(umbrales,
ncol=2,
byrow = TRUE)
clase4 <- classify(knr, matriz, right=TRUE)
plot(clase4)
A continuación, podemos superponer ambas imágenes para precisar el tipo de superficie representada en esa categoría.
imagenRGB <- imagen[[c(4,3,2)]]
plotRGB(imagenRGB, axes=TRUE, stretch="lin")
plot(clase4, add=TRUE, legend=FALSE, alpha = 0.4, col = "blue")
ANÁLISIS ESTADÍSTICO DE LA FRECUENCIA DE CATEGORÍAS
Una vez obtenidas las categorías, podría ser interesante extraer algún tipo de información estadística sobre la frecuencia con la que aparecen en la escena analizada. Para realizarlo, es necesario transformar las categorías en un factores, para posteriormente representarlas gráficamente y calcular algún valor adicional.
1er paso: en transformar los valores correspondientes a cada
categoría en un factor, con la función
terra::as.factor.
knr_factor <- as.factor(knr) # Define los valores del ráster como factores (variable categórica)
knr_factor## class : SpatRaster
## dimensions : 678, 772, 1 (nrow, ncol, nlyr)
## resolution : 30, 30 (x, y)
## extent : 432345, 455505, 2905185, 2925525 (xmin, xmax, ymin, ymax)
## coord. ref. : WGS 84 / UTM zone 44N (EPSG:32644)
## source(s) : memory
## varname : ejemplo_clas_no_supervisada
## categories : clase
## name : clase
## min value : 1
## max value : 52º paso: creación de la rat o tabla de atributos del
ráster (“Raster Attribute Table”) utilizando para ello la función
levels()
rat <- levels(knr_factor)[[1]]
rat## ID clase
## 1 1 1
## 2 2 2
## 3 3 3
## 4 4 4
## 5 5 53er paso: se añade a cada categoría el siguiente texto, y dicho texto se extiende al raster con las categorías numéricas.
rat$legend <- c("Clase 1", "Clase 2", "Clase 3", "Clase 4", "Clase 5")
levels(knr_factor) <- rat 4º paso: creación de los colores que identifican cada categoría en el diagrama de barras
col <- c("wheat", "white", "skyblue","green","blue")5º paso: representación gráfica de la frecuencia mediante un diagrama de barras
options(scipen=999) #Sustituye la notación científica por 0.
barplot(knr_factor,
ylim = c(0, 200000),
main = "Distribución de categorías de superficie",
col = col,
names.arg = c("Clase 1", "Clase 2", "Clase 3", "Clase 4", "Clase 5"),
las=0)
El gráfico anterior se puede mejorar, añadiendo una etiqueta con el número de casos por categoría.
vals <- values(knr_factor)
pix_1 <- length(subset(vals, vals == 1))
pix_2 <- length(subset(vals, vals == 2))
pix_3 <- length(subset(vals, vals == 3))
pix_4 <- length(subset(vals, vals == 4))
pix_5 <- length(subset(vals, vals == 5))
pix_values <- c(pix_1,pix_2,pix_3,pix_4, pix_5)
bp <- barplot(knr_factor,
ylim = c(0, 250000),
main = "Número de píxeles según categoría",
col = col,
names.arg = c("Clase 1", "Clase 2", "Clase 3", "Clase 4", "Clase 5"),
las=0)
text(bp, y=240000, labels=pix_values) #El valor de y debe ser elegido manualmente
Otra variante es incluir el porcentaje respecto al número de píxeles de la imagen.
pix_total <- length(vals)
area_1_perc <- round(pix_1/pix_total*100, digits = 5)
area_2_perc <- round(pix_2/pix_total*100, digits = 5)
area_3_perc <- round(pix_3/pix_total*100, digits = 5)
area_4_perc <- round(pix_4/pix_total*100, digits = 5)
area_5_perc <- round(pix_5/pix_total*100, digits = 5)
area_perc_covered <- c(area_1_perc,area_2_perc,area_3_perc,area_4_perc, area_5_perc)
bp <- barplot(knr_factor,
ylim = c(0, 250000),
main = "Porcentaje de píxeles ocupados por cada categoría",
col = col,
names.arg = c("Clase 1", "Clase 2", "Clase 3", "Clase 4", "Clase 5"),
las=0)
text(bp, y=240000, labels=area_perc_covered) 
Finalmente, es posible calcular ese porcentaje, pero respecto a la extensión del ráster en kilómetros cuadrados
ex <- ext(knr_factor)
x <- (ex[2]-ex[1])/1000
y <- (ex[4]-ex[3])/1000
area_total <- x*y
area_1 <- round((area_1_perc/100)*area_total, digits = 5)
area_2 <- round((area_2_perc/100)*area_total, digits = 5)
area_3 <- round((area_3_perc/100)*area_total, digits = 5)
area_4 <- round((area_4_perc/100)*area_total, digits = 5)
area_5 <- round((area_5_perc/100)*area_total, digits = 5)
area_covered <- c(area_1,area_2,area_3,area_4, area_5)
bp <- barplot(knr_factor,
main = "Área de cada categoría",
col = col,
names.arg = c("Clase 1", "Clase 2", "Clase 3", "Clase 4", "Clase 5"),
las=0,
ylim = c(0, 250000))
text(bp, y=240000, labels=area_covered) #y chosen manually
Finalmente,
tabla <- as.data.frame(rbind(pix_values, area_perc_covered, area_covered))
names(tabla) <- c("Cl1", "Cl2", "Cl3", "Cl4", "Cl5")