1 Introducción

Este cuaderno busca ilustrar los atributos geomorfológicos del departamento de Sucre, Colombia, usando modelos de elevación digital (DEM). Se utiliza el paquete Multiscale DTM de R para el procesamiento y análisis de los datos a fin de obtener variables que describan el terreno para comprender de mejor manera la variabilidad topográfica del departamento.

2 Modelos de elevación digital (DEMs)

Puede definirse como la representación digital del relieve de la tierra, donde a cada punto (Vectorial) o cuadrícula (Raster) le corresponde un valor de elevación. Estos modelos se utilizan para visualizar y analizar la topografía de un terreno, puede usarse en distintas disciplinas, entre ellas la agronomía, ya que permite tomar decisiones informadas según las necesidades del usuario.

##Limpieza de memoria e instalación de librerías

Inicialmente, se deben instalar los siguientes paquetes, que serán usados a lo largo del código para el análisis de datos geomorfológicos.

## SETUP
# INSTALL THIS PACKAGE FROM THE CONSOLE, NOT FROM THIS CHUNK
# paquetes = c("MultiscaleDTM", "exactextractr")
# install.packages(paquetes)
rm(list=ls())
library(elevatr)
## elevatr v0.99.0 NOTE: Version 0.99.0 of 'elevatr' uses 'sf' and 'terra'.  Use 
## of the 'sp', 'raster', and underlying 'rgdal' packages by 'elevatr' is being 
## deprecated; however, get_elev_raster continues to return a RasterLayer.  This 
## will be dropped in future versions, so please plan accordingly.
library(sf)
## Linking to GEOS 3.13.1, GDAL 3.11.0, PROJ 9.6.0; sf_use_s2() is TRUE
library(leaflet)
library(terra)
## terra 1.8.54
library(MultiscaleDTM)
library(exactextractr)

3 Datos de elevación del departamento de Sucre

A continuación, usaremos el DEM de Sucre y el shapefile del departamento.

se usa list.files para ubicar nuestro archivo

list.files("ELEVACION.tif")
## character(0)

Inicialmente, se carga el DEM raster de Sucre. Se lee el DEM usando la librería terra. Este archivo se llamará “dem”.

# change the path & DEM filename as needed
(dem = terra::rast("ELEVACION.tif"))
## class       : SpatRaster 
## size        : 223, 140, 1  (nrow, ncol, nlyr)
## resolution  : 0.008333333, 0.008333333  (x, y)
## extent      : -75.7, -74.53333, 8.283333, 10.14167  (xmin, xmax, ymin, ymax)
## coord. ref. : lon/lat WGS 84 (EPSG:4326) 
## source      : ELEVACION.tif 
## name        : ELEVACION 
## min value   :         0 
## max value   :       630

Propiedades del documento “dem”: Tiene clase SpatRaster, ideal para trabajar con datos de tipo espacial raster. Posee una tabla con un total de (número de filas), (número de columnas) y una sola capa. Cada celda tiene una resolución de aproximadamente (resolución en x,y) grados, es decir, se tiene un nivel de detalle bastante alto.

El sistema de referencia espacial es geográfico, con un elipsoide GRS80 (o el que corresponda a tu DEM).

Reducir la resolución del DEM, para evitar problemáticas en la memoria, ya que se agrupan las celdas del archivo DEM original en bloques de 2*2, promediando valores, así, se reduce la resolución y por lo tanto, la carga de la memoria del sistema:

dem2 = terra::aggregate(dem,2, "mean")

Se incluye el geopackage del departamento con las limitaciones del departamento. Se usa la librería sf para esta acción.

