Este notebook en R tiene como objetivo realizar un análisis espacial
de los municipios del departamento de Vichada, Colombia. Para ello, se
cargan datos geoespaciales en formato .gpkg que contienen
las características de los municipios y se complementan con datos de
elevación extraídos de un archivo raster. A lo largo del análisis, se
calculan variables geoespaciales como la pendiente y la orientación del
terreno, que son fundamentales para el análisis de la topografía de la
región. También se emplean herramientas de visualización interactivas,
como leaflet, para crear mapas interactivos que permiten
explorar los datos geoespaciales.
Aquí se instala el paquete terra, utilizado para
manipular datos raster, como los datos de elevación. Otros paquetes como
sf y leaflet se cargarán más adelante para
trabajar con datos vectoriales y crear mapas interactivos.
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.12.2, GDAL 3.9.3, PROJ 9.4.1; sf_use_s2() is TRUE
library(leaflet)
library(terra)
## terra 1.8.10
library(dplyr)
##
## Adjuntando el paquete: 'dplyr'
## The following objects are masked from 'package:terra':
##
## intersect, union
## The following objects are masked from 'package:stats':
##
## filter, lag
## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, union
library(ggplot2)
library(magrittr)
##
## Adjuntando el paquete: 'magrittr'
## The following objects are masked from 'package:terra':
##
## extract, inset
library(rmarkdown)
list.files("C:/Users/Danna Abril/Desktop/proyecto2/RStudio/Proyecto4/datos/mun_vichada.gpkg")
## character(0)
Aquí se muestra específicamente el archivo
mun_vichada.gpkg, que contiene la información geoespacial
de los municipios de Vichada.
El archivo .gpkg es leído y cargado en un objeto
sf que almacena los datos espaciales en formato
vectorial.
munic <- sf::st_read("C:/Users/Danna Abril/Desktop/proyecto2/RStudio/Proyecto4/datos/mun_vichada.gpkg")
## Multiple layers are present in data source C:\Users\Danna Abril\Desktop\proyecto2\RStudio\Proyecto4\datos\mun_vichada.gpkg, reading layer `municipios_colombia'.
## Use `st_layers' to list all layer names and their type in a data source.
## Set the `layer' argument in `st_read' to read a particular layer.
## Warning in CPL_read_ogr(dsn, layer, query, as.character(options), quiet, :
## automatically selected the first layer in a data source containing more than
## one.
## Reading layer `municipios_colombia' from data source
## `C:\Users\Danna Abril\Desktop\proyecto2\RStudio\Proyecto4\datos\mun_vichada.gpkg'
## using driver `GPKG'
## Simple feature collection with 4 features and 11 fields
## Geometry type: MULTIPOLYGON
## Dimension: XY
## Bounding box: xmin: -71.07793 ymin: 2.737109 xmax: -67.4098 ymax: 6.324317
## Geodetic CRS: MAGNA-SIRGAS
El archivo contiene las geometrías y los atributos de los municipios, que se utilizarán para calcular puntos centrales y elevarlos para análisis de elevación.
En este bloque, se calculan los centros geométricos (centroides) de los municipios para obtener una representación simplificada de su ubicación.
(centers <- st_centroid(munic))
## Warning: st_centroid assumes attributes are constant over geometries
Luego, se extraen las coordenadas (longitud y latitud) de los centroides:
centers$x = st_coordinates(centers)[,1]
centers$y = st_coordinates(centers)[,2]
Se usa el paquete elevatr para obtener datos de
elevación a partir de un raster de elevación.
elevation <- get_elev_raster(munic, z = 10)
## Mosaicing & Projecting
## Note: Elevation units are in meters.
El parámetro z = 10 establece la resolución de los datos
de elevación, que en este caso se elige como nivel 10 de resolución.
##5.1.Cargar librerías necesarias:
Aquí se cargan las librerías terra y rasterVis. La librería terra se usa para manejar y procesar datos espaciales raster (como imágenes geográficas, por ejemplo, mapas de elevación), mientras que rasterVis permite visualizar estos datos de manera más avanzada, como a través de mapas temáticos.
library(terra)
library(rasterVis)
## Cargando paquete requerido: lattice
Este bloque define un vector paquetes que contiene el nombre de las librerías necesarias para el análisis: terra (procesamiento de datos raster), elevatr (obtención de datos de elevación), sf (trabajo con objetos espaciales), y leaflet (para la visualización interactiva de mapas). Luego se cargan todas las librerías.
paquetes = c("terra","elevatr","sf", "leaflet")
library(elevatr)
library(sf)
library(leaflet)
library(terra)
Este comando obtiene el sistema de coordenadas espaciales (CRS, por sus siglas en inglés) del objeto munic (que contiene los municipios), lo cual es importante para garantizar que los datos estén correctamente alineados geográficamente.
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]]
En esta línea, se calcula el centroide de cada polígono en el objeto munic (los municipios), y el resultado se guarda en una nueva columna llamada centroid. Los centroides son los puntos centrales de los polígonos de los municipios y se utilizan para realizar análisis espaciales más detallados.
munic$centroid <- sf::st_centroid(munic$geom)
Se extraen las coordenadas de los centroides
coords <- st_coordinates(centers)
Aquí se extraen las coordenadas (x, y) de los centroides (almacenadas en coords) y se crean un nuevo dataframe llamado puntos. Este dataframe contendrá las ubicaciones específicas de los centroides de cada municipio, que más tarde serán utilizados para obtener los datos de elevación de cada punto.
puntos <- data.frame(x = coords[, 1], y = coords[, 2])
Utilizando la función get_elev_point() de la librería elevatr, se obtiene la elevación para cada uno de los puntos definidos en puntos. El parámetro locations toma las coordenadas de los centroides, prj asegura que las coordenadas se usen en el sistema de referencia correcto, y src = “aws”
elev_data <- get_elev_point(locations = puntos, prj = sf::st_crs(munic), src = "aws")
## Mosaicing & Projecting
## Warning: [extract] transforming vector data to the CRS of the raster
## Note: Elevation units are in meters
library(sf)
library(tidyr)
##
## Adjuntando el paquete: 'tidyr'
## The following object is masked from 'package:magrittr':
##
## extract
## The following object is masked from 'package:terra':
##
## extract
munic$elevacion <- elev_data$elevation
coords_elev <- st_coordinates(elev_data)
elev_data <- cbind(elev_data, coords_elev)
munic_clean <- na.omit(munic)
munic_clean <- munic %>%
drop_na(geom, elevacion)
munic_clean <- munic[sapply(munic$geom, function(x) !is.null(x)), ]
munic <- st_transform(munic, crs = 4326)
munic <- munic[!st_is_empty(munic$geom), ]
munic <- munic[st_is_valid(munic$geom), ]
munic <- munic[!st_is_empty(munic$geom) & st_is_valid(munic$geom), ]
Se crea un mapa interactivo usando leaflet para
visualizar los municipios de Vichada con un gradiente de color en
función de la elevación.
library(leaflet)
leaflet(munic) %>%
addTiles() %>%
addPolygons(fillColor = ~colorQuantile("YlOrRd", elevacion)(elevacion),
fillOpacity = 0.7, weight = 1) %>%
addMarkers(data = elev_data, ~X, ~Y, popup = ~paste("Elevación:", elevation, "m"))
library(leaflet)
map_interactivo <- leaflet(munic) %>%
addTiles() %>%
addPolygons(fillColor = ~colorQuantile("YlOrRd", elevacion)(elevacion),
fillOpacity = 0.7, weight = 1) %>%
addMarkers(data = elev_data, ~X, ~Y, popup = ~paste("Elevación:", elevation, "m"))
library(htmlwidgets)
saveWidget(map_interactivo, "map_interactivo.html", selfcontained = TRUE)
Este código genera un mapa con polígonos de los municipios de Vichada, coloreados según su elevación. También agrega marcadores con información sobre la elevación de cada municipio.
Aquí se genera un mapa interactivo de los centroides de los municipios de Vichada con un popup que muestra el nombre de cada municipio.
library(leaflet)
library(htmlwidgets)
# Crear el segundo mapa interactivo de los centroides
map_interactivo2 <- leaflet(data = centers) %>%
addTiles() %>%
addMarkers(
lng = ~x, # Longitud
lat = ~y, # Latitud
popup = ~mpio_cnmbr # Popup con el nombre del municipio
)
# Guardar el segundo mapa interactivo en un archivo HTML
saveWidget(map_interactivo2, "map_interactivo2.html", selfcontained = TRUE)
# Mostrar el mapa en la previsualización
map_interactivo2
Para crear una visualización estática más detallada, se utiliza
ggplot2 para generar un gráfico donde los municipios se
colorean según su elevación.
ggplot(munic) +
geom_sf(aes(fill = elevacion)) +
scale_fill_viridis_c() +
labs(title = "Elevación de los municipios del Vichada",
fill = "Elevación (m)") +
theme_minimal()
Este gráfico permite ver la distribución espacial de la elevación en el departamento de Vichada.
Aquí se carga el archivo raster de elevación
(elevancion_vichada.tif) y se realizan varias
transformaciones espaciales y de resolución.
dem = terra::rast("C:/Users/Danna Abril/Desktop/Proyecto 1/datos/elevancion_vichada.tif")
A continuación, se realiza un agregado del raster para reducir la resolución:
dem2 = terra::aggregate(dem, 2, "mean")
munic = sf::st_read("C:/Users/Danna Abril/Desktop/proyecto2/RStudio/Proyecto4/datos/mun_vichada.gpkg")
## Multiple layers are present in data source C:\Users\Danna Abril\Desktop\proyecto2\RStudio\Proyecto4\datos\mun_vichada.gpkg, reading layer `municipios_colombia'.
## Use `st_layers' to list all layer names and their type in a data source.
## Set the `layer' argument in `st_read' to read a particular layer.
## Warning in CPL_read_ogr(dsn, layer, query, as.character(options), quiet, :
## automatically selected the first layer in a data source containing more than
## one.
## Reading layer `municipios_colombia' from data source
## `C:\Users\Danna Abril\Desktop\proyecto2\RStudio\Proyecto4\datos\mun_vichada.gpkg'
## using driver `GPKG'
## Simple feature collection with 4 features and 11 fields
## Geometry type: MULTIPOLYGON
## Dimension: XY
## Bounding box: xmin: -71.07793 ymin: 2.737109 xmax: -67.4098 ymax: 6.324317
## Geodetic CRS: MAGNA-SIRGAS
print(munic)
## Simple feature collection with 4 features and 11 fields
## Geometry type: MULTIPOLYGON
## Dimension: XY
## Bounding box: xmin: -71.07793 ymin: 2.737109 xmax: -67.4098 ymax: 6.324317
## Geodetic CRS: MAGNA-SIRGAS
## dpto_ccdgo mpio_ccdgo mpio_cdpmp dpto_cnmbr mpio_cnmbr
## 1 99 001 99001 VICHADA PUERTO CARREÑO
## 2 99 524 99524 VICHADA LA PRIMAVERA
## 3 99 624 99624 VICHADA SANTA ROSALÍA
## 4 99 773 99773 VICHADA CUMARIBO
## mpio_crslc mpio_tipo mpio_narea mpio_nano
## 1 Decreto 1594 de Ago 5 de 1974 MUNICIPIO 12204.913 2023
## 2 Decreto 676 de Abril13 de 1987 MUNICIPIO 18569.338 2023
## 3 Ordenanza 19 de Noviembre 26 de 1993 MUNICIPIO 3691.869 2023
## 4 Ordenanza 66 de Noviembre 22 de 1996 MUNICIPIO 65597.251 2023
## shape_Leng shape_Area geom
## 1 5.475145 0.9859138 MULTIPOLYGON (((-67.80972 6...
## 2 8.080061 1.5062235 MULTIPOLYGON (((-69.03359 6...
## 3 3.805847 0.2999597 MULTIPOLYGON (((-70.65378 5...
## 4 18.794677 5.3085834 MULTIPOLYGON (((-68.47074 5...
El DEM es recortado para ajustarse a los límites de los municipios de Vichada:
dem3 = terra::crop(dem2, munic, mask = TRUE)
dem_plane = terra::project(dem3, "EPSG:9377")
print(dem_plane)
## class : SpatRaster
## dimensions : 216, 221, 1 (nrow, ncol, nlyr)
## resolution : 1849.053, 1849.053 (x, y)
## extent : 5212867, 5621507, 1859917, 2259312 (xmin, xmax, ymin, ymax)
## coord. ref. : MAGNA-SIRGAS 2018 / Origen-Nacional (EPSG:9377)
## source(s) : memory
## name : elevancion_vichada
## min value : 41.12041
## max value : 312.12048
munic_plane = sf::st_transform(munic, "EPSG:9377")
print(munic_plane)
## Simple feature collection with 4 features and 11 fields
## Geometry type: MULTIPOLYGON
## Dimension: XY
## Bounding box: xmin: 5213056 ymin: 1860707 xmax: 5619388 ymax: 2259698
## Projected CRS: MAGNA-SIRGAS 2018 / Origen-Nacional
## dpto_ccdgo mpio_ccdgo mpio_cdpmp dpto_cnmbr mpio_cnmbr
## 1 99 001 99001 VICHADA PUERTO CARREÑO
## 2 99 524 99524 VICHADA LA PRIMAVERA
## 3 99 624 99624 VICHADA SANTA ROSALÍA
## 4 99 773 99773 VICHADA CUMARIBO
## mpio_crslc mpio_tipo mpio_narea mpio_nano
## 1 Decreto 1594 de Ago 5 de 1974 MUNICIPIO 12204.913 2023
## 2 Decreto 676 de Abril13 de 1987 MUNICIPIO 18569.338 2023
## 3 Ordenanza 19 de Noviembre 26 de 1993 MUNICIPIO 3691.869 2023
## 4 Ordenanza 66 de Noviembre 22 de 1996 MUNICIPIO 65597.251 2023
## shape_Leng shape_Area geom
## 1 5.475145 0.9859138 MULTIPOLYGON (((5574599 225...
## 2 8.080061 1.5062235 MULTIPOLYGON (((5438952 224...
## 3 3.805847 0.2999597 MULTIPOLYGON (((5259904 215...
## 4 18.794677 5.3085834 MULTIPOLYGON (((5501962 217...
Se utilizan las funciones del paquete MultiscaleDTM para
calcular la pendiente y la orientación del terreno. Estos son dos
índices importantes para el análisis de la topografía.
slp_asp = MultiscaleDTM::SlpAsp(
dem_plane,
w = c(3, 3), # Tamaño de la ventana de análisis
unit = "degrees", # Unidad de medida en grados
method = "queen", # Método de cálculo
metrics = c("slope", "aspect"),
na.rm = TRUE,
include_scale = FALSE,
mask_aspect = TRUE
)
slope = subset(slp_asp, 1)
aspect = subset(slp_asp, 2)
slope_mean = exactextractr::exact_extract(slope, munic_plane, "mean")
## | | | 0% | |================== | 25% | |=================================== | 50% | |==================================================== | 75% | |======================================================================| 100%
print(slope_mean)
## [1] 0.12927212 0.12348688 0.09315646 0.17419241
slp_asp = MultiscaleDTM::SlpAsp(
dem_plane,
w = c(3, 3), # Tamaño de la ventana de análisis
unit = "degrees", # Unidad de medida en grados
method = "queen", # Método de cálculo
metrics = c("slope", "aspect"),
na.rm = TRUE,
include_scale = FALSE,
mask_aspect = TRUE
)
Se generan mapas de la pendiente y la orientación utilizando tanto
plot como leaflet para visualizaciones
interactivas.
plot(slope, main = "Mapa de Pendiente del Terreno de Vichada", col = terrain.colors(20))
Además, se genera un mapa interactivo de la pendiente:
leaflet() %>%
addTiles() %>%
addRasterImage(slope, colors = terrain.colors(20), opacity = 0.8) %>%
addLegend(pal = colorNumeric(terrain.colors(20), values(slope)), values = values(slope), title = "Pendiente de Vichada")
Finalmente, se genera el mapa de orientación del terreno.
library(leaflet)
custom_palette <- c("#ffffff", "#8B4513", "#DAA520", "#008000")
plot(aspect, main = "Orientación del Terreno", col = custom_palette, legend = TRUE)
Este código proporciona una completa visualización de los datos geoespaciales y una serie de análisis topográficos útiles para estudios de planificación territorial, manejo ambiental y desarrollo de infraestructuras en el departamento de Vichada.
en el siguiente chunk: Se carga la librería terra, que es útil para el procesamiento de datos raster. Se lista el contenido del paquete terra para verificar que se cargó correctamente, luego cargan dos archivos raster (elevancion_vichada.tif y elevacion_vichada.tif).calcula la pendiente (inclinación del terreno) en grados usando la función terrain(),extraen los valores del raster de pendiente con values(), y se usa ggplot2 para generar un histograma que muestra la distribución de la pendiente en grados.
# Cargar los paquetes necesarios
library(terra)
ls("package:terra")
## [1] "%in%" "activeCat" "activeCat<-"
## [4] "add_box" "add_grid" "add_legend"
## [7] "add_mtext" "add<-" "addCats"
## [10] "adjacent" "aggregate" "align"
## [13] "all.equal" "allNA" "animate"
## [16] "app" "approximate" "area"
## [19] "Arith" "as.array" "as.bool"
## [22] "as.contour" "as.data.frame" "as.factor"
## [25] "as.int" "as.lines" "as.list"
## [28] "as.matrix" "as.points" "as.polygons"
## [31] "as.raster" "atan_2" "atan2"
## [34] "autocor" "barplot" "bestMatch"
## [37] "blocks" "boundaries" "boxplot"
## [40] "buffer" "cartogram" "catalyze"
## [43] "categories" "cats" "cbind2"
## [46] "cellFromRowCol" "cellFromRowColCombine" "cellFromXY"
## [49] "cells" "cellSize" "centroids"
## [52] "clamp" "clamp_ts" "classify"
## [55] "clearance" "click" "colFromCell"
## [58] "colFromX" "colMeans" "colorize"
## [61] "colSums" "coltab" "coltab<-"
## [64] "combineGeoms" "combineLevels" "compare"
## [67] "Compare" "compareGeom" "concats"
## [70] "contour" "convHull" "costDist"
## [73] "countNA" "cover" "crds"
## [76] "crop" "crosstab" "crs"
## [79] "crs<-" "datatype" "deepcopy"
## [82] "delaunay" "densify" "density"
## [85] "depth" "depth<-" "describe"
## [88] "diff" "direction" "disagg"
## [91] "distance" "divide" "dots"
## [94] "draw" "droplevels" "elongate"
## [97] "emptyGeoms" "erase" "expanse"
## [100] "ext" "ext<-" "extend"
## [103] "extract" "extractAlong" "extractRange"
## [106] "fileBlocksize" "fillHoles" "fillTime"
## [109] "flip" "flowAccumulation" "focal"
## [112] "focal3D" "focalCpp" "focalMat"
## [115] "focalPairs" "focalReg" "focalValues"
## [118] "forceCCW" "free_RAM" "freq"
## [121] "gaps" "gdal" "gdalCache"
## [124] "geom" "geomtype" "getGDALconfig"
## [127] "getTileExtents" "global" "graticule"
## [130] "gridDist" "gridDistance" "halo"
## [133] "has.colors" "has.RGB" "has.time"
## [136] "hasMinMax" "hasValues" "head"
## [139] "hist" "hull" "identical"
## [142] "ifel" "image" "impose"
## [145] "inext" "init" "inMemory"
## [148] "inset" "interpIDW" "interpNear"
## [151] "interpolate" "intersect" "is.bool"
## [154] "is.empty" "is.factor" "is.int"
## [157] "is.lines" "is.lonlat" "is.points"
## [160] "is.polygons" "is.related" "is.rotated"
## [163] "is.valid" "isFALSE" "isTRUE"
## [166] "k_means" "lapp" "layerCor"
## [169] "levels" "libVersion" "linearUnits"
## [172] "lines" "logic" "Logic"
## [175] "longnames" "longnames<-" "makeNodes"
## [178] "makeTiles" "makeValid" "makeVRT"
## [181] "map.pal" "map_extent" "mask"
## [184] "match" "math" "Math"
## [187] "Math2" "mean" "median"
## [190] "mem_info" "merge" "mergeLines"
## [193] "mergeTime" "meta" "metags"
## [196] "metags<-" "minCircle" "minmax"
## [199] "minRect" "modal" "mosaic"
## [202] "na.omit" "NAflag" "NAflag<-"
## [205] "names" "ncell" "ncol"
## [208] "ncol<-" "nearby" "nearest"
## [211] "NIDP" "nlyr" "nlyr<-"
## [214] "noNA" "normalize.longitude" "north"
## [217] "not.na" "nrow" "nrow<-"
## [220] "nseg" "nsrc" "origin"
## [223] "origin<-" "pairs" "panel"
## [226] "patches" "perim" "persp"
## [229] "pitfinder" "plet" "plot"
## [232] "plotRGB" "points" "polys"
## [235] "prcomp" "predict" "princomp"
## [238] "project" "quantile" "query"
## [241] "rangeFill" "rapp" "rast"
## [244] "rasterize" "rasterizeGeom" "rasterizeWin"
## [247] "rcl" "readRDS" "readStart"
## [250] "readStop" "readValues" "rectify"
## [253] "regress" "relate" "removeDupNodes"
## [256] "res" "res<-" "resample"
## [259] "rescale" "rev" "RGB"
## [262] "RGB<-" "roll" "rotate"
## [265] "round" "rowColCombine" "rowColFromCell"
## [268] "rowFromCell" "rowFromY" "rowMeans"
## [271] "rowSums" "same.crs" "sapp"
## [274] "saveRDS" "sbar" "scale"
## [277] "scale_linear" "scoff" "scoff<-"
## [280] "sds" "segregate" "sel"
## [283] "selectHighest" "selectRange" "serialize"
## [286] "set.cats" "set.crs" "set.ext"
## [289] "set.names" "set.RGB" "set.values"
## [292] "setGDALconfig" "setMinMax" "setValues"
## [295] "shade" "sharedPaths" "shift"
## [298] "sieve" "simplifyGeom" "size"
## [301] "snap" "sort" "sources"
## [304] "spatSample" "spin" "split"
## [307] "sprc" "stdev" "stretch"
## [310] "subset" "subst" "summary"
## [313] "Summary" "surfArea" "svc"
## [316] "symdif" "t" "tail"
## [319] "tapp" "terrain" "terraOptions"
## [322] "text" "tighten" "time"
## [325] "time<-" "timeInfo" "tmpFiles"
## [328] "toMemory" "trans" "trim"
## [331] "union" "unique" "units"
## [334] "units<-" "unserialize" "unwrap"
## [337] "update" "values" "values<-"
## [340] "varnames" "varnames<-" "vect"
## [343] "vector_layers" "viewshed" "voronoi"
## [346] "vrt" "vrt_tiles" "watershed"
## [349] "weighted.mean" "where.max" "where.min"
## [352] "which.lyr" "which.max" "which.min"
## [355] "width" "window" "window<-"
## [358] "wrap" "wrapCache" "writeCDF"
## [361] "writeRaster" "writeStart" "writeStop"
## [364] "writeValues" "writeVector" "xapp"
## [367] "xFromCell" "xFromCol" "xmax"
## [370] "xmax<-" "xmin" "xmin<-"
## [373] "xres" "xyFromCell" "yFromCell"
## [376] "yFromRow" "ymax" "ymax<-"
## [379] "ymin" "ymin<-" "yres"
## [382] "zonal" "zoom"
dem <- rast("C:/Users/Danna Abril/Desktop/Proyecto 1/datos/elevancion_vichada.tif")
elevacion <- rast("C:/Users/Danna Abril/Desktop/Proyecto 1/datos/elevacion_vichada.tif")
slope_raster <- terrain(elevacion, v = "slope", unit = "degrees")
slope_values <- values(slope_raster)
library(ggplot2)
ggplot(data.frame(slope_values), aes(x = slope_values)) +
geom_histogram(binwidth = 1, fill = "blue", color = "black", alpha = 0.7) +
labs(title = "Distribución de la pendiente", x = "Pendiente (°)", y = "Frecuencia") +
theme_minimal()
## Warning: Removed 74616 rows containing non-finite outside the scale range
## (`stat_bin()`).
Se calcula el aspecto del terreno (orientación de la pendiente) y se genera un histograma para visualizar su distribución.
aspect_raster <- terrain(elevacion, v = "aspect")
aspect_values <- values(aspect_raster)
ggplot(data.frame(aspect_values), aes(x = aspect_values)) +
geom_histogram(binwidth = 10, fill = "green", color = "black", alpha = 0.7) +
labs(title = "Distribución del aspecto", x = "Aspecto (°)", y = "Frecuencia") +
theme_minimal()
## Warning: Removed 74616 rows containing non-finite outside the scale range
## (`stat_bin()`).
Se convierte la pendiente de grados a porcentaje utilizando la tangente
y se genera un histograma para analizar su distribución.
slope_raster <- terrain(elevacion, v = "slope", unit = "degrees")
slope_percent_raster <- tan(slope_raster * pi / 180) * 100
slope_percent_values <- values(slope_percent_raster)
ggplot(data.frame(slope_percent_values), aes(x = slope_percent_values)) +
geom_histogram(binwidth = 5, fill = "red", color = "black", alpha = 0.7) +
labs(title = "Distribución de la pendiente en porcentaje", x = "Pendiente (%)", y = "Frecuencia") +
theme_minimal()
## Warning: Removed 74616 rows containing non-finite outside the scale range
## (`stat_bin()`).
Se obtiene el valor medio de la pendiente y se visualiza en un histograma, marcando la media con una línea roja discontinua.
mean_slope <- global(slope_raster, fun = "mean")
slope_values <- values(slope_raster)
slope_values <- na.omit(slope_values)
slope_df <- data.frame(slope_values = slope_values)
mean_slope <- mean(slope_values)
ggplot(slope_df, aes(x = slope_values)) +
geom_histogram(binwidth = 1, fill = "blue", color = "black", alpha = 0.7) +
geom_vline(aes(xintercept = mean_slope), color = "red", linetype = "dashed", linewidth = 1) +
labs(title = "Distribución de la pendiente con pendiente media", x = "Pendiente (°)", y = "Frecuencia") +
theme_minimal()
Se agrupan los valores de pendiente en categorías definidas en una matriz de reclasificación. Luego, se genera un histograma para visualizar la frecuencia de cada categoría de pendiente.
rcl_matrix <- matrix(c(0, 5, 1, # valores entre 0 y 5 se recodifican como 1
5, 15, 2, # valores entre 5 y 15 se recodifican como 2
15, 30, 3, # valores entre 15 y 30 se recodifican como 3
30, 45, 4, # valores entre 30 y 45 se recodifican como 4
45, 90, 5, # valores entre 45 y 90 se recodifican como 5
90, 180, 6), # valores entre 90 y 180 se recodifican como 6
ncol = 3, byrow = TRUE)
slope_reclass <- classify(slope_raster, rcl = rcl_matrix)
slope_reclass_values <- values(slope_reclass)
ggplot(data.frame(slope_reclass_values), aes(x = slope_reclass_values)) +
geom_histogram(binwidth = 1, fill = "purple", color = "black", alpha = 0.7) +
labs(title = "Pendiente reclasificada", x = "Categorías de pendiente", y = "Frecuencia") +
theme_minimal()
## Warning: Removed 74616 rows containing non-finite outside the scale range
## (`stat_bin()`).
Este código realiza un análisis espacial del departamento de Vichada
utilizando R. Se cargan y procesan datos geoespaciales en formato
.gpkg y un archivo raster de elevación. A través de
herramientas como sf, terra, y
elevatr, se calculan variables geoespaciales clave, como la
pendiente y la orientación del terreno, y se visualizan en mapas
interactivos con leaflet y gráficos estáticos con
ggplot2. El análisis incluye la obtención de datos de
elevación, la creación de mapas temáticos, y el cálculo de indicadores
topográficos importantes para entender las características geográficas
de la región. Además, se utilizan técnicas de agregación y proyección de
datos raster para ajustarlos a la escala y forma de los municipios de
Vichada.
Lizarazo, I., 2025. Geomorphometric terrain attributes in R. Available at https://rpubs.com/ials2un/geomorphometric