title: “Análisis de Pendiente y Orientación en Guaviare usando R” |
df_print: paged |
Autor: “Jesus Manuel Riaño” |
Date: “20/06/2025” |
output: html_notebook |
Este trabajo tiene como objetivo analizar la topografía del departamento de Guaviare (Colombia), enfocándose en dos características fundamentales del relieve: la pendiente y la orientación del terren,para ello, se utilizó el lenguaje de programación R junto con paquetes especializados en análisis geoespacial y manejo de datos ráster.
Primero se cargargan las bibliotecas necesarias como elevatr, terra, sf, MultiscaleDTM y exactextractr, que permiten trabajar con archivos geográficos, cálculos de pendiente y estadísticas espaciales.
rm(list=ls())
Se inicia a cargar un archivo ráster tipo (.tif) que contiene datos de elevación del terreno de Guaviare, continunado se reduce la resolución para eso usando un promedio por bloques para que el procesamiento sea más eficiente,el archivo contiene 4094 filas y 5633 columnas ,con una resoluccion de 0.00068651 grados en ambos ejes la cubre un área geográfica entre -73.83 y -69.96 de longitud, y 0.35 a 3.16 de latitud esto proveniente del archivo Guaviare.tif con valores desde -553 hasta 977 .
dem <- rast("C:/Users/yorhl/Downloads/DATOS/Guaviare.tif")
print(dem)
## class : SpatRaster
## size : 4094, 5633, 1 (nrow, ncol, nlyr)
## resolution : 0.00068651, 0.00068651 (x, y)
## extent : -73.82812, -69.96101, 0.3518836, 3.162456 (xmin, xmax, ymin, ymax)
## coord. ref. : +proj=longlat +ellps=GRS80 +no_defs
## source : Guaviare.tif
## name : Guaviare
## min value : -553
## max value : 977
En esta caso el raster dem2 tendrá una resolución más baja, ya que se ha hecho un promedio de celdas vecinas para reducir su tamaño y resolución, Esto se usa comúnmente para hacer que los cálculos u visualizaciones sean más rápidos para trabajar con datos menos detallados.
dem2 = terra::aggregate(dem,2, "mean")
## |---------|---------|---------|---------|=========================================
Se importó un archivo vectorial (.gpkg) con la división política del departamento del Guaviare , lo cual permite delimitar el área de estudio este archivo GeoPackage (Guaviare.gpkg), contiene 4 características geoespaciales de tipo MULTIPOLYGON con 11 atributos asociados. La capa utiliza el sistema de coordenadas MAGNA-SIRGAS y cubre un área geográfica delimitada por las coordenadas de la caja delimitadora: xmin: -73.66, ymin: 0.66, xmax: -69.99, ymax: 2.92.
library("sf")
munic <- sf::st_read("C:/Users/yorhl/Downloads/DATOS/Guaviare.gpkg")
## Reading layer `mipio' from data source
## `C:\Users\yorhl\Downloads\DATOS\Guaviare.gpkg' using driver `GPKG'
## Simple feature collection with 4 features and 11 fields
## Geometry type: MULTIPOLYGON
## Dimension: XY
## Bounding box: xmin: -73.66391 ymin: 0.6554126 xmax: -69.99511 ymax: 2.924987
## Geodetic CRS: MAGNA-SIRGAS
Se recortó el DEM para que coincida con los límites de Guaviare. Esto asegura que el análisis se realice solo sobre la zona de interés.
este codigo recorta el raster dem2 según los límites del objeto munic, las áreas fuera de estos límites se enmascaran y se asignan valores NA, limitando el análisis a la zona deseada.
dem3 = terra::crop(dem2,munic, mask=TRUE)
## Warning: [crop] CRS do not match
Esta liena dem_plane es una versión de dem3 reproyectada al sistema de coordenadas EPSG:9377, lo cual es útil para realizar cálculos precisos en ese sistema y facilita la comparación con otros datos en el mismo CRS , este poligono asegura que las ubicaciones de los vértices que sean precisos y estén alineadas con otros datos geoespaciales.
(dem_plane = project(dem3, "EPSG:9377"))
## class : SpatRaster
## size : 1648, 2678, 1 (nrow, ncol, nlyr)
## resolution : 152.4372, 152.4372 (x, y)
## extent : 4926217, 5334444, 1630446, 1881662 (xmin, xmax, ymin, ymax)
## coord. ref. : MAGNA-SIRGAS 2018 / Origen-Nacional (EPSG:9377)
## source(s) : memory
## name : Guaviare
## min value : -59.91631
## max value : 812.65155
Se reproyectaron tanto el DEM como el polígono del municipio a un sistema de coordenadas planas (EPSG:9377), munic a EPSG:9377 (MAGNA-SIRGAS 2018) para que coincida con el sistema de coordenadas del raster dem3. La advertencia “CRS do not match” indica que las capas tienen CRS diferentes y necesitan ser alineadas antes de realizar operaciones como el recorte o análisis conjunto. La transformación asegura que ambas capas puedan ser comparadas o utilizadas juntas.
(munic_plane = sf::st_transform(munic, "EPSG:9377"))
## Simple feature collection with 4 features and 11 fields
## Geometry type: MULTIPOLYGON
## Dimension: XY
## Bounding box: xmin: 4926180 ymin: 1630473 xmax: 5334140 ymax: 1881362
## Projected CRS: MAGNA-SIRGAS 2018 / Origen-Nacional
## dpto_ccdgo mpio_ccdgo mpio_cdpmp dpto_cnmbr mpio_cnmbr
## 1 95 001 95001 GUAVIARE SAN JOSÉ DEL GUAVIARE
## 2 95 015 95015 GUAVIARE CALAMAR
## 3 95 025 95025 GUAVIARE EL RETORNO
## 4 95 200 95200 GUAVIARE MIRAFLORES
## mpio_crslc mpio_tipo mpio_narea mpio_nano
## 1 Decreto Nal 1165 del 7 de Junio de 1976 MUNICIPIO 16757.08 2023
## 2 Ordenanza 001 del 7 de Agosto de 1992 MUNICIPIO 13069.59 2023
## 3 Ordenanza 001 del 7 de Agosto de 1992 MUNICIPIO 12909.12 2023
## 4 Ordenanza 001 del 7 de Agosto de 1992 MUNICIPIO 12808.08 2023
## shape_Leng shape_Area geom
## 1 14.085997 1.360393 MULTIPOLYGON (((5187468 188...
## 2 7.529464 1.061812 MULTIPOLYGON (((5036504 179...
## 3 12.315742 1.047240 MULTIPOLYGON (((5029855 182...
## 4 7.385424 1.039462 MULTIPOLYGON (((5109660 177...
Se calcula dos características del terreno un es la pendiente (cuánto se inclina el terreno) y la orientación (dirección en la que se inclina). Para hacerlo se usa una ventana de 3x3 celdas y expresa los resultados en grados esta salida muestra un raster con dos capas (una para la pendiente y otra para la orientación) la cual tiene una resolución de 152.4372. Los valores de pendiente varían desde un valor muy bajo (prácticamente plano) hasta un máximo de 42.14, y los valores de orientación cubren todo el rango de direcciones posibles, de 0 a 360 grados.
(slp_asp = MultiscaleDTM::SlpAsp(
dem_plane,
w = c(3, 3),
unit = "degrees",
method = "queen",
metrics = c("slope", "aspect"),
na.rm = TRUE,
include_scale = FALSE,
mask_aspect = TRUE
))
## class : SpatRaster
## size : 1648, 2678, 2 (nrow, ncol, nlyr)
## resolution : 152.4372, 152.4372 (x, y)
## extent : 4926217, 5334444, 1630446, 1881662 (xmin, xmax, ymin, ymax)
## coord. ref. : MAGNA-SIRGAS 2018 / Origen-Nacional (EPSG:9377)
## source(s) : memory
## names : slope, aspect
## min values : 6.291552e-04, 1.815208e-04
## max values : 4.214191e+01, 3.599999e+02
esta capa especifica solo la capa de pendiente, excluyendo la orientación y extrae la capa de pendiente del objeto slp_asp, que contiene tanto la pendiente como la orientación del terreno esta salida muestra un raster de una capa con resolución de 152.4372 y una extensión que abarca el área definida en el sistema de coordenadas MAGNA-SIRGAS 2018 (EPSG:9377). Los valores de pendiente varían desde un valor muy bajo de 0.000629 (casi plano) hasta un máximo de 42.14 (una pendiente considerable), representando las variaciones en la inclinación del terreno.
(slope = subset(slp_asp, 1))
## class : SpatRaster
## size : 1648, 2678, 1 (nrow, ncol, nlyr)
## resolution : 152.4372, 152.4372 (x, y)
## extent : 4926217, 5334444, 1630446, 1881662 (xmin, xmax, ymin, ymax)
## coord. ref. : MAGNA-SIRGAS 2018 / Origen-Nacional (EPSG:9377)
## source(s) : memory
## name : slope
## min value : 6.291552e-04
## max value : 4.214191e+01
Este genera un histograma de los valores de pendiente (inclinación del terreno) de la capa slope permite visualizar la distribución de las pendientes en grados en el área de Guaviare. El título es “Pendiente de Guaviare” y el eje X está etiquetado como “Pendiente (en grados)”, mostrando cómo varía la inclinación del terreno en esa región dond epredomina las zonas planas.
library(terra)
terra::hist(slope,
main = "Pendiente de Guaviare",
xlab = "Pendiente (en grados)")
## Warning: [hist] a sample of 23% of the cells was used (of which 46% was NA)
Se muestra la orientacion del terreno
(aspect =subset(slp_asp, 2))
## class : SpatRaster
## size : 1648, 2678, 1 (nrow, ncol, nlyr)
## resolution : 152.4372, 152.4372 (x, y)
## extent : 4926217, 5334444, 1630446, 1881662 (xmin, xmax, ymin, ymax)
## coord. ref. : MAGNA-SIRGAS 2018 / Origen-Nacional (EPSG:9377)
## source(s) : memory
## name : aspect
## min value : 1.815208e-04
## max value : 3.599999e+02
Este codigo terra::hist(aspect).Genera un histograma de la capa de orientación (aspect) del terreno, representando la distribución de los valores de orientación en el área de Guaviare. Los valores de aspecto indican la dirección de la inclinación del terreno en grados, donde 0° y 360° son el norte, y los demás valores corresponden a otras direcciones (este, sur, oeste).
terra::hist(aspect,
main = "Aspecto de Guaviare",
xlab = " Aspecto (en grados)")
## Warning: [hist] a sample of 23% of the cells was used (of which 46% was NA)
Este código calcula la pendiente en porcentaje de una capa raster utilizando la fórmula:slope_perc = tan(slope(pi/180))100),la capa tiene una resolución de 152.44 metros por celda, cubriendo un área geográfica definida. Los valores de pendiente varían entre 0.0011% (prácticamente plana) y 90.49% (casi vertical), y se usa el sistema de referencia MAGNA-SIRGAS 2018 / Origen-Nacional (EPSG:9377).
(slope_perc = tan(slope*(pi/180))*100)
## class : SpatRaster
## size : 1648, 2678, 1 (nrow, ncol, nlyr)
## resolution : 152.4372, 152.4372 (x, y)
## extent : 4926217, 5334444, 1630446, 1881662 (xmin, xmax, ymin, ymax)
## coord. ref. : MAGNA-SIRGAS 2018 / Origen-Nacional (EPSG:9377)
## source(s) : memory
## name : slope
## min value : 0.001098083
## max value : 90.489891938
terra::hist(slope_perc, main = “Pendiente de Guaviare”, xlab = “Pendiente (en porcentaje)”) este codigo genera un histograma de la capa raster slope_perc que muestra la distribución de la pendiente en porcentaje en la región de Guaviare. El título del gráfico es “Pendiente de Guaviare” y el eje X está etiquetado como “Pendiente (en porcentaje)”, lo que permite visualizar la frecuencia de los diferentes valores de pendiente en esa área.
terra::hist(slope_perc,
main = "Pendiente de Guaviare",
xlab = "Pendiente (en porcentaje)")
## Warning: [hist] a sample of 23% of the cells was used (of which 46% was NA)
Este codigo permite slope_perc utilizando la función classify() del paquete terra,la cual esta define en una matriz de rangos de pendiente y clases asociadas,por ejemplo, pendientes entre 0 y 3 se asignan a la clase 1, entre 3 y 7 a la clase 2, y así sucesivamente con el rrsto de datos por otro lado el parámetro right=TRUE asegura que el valor máximo de cada intervalo sea exclusivo en el resultado es un nuevo raster con las clases asignadas según los intervalos de pendiente.
#rc <- classify(slope_perc, c(0, 3, 7, 12, 25,50, 75), include.lowest=TRUE, brackets=TRUE)
m <- c(0, 3, 1,
3, 7, 2,
7, 12, 3,
12, 25, 4,
25, 50, 5,
50, 75, 6,
75, 160, 7)
m <- matrix(m, ncol=3, byrow = TRUE)
rc <- classify(slope_perc, m, right=TRUE)
Mediante la función exact_extract() extrae el valor promedio de la capa raster slope_perc dentro de las áreas definidas por el objeto minic el cual representa los limites de los municipios y el resultado se asigna a la nueva columna mean_slope del dataframe munic los primeros valores extraídos en este caso son 2.275636, 2.279084, 2.531816 y 1.781586, la cual corresponde a la pendiente media por municipio.
(munic$mean_slope <- exactextractr::exact_extract(slope_perc, munic, 'mean'))
## Warning in .local(x, y, ...): Polygons transformed to raster CRS (EPSG:9377)
## | | | 0% | |================== | 25% | |=================================== | 50% | |==================================================== | 75% | |======================================================================| 100%
## [1] 2.275636 2.279084 2.531816 1.781586
El histograma generado con el código muestra la distribución de las pendientes medias de los municipios en Guaviare. En el eje X están los valores de pendiente en porcentaje, y en el eje Y se muestra la frecuencia. Si se observan picos en 1.6, 22-23 y cerca de 1, significa que la mayoría de los municipios tienen pendientes en esos rangos específicos, con mayor concentración en los valores de 1.6 y 2.3, lo que revela tendencias en la pendiente de la región.
hist(munic$mean_slope,
main = "Pendiente media de Guaviare",
xlab = " Pendiente (en porcentaje)")
este codigo ejecutado extrae el valor más frecuente (modo) de la capa
raster rc para cada municipio en munic. El resultado 11111 significa
que, para cada municipio, el valor más común de la capa reclasificada es
11111.
(munic$class <- exactextractr::exact_extract(rc, munic, 'mode'))
## Warning in .local(x, y, ...): Polygons transformed to raster CRS (EPSG:9377)
## | | | 0% | |================== | 25% | |=================================== | 50% | |==================================================== | 75% | |======================================================================| 100%
## [1] 1 1 1 1
El histograma ejecutado muestra la distribución de las clases reclasificadas de la pendiente en Guaviare. La columna class contiene los valores de clase asignados a cada municipio, donde cada clase representa un rango específico de pendiente.
hist(munic$class,
main = "Pendiente reclasificada de Guaviare",
xlab = "Pendiente (como categoría)")
realiza la proyección de la capa raster rc a un nuevo sistema de coordenadas, en este caso, EPSG:4326.como resultado 1657 filas y 2675 columnas, con una resolución de 0.00137 grados. Cubre un área geográfica con coordenadas de longitud entre -73.66418 y -69.99133 y latitud entre 0.6538653 y 2.928971. Usa el sistema de coordenadas WGS 84 (EPSG:4326).
(rc.geo = project(rc, "EPSG:4326"))
## class : SpatRaster
## size : 1657, 2675, 1 (nrow, ncol, nlyr)
## resolution : 0.001373027, 0.001373027 (x, y)
## extent : -73.66418, -69.99133, 0.6538653, 2.928971 (xmin, xmax, ymin, ymax)
## coord. ref. : lon/lat WGS 84 (EPSG:4326)
## source(s) : memory
## name : slope
## min value : 1
## max value : 7
(slope.geo = project(slope_perc, "EPSG:4326"))
## class : SpatRaster
## size : 1657, 2675, 1 (nrow, ncol, nlyr)
## resolution : 0.001373027, 0.001373027 (x, y)
## extent : -73.66418, -69.99133, 0.6538653, 2.928971 (xmin, xmax, ymin, ymax)
## coord. ref. : lon/lat WGS 84 (EPSG:4326)
## source(s) : memory
## name : slope
## min value : 0.02170257
## max value : 80.95211029
Con este codigo se reduce la resolución del raster slope.geo a la mitad utilizando la función aggregate(), promediando los valores en un factor de 2. Luego, crea una paleta de colores para visualizar la pendiente con la función colorNumeric(). Usando leaflet, se genera un mapa interactivo con los valores del raster reducido, aplicando la paleta de colores. También se añade una leyenda para mostrar los valores de pendiente en porcentaje, con una opacidad del 80% en la capa raster.
# Reducir la resolución del raster (por ejemplo, a la mitad)
slope.lowres <- aggregate(slope.geo, fact = 2, fun = mean)
# Usar el raster reducido en el mapa
palblue <- colorNumeric(
c("lightblue", "purple", "yellow2", "orange", "brown2"),
values(slope.lowres),
na.color = "transparent"
)
leaflet() %>%
addTiles() %>%
addRasterImage(slope.lowres, colors = palblue, opacity = 0.8) %>%
addLegend(pal = palblue, values = values(slope.lowres), title = "Pendiente (%)")
En este caso muestra información detallada sobre el entorno de trabajo actual, incluyendo la versión de R, los paquetes cargados, el sistema operativo y otros detalles sobre el entorno de la sesión. Es útil para verificar configuraciones, versiones de librerías y su replicacion.
sessionInfo()
## R version 4.4.1 (2024-06-14 ucrt)
## Platform: x86_64-w64-mingw32/x64
## Running under: Windows 11 x64 (build 22631)
##
## Matrix products: default
##
##
## locale:
## [1] LC_COLLATE=Spanish_Colombia.utf8 LC_CTYPE=Spanish_Colombia.utf8
## [3] LC_MONETARY=Spanish_Colombia.utf8 LC_NUMERIC=C
## [5] LC_TIME=Spanish_Colombia.utf8
##
## time zone: America/Bogota
## tzcode source: internal
##
## attached base packages:
## [1] stats graphics grDevices utils datasets methods base
##
## other attached packages:
## [1] exactextractr_0.10.0 MultiscaleDTM_1.0 terra_1.8-54
## [4] leaflet_2.2.2 sf_1.0-21 elevatr_0.99.0
##
## loaded via a namespace (and not attached):
## [1] sass_0.4.9 utf8_1.2.4 generics_0.1.3 class_7.3-22
## [5] KernSmooth_2.23-24 lattice_0.22-7 digest_0.6.37 magrittr_2.0.3
## [9] rgl_1.3.18 evaluate_1.0.1 grid_4.4.1 fastmap_1.2.0
## [13] jsonlite_1.8.9 e1071_1.7-16 DBI_1.2.3 promises_1.3.2
## [17] purrr_1.0.2 fansi_1.0.6 scales_1.3.0 crosstalk_1.2.1
## [21] codetools_0.2-20 jquerylib_0.1.4 cli_3.6.3 shiny_1.10.0
## [25] rlang_1.1.4 units_0.8-5 munsell_0.5.1 base64enc_0.1-3
## [29] cachem_1.1.0 raster_3.6-32 tools_4.4.1 dplyr_1.1.4
## [33] colorspace_2.1-1 httpuv_1.6.15 png_0.1-8 vctrs_0.6.5
## [37] R6_2.5.1 mime_0.12 proxy_0.4-27 lifecycle_1.0.4
## [41] classInt_0.4-10 htmlwidgets_1.6.4 pkgconfig_2.0.3 progressr_0.15.1
## [45] bslib_0.8.0 pillar_1.9.0 later_1.4.1 glue_1.8.0
## [49] Rcpp_1.0.13 tidyselect_1.2.1 xfun_0.49 tibble_3.2.1
## [53] rstudioapi_0.17.1 knitr_1.49 farver_2.1.2 xtable_1.8-4
## [57] htmltools_0.5.8.1 rmarkdown_2.29 compiler_4.4.1 sp_2.2-0