(munic <-  sf::st_read("SUCRE2.gpkg"))
## Reading layer `sucre' from data source `D:\GB2\P4\SUCRE2.gpkg' using driver `GPKG'
## Simple feature collection with 26 features and 12 fields
## Geometry type: MULTIPOLYGON
## Dimension:     XY
## Bounding box:  xmin: -75.70603 ymin: 8.277882 xmax: -74.53271 ymax: 10.14548
## Geodetic CRS:  MAGNA-SIRGAS
## Simple feature collection with 26 features and 12 fields
## Geometry type: MULTIPOLYGON
## Dimension:     XY
## Bounding box:  xmin: -75.70603 ymin: 8.277882 xmax: -74.53271 ymax: 10.14548
## Geodetic CRS:  MAGNA-SIRGAS
## First 10 features:
##    dpto_ccdgo mpio_ccdgo mpio_cdpmp dpto_cnmbr mpio_cnmbr
## 1          70        001      70001      SUCRE  SINCELEJO
## 2          70        110      70110      SUCRE BUENAVISTA
## 3          70        124      70124      SUCRE    CAIMITO
## 4          70        204      70204      SUCRE     COLOSO
## 5          70        215      70215      SUCRE    COROZAL
## 6          70        221      70221      SUCRE   COVENIAS
## 7          70        230      70230      SUCRE     CHALAN
## 8          70        233      70233      SUCRE   EL ROBLE
## 9          70        235      70235      SUCRE    GALERAS
## 10         70        265      70265      SUCRE   GUARANDA
##                                                     mpio_crslc mpio_tipo
## 1                                  Ley 47 de Agosto 30 de 1966 MUNICIPIO
## 2                       Ordenanza 9 del 5 de Noviembre de 1968 MUNICIPIO
## 3                             Ordenanza 58 de Abril 28 de 1914 MUNICIPIO
## 4                                                         1907 MUNICIPIO
## 5                                               Ley 44 de 1876 MUNICIPIO
## 6                           Decreto 0063 de febrero 07 de 2002 MUNICIPIO
## 7                            Ordenanza 6 de Octubre 30 de 1968 MUNICIPIO
## 8  Decreto 356 de 98, fallo 5 de Mayo de 99 y Sentencia del 17 MUNICIPIO
## 9                            Ordenanza 4 de Octubre 15 de 1968 MUNICIPIO
## 10                           Ordenanza 2 de Octubre 31 de 1984 MUNICIPIO
##    mpio_narea mpio_nano shape_Leng  shape_Area     NOMBRE
## 1   279.69771      2023  1.1431571 0.023005955  SINCELEJO
## 2   146.47566      2023  0.8805553 0.012052022 BUENAVISTA
## 3   412.46981      2023  1.1794901 0.033885959    CAIMITO
## 4   131.16269      2023  0.7882677 0.010795532     COLOSO
## 5   286.11257      2023  1.4119720 0.023529464    COROZAL
## 6    56.55515      2023  0.4268842 0.004651922   COVENIAS
## 7    83.71001      2023  0.5720244 0.006890902     CHALAN
## 8   198.43440      2023  1.2342195 0.016314220   EL ROBLE
## 9   321.85643      2023  1.1237950 0.026468388    GALERAS
## 10  361.86764      2023  1.4238869 0.029704319   GUARANDA
##                              geom
## 1  MULTIPOLYGON (((-75.49687 9...
## 2  MULTIPOLYGON (((-74.93518 9...
## 3  MULTIPOLYGON (((-75.19855 8...
## 4  MULTIPOLYGON (((-75.3837 9....
## 5  MULTIPOLYGON (((-75.28363 9...
## 6  MULTIPOLYGON (((-75.60529 9...
## 7  MULTIPOLYGON (((-75.35705 9...
## 8  MULTIPOLYGON (((-75.13801 9...
## 9  MULTIPOLYGON (((-74.93398 9...
## 10 MULTIPOLYGON (((-74.56933 8...

Preguardado de los archivos en un nuevo atributo. Esto recorta el DEM al límite del departamento de Sucre.

dem3 = terra::crop(dem2,munic, mask=TRUE)
## Warning: [crop] CRS do not match

3.1 Transformación de coordenadas

Para calcular correctamente la geomorfología del departamento, el DEM debe tener coordenadas planas. Para el año 2020, el Instituto Geográfico Agustín Codazzi (IGAC) adopta el sistema de coordenadas plana con origen único nacional en Colombia, con código EPSG 9377. Se usa el siguiente código para la transformación de coordenadas:

# check the parameters of this terra function here
# https://rdrr.io/github/rspatial/terra/man/project.html
(dem_plane = project(dem3, "EPSG:9377"))
## class       : SpatRaster 
## size        : 113, 70, 1  (nrow, ncol, nlyr)
## resolution  : 1839.834, 1839.834  (x, y)
## extent      : 4702678, 4831466, 2472003, 2679904  (xmin, xmax, ymin, ymax)
## coord. ref. : MAGNA-SIRGAS 2018 / Origen-Nacional (EPSG:9377) 
## source(s)   : memory
## name        : ELEVACION 
## min value   :   1.59552 
## max value   : 530.25793

Como se observa anteriormente, la transformación de coordenadas fue exitosa. Se transforma el objeto vectorial (límites municipales) al mismo sistema de coordenadas planas:

# check the parameters of this sf function here
# https://github.com/rstudio/cheatsheets/blob/main/sf.pdf
(munic_plane = sf::st_transform(munic, "EPSG:9377"))
## Simple feature collection with 26 features and 12 fields
## Geometry type: MULTIPOLYGON
## Dimension:     XY
## Bounding box:  xmin: 4702923 ymin: 2473141 xmax: 4831339 ymax: 2680149
## Projected CRS: MAGNA-SIRGAS 2018 / Origen-Nacional
## First 10 features:
##    dpto_ccdgo mpio_ccdgo mpio_cdpmp dpto_cnmbr mpio_cnmbr
## 1          70        001      70001      SUCRE  SINCELEJO
## 2          70        110      70110      SUCRE BUENAVISTA
## 3          70        124      70124      SUCRE    CAIMITO
## 4          70        204      70204      SUCRE     COLOSO
## 5          70        215      70215      SUCRE    COROZAL
## 6          70        221      70221      SUCRE   COVENIAS
## 7          70        230      70230      SUCRE     CHALAN
## 8          70        233      70233      SUCRE   EL ROBLE
## 9          70        235      70235      SUCRE    GALERAS
## 10         70        265      70265      SUCRE   GUARANDA
##                                                     mpio_crslc mpio_tipo
## 1                                  Ley 47 de Agosto 30 de 1966 MUNICIPIO
## 2                       Ordenanza 9 del 5 de Noviembre de 1968 MUNICIPIO
## 3                             Ordenanza 58 de Abril 28 de 1914 MUNICIPIO
## 4                                                         1907 MUNICIPIO
## 5                                               Ley 44 de 1876 MUNICIPIO
## 6                           Decreto 0063 de febrero 07 de 2002 MUNICIPIO
## 7                            Ordenanza 6 de Octubre 30 de 1968 MUNICIPIO
## 8  Decreto 356 de 98, fallo 5 de Mayo de 99 y Sentencia del 17 MUNICIPIO
## 9                            Ordenanza 4 de Octubre 15 de 1968 MUNICIPIO
## 10                           Ordenanza 2 de Octubre 31 de 1984 MUNICIPIO
##    mpio_narea mpio_nano shape_Leng  shape_Area     NOMBRE
## 1   279.69771      2023  1.1431571 0.023005955  SINCELEJO
## 2   146.47566      2023  0.8805553 0.012052022 BUENAVISTA
## 3   412.46981      2023  1.1794901 0.033885959    CAIMITO
## 4   131.16269      2023  0.7882677 0.010795532     COLOSO
## 5   286.11257      2023  1.4119720 0.023529464    COROZAL
## 6    56.55515      2023  0.4268842 0.004651922   COVENIAS
## 7    83.71001      2023  0.5720244 0.006890902     CHALAN
## 8   198.43440      2023  1.2342195 0.016314220   EL ROBLE
## 9   321.85643      2023  1.1237950 0.026468388    GALERAS
## 10  361.86764      2023  1.4238869 0.029704319   GUARANDA
##                              geom
## 1  MULTIPOLYGON (((4725807 258...
## 2  MULTIPOLYGON (((4787611 260...
## 3  MULTIPOLYGON (((4758336 254...
## 4  MULTIPOLYGON (((4738508 262...
## 5  MULTIPOLYGON (((4749282 259...
## 6  MULTIPOLYGON (((4714038 260...
## 7  MULTIPOLYGON (((4741438 262...
## 8  MULTIPOLYGON (((4765148 257...
## 9  MULTIPOLYGON (((4787599 257...
## 10 MULTIPOLYGON (((4827330 249...

