El paquete sf de R fue diseñado para ofrecer soporte al
estándar simple features (sf) —o «objetos simples»—, utilizado
para la representación de datos geográficos vectoriales dentro del
ecosistema R. El paquete sf se lanzó por primera vez en
2016 y, en pocos años, superó con creces al paquete sp
(lanzado en 2005) como la opción más utilizada en R para representar
datos geográficos vectoriales. Además de brindar soporte al estándar
simple features, sf incluye una amplia gama de
funciones para trabajar con datos espaciales vectoriales en R. Si usted
está acostumbrado a trabajar con datos espaciales mediante software SIG
de escritorio tradicional, como QGIS o ArcGIS, encontrará en el paquete
sf —potenciado por los paquetes del tidyverse—
equivalentes a las herramientas SIG habituales para el análisis, la
consulta, la selección, el geoprocesamiento y la visualización.
«Simple features» (sf), o acceso a «simple features», hace referencia a un estándar formal (ISO 19125-1:2004) que describe cómo pueden representarse en sistemas informáticos los objetos del mundo real, poniendo énfasis en la geometría espacial de dichos objetos. Asimismo, describe cómo es posible almacenar estos objetos en bases de datos y recuperarlos de ellas, además de definir qué operaciones geométricas deben aplicarse a los mismos. Este estándar cuenta con una amplia implementación en bases de datos espaciales (como PostGIS) y sistemas de información geográfica (SIG) comerciales (por ejemplo, ESRI ArcGIS), y constituye la base de datos vectoriales para bibliotecas como GDAL.
Este cuaderno te presentará el paquete sf y el trabajo con datos vectoriales geográficos en el estándar de datos de objetos simples (sf). Este cuaderno solo pretende ser una introducción rápida al paquete sf, ya que muchas de sus funciones no se abordarán aquí.
Tener una copia del archivo shapefile MGN2023_MPIO_POLITICO disponible en DANE (y utilizado desde el inicio del curso de Geomática Básica).
https://www.jessesadler.com/post/simple-feature-objects/
https://learning.nceas.ucsb.edu/2023-04-coreR/session_15.html
Configuración
Lectura y escritura de datos espaciales
Creación y visualización de objetos espaciales
Este notebook debe ubicarse en tu directorio GB2, cerca del directorio donde guardaste tus datos.
Este notebook requiere los siguientes paquetes y funciones de R.
st_as_sf· convertir un objeto externo en un objeto sf
st_bbox · devolver el cuadro delimitador (bounding box) de una objetos simple o de un conjunto de objetos simples
st_geometry · obtener, establecer, reemplazar o renombrar la geometría de un objeto sf
st_read · leer entidades simples o capas desde un archivo o una base de datos
st_write · escribir un objeto de entidades simples en un archivo o una base de datos
https://dplyr.tidyverse.org/ - dplyr es una gramática para la manipulación de datos que ofrece un conjunto coherente de verbos, los cuales ayudan a resolver los desafíos más habituales en la manipulación de datos:
mutate() añade nuevas variables que son funciones de variables ya existentes.
select() selecciona variables basándose en sus nombres.
filter() selecciona casos basándose en sus valores.
summarise() reduce múltiples valores a un único resumen.
arrange() modifica el orden de las filas.
https://ggplot2.tidyverse.org/ ggplot2 es un sistema para crear gráficos de manera declarativa, basado en «La gramática de los gráficos» (The Grammar of Graphics). Tú proporcionas los datos, indicas a ggplot2 cómo asignar las variables a los elementos estéticos y qué primitivas gráficas utilizar; el sistema se encarga del resto de los detalles.
https://paleolimbot.github.io/ggspatial/ - Los datos espaciales, sumados a la potencia del marco de trabajo ggplot2, facilitan la creación de mapas.
library(sf)
## Linking to GEOS 3.14.1, GDAL 3.12.1, PROJ 9.7.1; sf_use_s2() is TRUE
library(dplyr)
##
## Adjuntando el paquete: 'dplyr'
## 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(ggspatial)
La función st_read, del paquete sf, lee
datos vectoriales espaciales a partir de una variedad de formatos
compatibles con GDAL (Geospatial Data Abstraction Library), incluidos
ESRI Shapefile (.shp), GeoJSON (.geojson), GeoPackage (.gpkg), KML
(.kml) y KMZ (.kmz), archivos CSV con WKT (well-known text) (.csv),
entre otros.
st_read simplifica el proceso de lectura de datos
espaciales al reconocer automáticamente el formato del archivo y
devolver los datos como un objeto de tipo «simple features» (sf).
#Get the filenames close to your notebook
list.files("D:/MI PC/GB02/Proyecto7/Mun_Cesar.gpkg")
## character(0)
Ahora puedes especificar la ruta relativa para leer el shapefile:
# read the *Municipios* shapefile as an sf object
Cesar <- st_read("D:/MI PC/GB02/Proyecto7/Mun_Cesar.gpkg")
## Multiple layers are present in data source D:\MI PC\GB02\Proyecto7\Mun_Cesar.gpkg, reading layer `Municipios'.
## 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' from data source
## `D:\MI PC\GB02\Proyecto7\Mun_Cesar.gpkg' using driver `GPKG'
## Simple feature collection with 24 features and 11 fields
## Geometry type: MULTIPOLYGON
## Dimension: XY
## Bounding box: xmin: -74.15633 ymin: 7.683902 xmax: -72.88663 ymax: 10.8524
## Geodetic CRS: WGS 84
Los datos se cargaron mediante el controlador de archivos shape de ESRI.
Los datos incluyen 24 objetos (el número de municipios en Cesar) y cuentan con 11 campos de atributos.
Cesar
## Simple feature collection with 24 features and 11 fields
## Geometry type: MULTIPOLYGON
## Dimension: XY
## Bounding box: xmin: -74.15633 ymin: 7.683902 xmax: -72.88663 ymax: 10.8524
## Geodetic CRS: WGS 84
## First 10 features:
## ID_0 ISO NAME_0 ID_1 NAME_1 ID_2 NAME_2 TYPE_2 ENGTYPE_2
## 1 53 COL Colombia 12 Cesar 445 Aguachica Municipio Municipality
## 2 53 COL Colombia 12 Cesar 446 Agustín Codazzi Municipio Municipality
## 3 53 COL Colombia 12 Cesar 447 Astrea Municipio Municipality
## 4 53 COL Colombia 12 Cesar 448 Becerril Municipio Municipality
## 5 53 COL Colombia 12 Cesar 449 Bosconia Municipio Municipality
## 6 53 COL Colombia 12 Cesar 450 Chimichagua Municipio Municipality
## 7 53 COL Colombia 12 Cesar 451 Chiriguaná Municipio Municipality
## 8 53 COL Colombia 12 Cesar 452 Curumaní Municipio Municipality
## 9 53 COL Colombia 12 Cesar 453 El Copey Municipio Municipality
## 10 53 COL Colombia 12 Cesar 454 El Paso Municipio Municipality
## NL_NAME_2 VARNAME_2 geom
## 1 <NA> <NA> MULTIPOLYGON (((-73.4802 8....
## 2 <NA> <NA> MULTIPOLYGON (((-72.95625 1...
## 3 <NA> <NA> MULTIPOLYGON (((-73.8072 9....
## 4 <NA> <NA> MULTIPOLYGON (((-72.95625 1...
## 5 <NA> <NA> MULTIPOLYGON (((-73.7618 9....
## 6 <NA> <NA> MULTIPOLYGON (((-73.6973 9....
## 7 <NA> <NA> MULTIPOLYGON (((-73.20842 9...
## 8 <NA> <NA> MULTIPOLYGON (((-73.38644 9...
## 9 <NA> <NA> MULTIPOLYGON (((-74.0645 9....
## 10 <NA> <NA> MULTIPOLYGON (((-73.8228 9....
Un MULTIPOLYGON (multipolígono) en el contexto de sistemas de información geográfica (SIG) y programación en R (por ejemplo, con el paquete sf) es una geometría vectorial que agrupa varios polígonos independientes dentro de un único objeto o registro.
st_bbox · devuelve el límite de un objeto simple o un conjunto de objetos simples.
MAGNA-SIRGAS utiliza el elipsoide GRS80 (Geodetic Reference System 1980), cuyos valores oficiales son:
• Semieje mayor (a): 6 378 137 m (Radio ecuatorial) • Semieje menor (b): 6 356 752,31414 m (Radio polar) • Aplanamiento (1/f): 298,257222101 • Unidad angular: Grados sexagesimales / Radianes • Meridiano de referencia: Greenwich (0°)
Mediante resoluciones oficiales del IGAC (como la Resolución 471 de 2020 y 370 de 2021), Colombia adoptó un Origen Único Nacional para la cartografía y el catastro multipropósito (identificado con el código internacional EPSG: 9377):
• Nombre del sistema: MAGNA-SIRGAS 2018 / Origen-Nacional (CTM12) • Tipo de proyección: Transversa de Mercator secante • Elipsoide: GRS 1980 • Falso Este (Coordenada X de origen): 5 000 000 m • Falso Norte (Coordenada Y de origen): 2 000 000 m • Meridiano central (λ₀): \(-73^\circ\,00^\prime\,00^{''} \,(73^\circ\text{ W})\) • Latitud de origen (φ₀): \(4^\circ\,00^\prime\,00^{''} \,(4^\circ\text{ N})\) • Factor de escala (k₀): 0,9992
head(Cesar)
## Simple feature collection with 6 features and 11 fields
## Geometry type: MULTIPOLYGON
## Dimension: XY
## Bounding box: xmin: -74.15633 ymin: 7.9975 xmax: -72.89888 ymax: 10.25736
## Geodetic CRS: WGS 84
## ID_0 ISO NAME_0 ID_1 NAME_1 ID_2 NAME_2 TYPE_2 ENGTYPE_2
## 1 53 COL Colombia 12 Cesar 445 Aguachica Municipio Municipality
## 2 53 COL Colombia 12 Cesar 446 Agustín Codazzi Municipio Municipality
## 3 53 COL Colombia 12 Cesar 447 Astrea Municipio Municipality
## 4 53 COL Colombia 12 Cesar 448 Becerril Municipio Municipality
## 5 53 COL Colombia 12 Cesar 449 Bosconia Municipio Municipality
## 6 53 COL Colombia 12 Cesar 450 Chimichagua Municipio Municipality
## NL_NAME_2 VARNAME_2 geom
## 1 <NA> <NA> MULTIPOLYGON (((-73.4802 8....
## 2 <NA> <NA> MULTIPOLYGON (((-72.95625 1...
## 3 <NA> <NA> MULTIPOLYGON (((-73.8072 9....
## 4 <NA> <NA> MULTIPOLYGON (((-72.95625 1...
## 5 <NA> <NA> MULTIPOLYGON (((-73.7618 9....
## 6 <NA> <NA> MULTIPOLYGON (((-73.6973 9....
Para exportar datos espaciales, la función st_write()
permite guardar un objeto sf en diversos formatos.
Puede especificar el formato de salida añadiendo la extensión de archivo correspondiente o indicando el controlador GDAL adecuado. Se pueden especificar opciones adicionales para controlar cómo se escriben los datos.
append = F permite sobrescribir
un archivo existente con el mismo nombre en el directorio de
trabajo.st_write(Cesar, "municipios.gpkg", driver = "GPKG", append = F)
## Deleting layer `municipios' using driver `GPKG'
## Writing layer `municipios' to data source `municipios.gpkg' using driver `GPKG'
## Writing 24 features with 11 fields and geometry type Multi Polygon.
# Check the new geopackage file was created
list.files(pattern="gpkg")
## [1] "cesar_munic.gpkg" "Mun_Cesar.gpkg" "municipios.gpkg" "stder_munic.gpkg"
# read the *Municipios* geopackage as an sf object
colombia2 <- st_read("./municipios.gpkg")
## Reading layer `municipios' from data source
## `D:\MI PC\GB02\Proyecto7\municipios.gpkg' using driver `GPKG'
## Simple feature collection with 24 features and 11 fields
## Geometry type: MULTIPOLYGON
## Dimension: XY
## Bounding box: xmin: -74.15633 ymin: 7.683902 xmax: -72.88663 ymax: 10.8524
## Geodetic CRS: WGS 84