4 Cálculos de los atributos geomorfométricos

Para calcular el aspecto y la pendiente, se debe usar la función SlpAsp de la librería MultiscaleDTM. Inicialmente, se define el parámetro w como el tamaño de la ventana usada para calcular los atributos. Para este caso, w=c(3,3), es decir, una ventana con celdas de 3x3. Continuando con el código, el parámetro method define el tipo de vecindad utilizado para los cálculos, donde “queen” es un tipo de vecindad que incluye aquellos datos que pueden llamarse vecinos directos o intermedios.

# Explain what is w 
# Explain what is method
# Change if needed
(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        : 113, 70, 2  (nrow, ncol, nlyr)
## resolution  : 1839.834, 1839.834  (x, y)
## extent      : 4702678, 4831466, 2472003, 2679904  (xmin, xmax, ymin, ymax)
## coord. ref. : MAGNA-SIRGAS 2018 / Origen-Nacional (EPSG:9377) 
## source(s)   : memory
## names       :        slope,       aspect 
## min values  : 0.0009308692,   0.01101528 
## max values  : 4.4821872689, 359.92884285

4.1 capa 1: pendiente

(slope = subset(slp_asp, 1))
## class       : SpatRaster 
## size        : 113, 70, 1  (nrow, ncol, nlyr)
## resolution  : 1839.834, 1839.834  (x, y)
## extent      : 4702678, 4831466, 2472003, 2679904  (xmin, xmax, ymin, ymax)
## coord. ref. : MAGNA-SIRGAS 2018 / Origen-Nacional (EPSG:9377) 
## source(s)   : memory
## name        :        slope 
## min value   : 0.0009308692 
## max value   : 4.4821872689

la capa de la pendiente puede visualizarse a continuación:

terra::hist(slope, 
     main = "Pendiente SUCRE", 
     xlab = "Pendiente (en grados)")

4.2 capa 2: Orientación (aspecto)

(aspect =subset(slp_asp, 2))
## class       : SpatRaster 
## size        : 113, 70, 1  (nrow, ncol, nlyr)
## resolution  : 1839.834, 1839.834  (x, y)
## extent      : 4702678, 4831466, 2472003, 2679904  (xmin, xmax, ymin, ymax)
## coord. ref. : MAGNA-SIRGAS 2018 / Origen-Nacional (EPSG:9377) 
## source(s)   : memory
## name        :       aspect 
## min value   :   0.01101528 
## max value   : 359.92884285

la capa orientacion puede visualizarse a continuación:

terra::hist(aspect, 
     main = "Aspecto de SUCRE", 
     xlab = "Aspecto (en grados)")

se debe convertir la pendiente de grados a porcentaje:

# 
(slope_perc = tan(slope*(pi/180))*100)
## class       : SpatRaster 
## size        : 113, 70, 1  (nrow, ncol, nlyr)
## resolution  : 1839.834, 1839.834  (x, y)
## extent      : 4702678, 4831466, 2472003, 2679904  (xmin, xmax, ymin, ymax)
## coord. ref. : MAGNA-SIRGAS 2018 / Origen-Nacional (EPSG:9377) 
## source(s)   : memory
## name        :       slope 
## min value   : 0.001624673 
## max value   : 7.838889802

4.3 gráfica con la nueva conversión

terra::hist(slope_perc, 
     main = "pendiente SUCRE", 
     xlab = "Pendiente (en porcentaje)")

#terra::hist(slope_perc)

4.4 Cálculo de estadísticas zonales de Sucre

Inicialmente, para analizar las pendientes en el departamento, se utiliza la clasificación diseñada por el IGAC con fines agrológicos, que puede observarse a continuación:

# Reclassify the slope raster
#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)

4.5 Cálculos de las características zonales usando axactextractr:

(munic$mean_slope <- exactextractr::exact_extract(slope_perc, munic, 'mean'))
## Warning in .local(x, y, ...): Polygons transformed to raster CRS (EPSG:9377)
##   |                                                                              |                                                                      |   0%  |                                                                              |===                                                                   |   4%  |                                                                              |=====                                                                 |   8%  |                                                                              |========                                                              |  12%  |                                                                              |===========                                                           |  15%  |                                                                              |=============                                                         |  19%  |                                                                              |================                                                      |  23%  |                                                                              |===================                                                   |  27%  |                                                                              |======================                                                |  31%  |                                                                              |========================                                              |  35%  |                                                                              |===========================                                           |  38%  |                                                                              |==============================                                        |  42%  |                                                                              |================================                                      |  46%  |                                                                              |===================================                                   |  50%  |                                                                              |======================================                                |  54%  |                                                                              |========================================                              |  58%  |                                                                              |===========================================                           |  62%  |                                                                              |==============================================                        |  65%  |                                                                              |================================================                      |  69%  |                                                                              |===================================================                   |  73%  |                                                                              |======================================================                |  77%  |                                                                              |=========================================================             |  81%  |                                                                              |===========================================================           |  85%  |                                                                              |==============================================================        |  88%  |                                                                              |=================================================================     |  92%  |                                                                              |===================================================================   |  96%  |                                                                              |======================================================================| 100%
##  [1] 0.99318403 0.38188815 0.17414238 3.35606432 0.46806931 0.32033587
##  [7] 4.02443218 0.34662846 0.33519977 0.04102257 0.38929290 0.72070783
## [13] 0.02801267 0.92830676 1.41653001 0.42902663 0.54026157 0.10614514
## [19] 0.55798924 0.17493942 0.70880634 0.65105468 0.45322528 0.02788872
## [25] 0.24941583 1.30437350

visualización gráfica mediante histograma de la pendiente media en Sucre:

hist(munic$mean_slope, 
     main = "Pendiente media de SUCRE", 
     xlab = "Pendiente (en porcentaje)")

4.6 Pendiente predominante de los municipios de Sucre:

(munic$class <- exactextractr::exact_extract(rc, munic, 'mode'))
## Warning in .local(x, y, ...): Polygons transformed to raster CRS (EPSG:9377)
##   |                                                                              |                                                                      |   0%  |                                                                              |===                                                                   |   4%  |                                                                              |=====                                                                 |   8%  |                                                                              |========                                                              |  12%  |                                                                              |===========                                                           |  15%  |                                                                              |=============                                                         |  19%  |                                                                              |================                                                      |  23%  |                                                                              |===================                                                   |  27%  |                                                                              |======================                                                |  31%  |                                                                              |========================                                              |  35%  |                                                                              |===========================                                           |  38%  |                                                                              |==============================                                        |  42%  |                                                                              |================================                                      |  46%  |                                                                              |===================================                                   |  50%  |                                                                              |======================================                                |  54%  |                                                                              |========================================                              |  58%  |                                                                              |===========================================                           |  62%  |                                                                              |==============================================                        |  65%  |                                                                              |================================================                      |  69%  |                                                                              |===================================================                   |  73%  |                                                                              |======================================================                |  77%  |                                                                              |=========================================================             |  81%  |                                                                              |===========================================================           |  85%  |                                                                              |==============================================================        |  88%  |                                                                              |=================================================================     |  92%  |                                                                              |===================================================================   |  96%  |                                                                              |======================================================================| 100%
##  [1] 1 1 1 1 1 1 2 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
hist(munic$class, 
     main = "Pendiente reclasificada de SUCRE", 
     xlab = "Pendiente (y su categoria)")

transformación de la pendiente a coordenadas geográficas:

#terra::hist(slope_perc)
# This is the reclassified slope
(rc.geo = project(rc, "EPSG:4326"))
## class       : SpatRaster 
## size        : 113, 71, 1  (nrow, ncol, nlyr)
## resolution  : 0.01666703, 0.01666703  (x, y)
## extent      : -75.71422, -74.53086, 8.265791, 10.14917  (xmin, xmax, ymin, ymax)
## coord. ref. : lon/lat WGS 84 (EPSG:4326) 
## source(s)   : memory
## name        : slope 
## min value   :     1 
## max value   :     3

Luego de obtener las coordenadas geográficas, estas se convierten a porcentaje de pendiente (ya lo teníamos en slope_perc, solo se proyecta a WGS84 para visualización en leaflet).

# This is the percentage slope
(slope.geo = project(slope_perc, "EPSG:4326"))
## class       : SpatRaster 
## size        : 113, 71, 1  (nrow, ncol, nlyr)
## resolution  : 0.01666703, 0.01666703  (x, y)
## extent      : -75.71422, -74.53086, 8.265791, 10.14917  (xmin, xmax, ymin, ymax)
## coord. ref. : lon/lat WGS 84 (EPSG:4326) 
## source(s)   : memory
## name        :       slope 
## min value   : 0.003078969 
## max value   : 7.526045322

5 Representación gráfica de la pendiente del terreno

Al querer graficar los procedimientos hechos anteriormente, se genera un “mapa” interactivo, con la intención de visualizar la pendiente de los municipios en Sucre. Es bastante útil ya que permite comprender las características topográficas de mejor manera y generar un material esencial para tomar decisiones relacionadas con el uso de los suelos y la planificación de actividades en el territorio.

Para la configuración del mapa interactivo, inicialmente se definen los colores para los diferentes niveles de la pendiente, usando el más bajo como aguamarina, y el más alto como rojo oscuro. Se identifican los límites municipales con un borde de color gris, y la opacidad del relleno es de un 40%, además, al interactuar el mapa haciendo clic en algún municipio, se desplegará el cuadro emergente con el nombre del municipio y la clase de pendiente más frecuente.

# expand the number of colors to improve your visualization
# https://r-graph-gallery.com/42-colors-names.html
palredgreen <- colorNumeric(c("darkseagreen3","yellow2", "orange", "brown2", "darkred"), values(slope.geo),
  na.color = "transparent")
leaflet(munic) %>% addTiles() %>% setView(-75.9, 8.7, 9) %>% 
    addPolygons(color = "gray", weight = 1.0, smoothFactor = 0.5,
    opacity = 0.4, fillOpacity = 0.10,
    popup = paste("Municipio: ", munic$mpio_cnmbr, "<br>",
                          "Slope class: ", munic$class, "<br>")) %>%
  addRasterImage(slope.geo, colors = palredgreen, opacity = 0.8)  %>%
  addLegend(pal = palredgreen, values = values(slope.geo),
    title = "Pendiente del terreno en Sucre(%)")
## Warning: sf layer has inconsistent datum (+proj=longlat +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +no_defs).
## Need '+proj=longlat +datum=WGS84'

6 Visualización alternativa de la pendiente en grados

Se propone visualizar la pendiente del terreno en grados, junto con las etiquetas de cada municipio con el siguiente código, manteniendo las ventanas emergentes. Inicialmente se cambia de nuevo el sistema de coordenadas para que no se generen errores en la geometría.

library(sf)
library(terra)
library(leaflet)
print(st_crs(munic))
## Coordinate Reference System:
##   User input: MAGNA-SIRGAS 
##   wkt:
## GEOGCRS["MAGNA-SIRGAS",
##     DATUM["Marco Geocentrico Nacional de Referencia",
##         ELLIPSOID["GRS 1980",6378137,298.257222101,
##             LENGTHUNIT["metre",1]]],
##     PRIMEM["Greenwich",0,
##         ANGLEUNIT["degree",0.0174532925199433]],
##     CS[ellipsoidal,2],
##         AXIS["geodetic latitude (Lat)",north,
##             ORDER[1],
##             ANGLEUNIT["degree",0.0174532925199433]],
##         AXIS["geodetic longitude (Lon)",east,
##             ORDER[2],
##             ANGLEUNIT["degree",0.0174532925199433]],
##     USAGE[
##         SCOPE["Horizontal component of 3D system."],
##         AREA["Colombia - onshore and offshore. Includes San Andres y Providencia, Malpelo Islands, Roncador Bank, Serrana Bank and Serranilla Bank."],
##         BBOX[-4.23,-84.77,15.51,-66.87]],
##     ID["EPSG",4686]]
if (st_crs(munic)$epsg != 4326) {
  munic <- st_transform(munic, crs = 4326)
}
print(terra::crs(slope.geo))
## [1] "GEOGCRS[\"WGS 84\",\n    ENSEMBLE[\"World Geodetic System 1984 ensemble\",\n        MEMBER[\"World Geodetic System 1984 (Transit)\"],\n        MEMBER[\"World Geodetic System 1984 (G730)\"],\n        MEMBER[\"World Geodetic System 1984 (G873)\"],\n        MEMBER[\"World Geodetic System 1984 (G1150)\"],\n        MEMBER[\"World Geodetic System 1984 (G1674)\"],\n        MEMBER[\"World Geodetic System 1984 (G1762)\"],\n        MEMBER[\"World Geodetic System 1984 (G2139)\"],\n        MEMBER[\"World Geodetic System 1984 (G2296)\"],\n        ELLIPSOID[\"WGS 84\",6378137,298.257223563,\n            LENGTHUNIT[\"metre\",1]],\n        ENSEMBLEACCURACY[2.0]],\n    PRIMEM[\"Greenwich\",0,\n        ANGLEUNIT[\"degree\",0.0174532925199433]],\n    CS[ellipsoidal,2],\n        AXIS[\"geodetic latitude (Lat)\",north,\n            ORDER[1],\n            ANGLEUNIT[\"degree\",0.0174532925199433]],\n        AXIS[\"geodetic longitude (Lon)\",east,\n            ORDER[2],\n            ANGLEUNIT[\"degree\",0.0174532925199433]],\n    USAGE[\n        SCOPE[\"Horizontal component of 3D system.\"],\n        AREA[\"World.\"],\n        BBOX[-90,-180,90,180]],\n    ID[\"EPSG\",4326]]"
if (terra::crs(slope.geo, describe = TRUE)$code != 4326) { names(slope.geo) <- "slope"
 slope.geo <- terra::project(slope.geo, "EPSG:4326")
}

Calcular los centroides de los municipios luego de la transformación de las coordenadas:

centers <- st_centroid(munic)
## Warning: st_centroid assumes attributes are constant over geometries

7 Representación de la visualización alternativa de la pendiente en grados

Se propone el siguiente código para visualizar la pendiente en grados de forma alternativa, además, se agrega información adicional como el nombre de cada municipio y los valores medios de la pendiente.

palredgreen <- colorNumeric(palette = c("aquamarine","chartreuse","darkseagreen3","yellow2", "orange", "brown2","chocolate", "darkred"),
domain = values(slope.geo, na.rm = TRUE))
leaflet(munic) %>%
  addTiles() %>%
  setView(lng = -76.99, lat =  6.31, zoom = 7) %>%
  addPolygons(color = "white", weight = 1.0, smoothFactor = 0.5,
              opacity = 0.25, fillOpacity = 0.15,
              popup = paste("Municipio: ", munic$mpio_cnmbr, "<br>",
                            "Km2: ", munic$area, "<br>",
                            "Pendiente media de SUCRE (en grados): ", round(munic$class, 2), "<br>")) %>%
  addLabelOnlyMarkers(data = centers,
                      lng = ~st_coordinates(centers)[,1], lat = ~st_coordinates(centers)[,2],
                      label = ~paste(mpio_cnmbr, " - ", round(munic$class, 2), "°"),
                      labelOptions = labelOptions(noHide = TRUE, direction = 'top',
                                                  textOnly = TRUE, textsize = "10px")) %>%
  addRasterImage(slope.geo, colors = palredgreen, opacity = 0.9, project = FALSE) %>%
  addLegend(pal = palredgreen, values = values(slope.geo) ,
            title = "Pendiente SUCRE")

El gráfico anterior representa la pendiente del departamento en grados, usando una paleta de colores que se maneja desde aguamarina para valores bajos y rojo oscuro para los valores más altos. Con la intención de dar mayor información, se agregan las etiquetas de los municipios, el valor promedio de la pendiente en grados y un cuadro emergente con detalles adicionales. Este tipo de representaciones ofrecen información útil para la toma de decisiones dentro de los territorios, enfocadas al análisis de riesgos, usos recomendados o posibles del terreno, planificación en temas de infraestructura vial, por ejemplo, para que se adapte integralmente a las condiciones del terreno. Y por último, pero no menos importante, la planeación y selección de áreas para usos agrícolas y pecuarios.

8 Visualización del aspecto de Sucre

A continuación, se presenta un código que permite visualizar el aspecto de Sucre, usando un método plano, se aplican colores llamativos que permitan identificar las variaciones de las pendientes, todo esto, usando la librería raster.

library(raster)
## Cargando paquete requerido: sp
elev_sucre <- raster("ELEVACION.tif")
exten_sucre <- extent( -75.7, -74.53333, 8.283333, 10.14167)
elev_sucre_final <- crop(elev_sucre, exten_sucre)
elev_sucre_final_reducida <- aggregate(elev_sucre, fact=1)
## Warning in .local(x, ...): all fact(s) were 1, nothing to aggregate
aspecto_final <- terrain(elev_sucre_final_reducida, opt="aspect", unit= "degrees", neighbors=8)

#Paleta de colores
colors <- colorRampPalette(c("magenta","orange","red", "yellow", "green","blue"))(255)

leaflet() %>%
  addTiles() %>%
  setView(lng = -76.99, lat =  6.31, zoom = 7) %>%
  addRasterImage(aspecto_final, colors = colors, opacity = 0.8) %>%
  addLegend(pal = colorNumeric(c("magenta", "red", "yellow", "orange", "green", "blue"), values(aspecto_final)), values = values(aspecto_final), title = "Aspecto del terreno")

Como se observa anteriormente, en el departamento la mayor parte del terreno tiene un aspecto que oscila entre los (rangos observados en tu mapa). La variabilidad observada gracias a la paleta de colores, permite entender la topografía de Sucre.

9 Bibliografía

Lizarazo, I., 2025. Geomorphometric terrain attributes in R. Available at https://rpubs.com/ials2un/geomorphometric

sessionInfo()
## R version 4.5.1 (2025-06-13 ucrt)
## Platform: x86_64-w64-mingw32/x64
## Running under: Windows 11 x64 (build 26100)
## 
## Matrix products: default
##   LAPACK version 3.12.1
## 
## 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] raster_3.6-32        sp_2.2-0             exactextractr_0.10.0
## [4] MultiscaleDTM_1.0    terra_1.8-54         leaflet_2.2.2       
## [7] sf_1.0-21            elevatr_0.99.0      
## 
## loaded via a namespace (and not attached):
##  [1] s2_1.1.9           sass_0.4.10        generics_0.1.4     class_7.3-23      
##  [5] KernSmooth_2.23-26 lattice_0.22-7     digest_0.6.37      magrittr_2.0.3    
##  [9] rgl_1.3.24         RColorBrewer_1.1-3 evaluate_1.0.4     grid_4.5.1        
## [13] fastmap_1.2.0      jsonlite_2.0.0     e1071_1.7-16       DBI_1.2.3         
## [17] promises_1.3.3     purrr_1.1.0        scales_1.4.0       crosstalk_1.2.1   
## [21] codetools_0.2-20   jquerylib_0.1.4    cli_3.6.5          shiny_1.11.1      
## [25] rlang_1.1.6        units_0.8-7        base64enc_0.1-3    cachem_1.1.0      
## [29] yaml_2.3.10        tools_4.5.1        dplyr_1.1.4        httpuv_1.6.16     
## [33] png_0.1-8          vctrs_0.6.5        R6_2.6.1           mime_0.13         
## [37] proxy_0.4-27       lifecycle_1.0.4    classInt_0.4-11    htmlwidgets_1.6.4 
## [41] pkgconfig_2.0.3    progressr_0.15.1   bslib_0.9.0        pillar_1.11.0     
## [45] later_1.4.2        glue_1.8.0         Rcpp_1.1.0         tidyselect_1.2.1  
## [49] xfun_0.52          tibble_3.3.0       rstudioapi_0.17.1  knitr_1.50        
## [53] farver_2.1.2       xtable_1.8-4       htmltools_0.5.8.1  rmarkdown_2.29    
## [57] wk_0.9.4           compiler_4.5.1