La herramienta “Seleccionar por atributo” en el software SIG permite introducir una consulta SQL para seleccionar las entidades que cumplen los criterios de selección.
Por otro lado, sf y dplyr ofrecen opciones bastante sencillas para realizar operaciones similares. Analicemos una de ellas.
Comprueba el número de entidades seleccionadas. ¿Coincide con la información previa sobre tu departamento?
A continuación, vamos a crear un gráfico muy básico de los municipios de Cesar:
plot(st_geometry(Cesar), col = sf.colors(12, categorical = TRUE), border = 'grey', axes = TRUE, las = 1)
plot(st_geometry(st_centroid(Cesar)), pch = 3, col = 'red', add = TRUE)
## Warning: st_centroid assumes attributes are constant over geometries
Nótese que hemos trazado los límites de los municipios, así como sus centroides.
Ahora, podemos guardar las entidades seleccionadas en un nuevo archivo geopackage para su uso posterior:
st_write(Cesar, "cesar_munic.gpkg", driver = "GPKG", append = F)
## Deleting layer `cesar_munic' using driver `GPKG'
## Writing layer `cesar_munic' to data source `cesar_munic.gpkg' using driver `GPKG'
## Writing 24 features with 11 fields and geometry type Multi Polygon.
Una operación muy común al trabajar con geodatos es la selección de entidades en función de su relación espacial con otro objeto (capa). En el software SIG, este proceso se suele realizar mediante la herramienta Seleccionar por ubicación. Esta herramienta permite seleccionar entidades según las siguientes relaciones:
Intersect: Las entidades de la capa de entrada se seleccionarán si intersecan una entidad de selección. Esta es la opción predeterminada.
Within_a_distance: Las entidades de la capa de entrada se
seleccionarán si se encuentran dentro de la distancia especificada
(mediante la distancia euclidiana) de una entidad de selección. Utilice
el parámetro search_distance para especificar la
distancia.
Contains: Las entidades de la capa de entrada se seleccionarán si contienen una entidad de selección.
Completely_contains: Las entidades de la capa de entrada se seleccionarán si contienen completamente una entidad de selección.
Within: Las entidades de la capa de entrada se seleccionarán si se encuentran dentro de una entidad de selección.
Completely_within: Las entidades de la capa de entrada se seleccionarán si están completamente dentro o contenidas por una entidad de selección.
Are_identical_to: Las entidades de la capa de entrada se seleccionarán si son idénticas (en geometría) a una entidad de selección.
boundary_touches: Las entidades de la capa de entrada se seleccionarán si su límite toca una entidad de selección. Cuando las entidades de entrada son líneas o polígonos, el límite de la entidad de entrada solo puede tocar el límite de la entidad de selección, y ninguna parte de la entidad de entrada puede cruzar el límite de la entidad de selección.
share_a_line_segment_with: Las entidades de la capa de entrada se seleccionarán si comparten un segmento de línea con una entidad de selección. Las entidades de entrada y de selección deben ser líneas o polígonos.
crossed_by_the_outline_of: Las entidades de la capa de entrada se seleccionarán si son cruzadas por el contorno de una entidad de selección. Las entidades de entrada y de selección deben ser líneas o polígonos. Si se utilizan polígonos para la capa de entrada o de selección, se utilizará el límite (línea) del polígono. Se seleccionarán las líneas que se cruzan en un punto; las líneas que comparten un segmento de línea no se seleccionarán.
have_their_center_in: Las entidades de la capa de entrada se seleccionarán si su centro se encuentra dentro de una entidad de selección. El centro de la entidad se calcula de la siguiente manera: para polígonos y multipuntos, se utiliza el centroide de la geometría; Para la entrada de líneas, se utiliza el punto medio de la geometría.
Vamos a ilustrar la funcionalidad de selección por ubicación.
Primero, leamos un archivo CSV con ciudades del mundo. Sí, es el mismo archivo que usamos hace varias semanas en este curso. ¿Recuerdas el sistema de coordenadas de referencia (CRS) de este archivo?
cities = read.csv('D:/MI PC/GB02/Proyecto7/Cities.csv') %>%
st_as_sf(coords=c("lng","lat"), crs=4326) # remember x=longitude and y=latitude
cities
## Simple feature collection with 105 features and 7 fields
## Geometry type: POINT
## Dimension: XY
## Bounding box: xmin: -81.7006 ymin: -4.215 xmax: -67.4837 ymax: 12.5847
## Geodetic CRS: WGS 84
## First 10 features:
## city country iso2 admin_name capital population
## 1 Bogotá Colombia CO Bogotá primary 7968095
## 2 Cali Colombia CO Valle del Cauca admin 2916790
## 3 Medellín Colombia CO Antioquia admin 2490164
## 4 Barranquilla Colombia CO Atlántico admin 1326588
## 5 Cartagena Colombia CO Bolívar admin 914552
## 6 Bucaramanga Colombia CO Santander admin 893040
## 7 Palermo Colombia CO Huila minor 800000
## 8 Cúcuta Colombia CO Norte de Santander admin 748948
## 9 Soledad Colombia CO Atlántico minor 676014
## 10 Soacha Colombia CO Cundinamarca minor 655025
## population_proper geometry
## 1 8034649 POINT (-74.0817 4.6097)
## 2 2471474 POINT (-76.5197 3.44)
## 3 2529403 POINT (-75.5676 6.2502)
## 4 1326588 POINT (-74.8019 10.9833)
## 5 914552 POINT (-75.5144 10.3997)
## 6 614860 POINT (-73.1161 7.1186)
## 7 800000 POINT (-75.4339 2.8883)
## 8 806378 POINT (-72.5047 7.9075)
## 9 676014 POINT (-74.7858 10.9097)
## 10 782647 POINT (-74.2144 4.5781)
Es posible que hayamos notado que el sistema de referencia de coordenadas (CRS) de los dos conjuntos de datos (es decir, Santander y ciudades) no coincide. Vamos a demostrarlo:
st_crs(cities)$epsg
## [1] 4326
st_crs(Cesar)$epsg
## [1] 4326
Antes de utilizar cualquier selección por ubicación, necesitamos transformar el sistema de coordenadas de referencia (SCR) de las ciudades para que coincida con el SCR de Santander. Podemos hacerlo utilizando la función st_transform:
# take note of this function sintax
ncities <- st_transform(cities, crs= st_crs(Cesar))
Ahora, estamos listos para partir:
stder_cities <- ncities[Cesar, , op = st_within]
¿Obtuvimos lo que esperábamos?
Comprobémoslo de nuevo con un gráfico rudimentario:
plot(st_geometry(Cesar), col = sf.colors(12, categorical = TRUE), border = 'grey', axes = TRUE, las = 1)
plot(st_geometry(stder_cities), pch = 20, col = 'red', add = TRUE)
Nótese la variedad de símbolos que se pueden utilizar para representar geometrías de puntos:

Usaremos ggplot() y llamaremos a geom_sf(), que le indica a R que estamos trabajando con datos espaciales. Al crear un mapa, también es fundamental considerar el orden de las capas. Si colocamos los datos de las ciudades (stder_cities) antes del polígono de municipios de Cesar (Cesar), las ciudades quedarán debajo del polígono de municipios y no podremos verlas.
ggplot() +
#Add municipalities
geom_sf(data = Cesar) +
#Add cities layer
geom_sf(data = stder_cities, aes(color = city, label = city), size = 3) +
#Add titles
labs(x = "Longitud", y = "Latitud", title = "Ciudades de Cesar") +
#Add theme
theme_bw()
# 4.2 Representación gráfica mediante ggplot y ggspatial
Ahora, agregaremos algunos elementos personalizados a nuestro mapa, como una flecha que indique el norte y una barra de escala, utilizando el paquete ggspatial de R.
ggplot() +
#Crop Virginia boundary to spatial extent of cities and add Virginia layer
geom_sf(data = Cesar) +
#Add cities layer
geom_sf(data = stder_cities, aes(color = city, label = city), size = 3) +
#Add scale bar to bottom left from ggspatial
annotation_scale(location = "tr", height = unit(.25, "cm"),
width = unit(1, "cm"), pad_x = unit(0.3, "in"),
pad_y = unit(0.5, "in")) +
#Add north arrow to bottom left from ggspatial
annotation_north_arrow(height = unit(1, "cm"), width = unit(1, "cm"),
which_north = "true", location = "tr",
pad_x = unit(0.5, "in"), pad_y = unit(0.05, "in")) +
#Add titles
labs(x = "Longitud", y = "Latitud", title = "Ciudades del Cesar") +
#Add theme
theme_bw()
Si obtuviste resultados similares a los mostrados anteriormente,
¡felicidades!
sessionInfo()
## R version 4.6.1 (2026-06-24 ucrt)
## Platform: x86_64-w64-mingw32/x64
## Running under: Windows 11 x64 (build 22000)
##
## 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] ggspatial_1.1.11 ggplot2_4.0.3 dplyr_1.2.1 sf_1.1-3
##
## loaded via a namespace (and not attached):
## [1] gtable_0.3.6 jsonlite_2.0.0 compiler_4.6.1 tidyselect_1.2.1
## [5] Rcpp_1.1.2 jquerylib_0.1.4 scales_1.4.0 yaml_2.3.12
## [9] fastmap_1.2.0 R6_2.6.1 generics_0.1.4 classInt_0.4-11
## [13] s2_1.1.13 knitr_1.52 tibble_3.3.1 units_1.0-1
## [17] DBI_1.3.0 bslib_0.12.0 pillar_1.11.1 RColorBrewer_1.1-3
## [21] rlang_1.3.0 cachem_1.1.0 xfun_0.62 sass_0.4.10
## [25] S7_0.2.2 otel_0.2.0 cli_3.6.6 withr_3.0.3
## [29] magrittr_2.0.5 wk_0.9.5 class_7.3-23 digest_0.6.39
## [33] grid_4.6.1 rstudioapi_0.19.0 lifecycle_1.0.5 vctrs_0.7.3
## [37] KernSmooth_2.23-26 proxy_0.4-29 evaluate_1.0.5 glue_1.8.1
## [41] farver_2.1.2 e1071_1.7-17 rmarkdown_2.32 tools_4.6.1
## [45] pkgconfig_2.0.3 htmltools_0.5.9