Mapas
Análisis espacial moderno en R: vector, raster y visualización cartográfica
Sobre análisis espacial en R
El ecosistema espacial de R ha vivido un cambio generacional en los últimos años. La pila histórica (sp para vectores, rgdal y rgeos como puentes a GDAL/GEOS, raster para rejillas) fue retirada oficialmente en 2023: rgdal, rgeos y maptools están archivados en CRAN desde octubre de 2023. El stack moderno es relativamente compacto:
sf: datos vectoriales, basado en el estándar OGC Simple Features. Sustituye asp+rgdal+rgeos.terra: datos raster (y vector secundario), reescritura moderna del antiguoraster. Más rápido y con menor consumo de memoria gracias a su backend en C++.tmap,mapview,leaflet, capas de visualización: tematización rigurosa, exploración interactiva rápida y mapas web embebidos respectivamente.
Tres principios prácticos que conviene tener interiorizados antes de operar:
- CRS y EPSG son obligatorios, no opcionales. Cualquier objeto espacial sin sistema de referencia coordenado declarado es una bomba de tiempo.
sf::st_crs()yterra::crs()deben devolver siempre un EPSG conocido (EPSG:4326para lon/lat WGS84,EPSG:3857para web Mercator, los EPSG nacionales, 25830 / 25831 en España, para análisis métricos). - Vector vs raster son dos modelos de datos distintos y requieren paquetes distintos.
sfno hace raster,terrahace ambos pero su API vectorial es menos idiomática que la desf. La interoperabilidad es directa:sf::st_as_sf(terra_vect)yterra::vect(sf_object). - Orden de coordenadas. GDAL ≥ 3 respeta el orden de eje declarado por la autoridad (lat/lon para EPSG:4326 según la EPSG database).
sfyterranormalizan internamente a lon/lat para mantener la convención GIS, pero al leer KML, GML o servicios WFS desactualizados aparecen ejes invertidos con regularidad. Verifica siempre el primer plot.
Esta página cataloga seis paquetes que cubren el flujo completo, desde I/O hasta cartografía publicable. El orden refleja jerarquía conceptual: primero los datos (vector → raster), después las capas de visualización (estática → exploratoria → web interactiva), y al final los paquetes legacy que aún aparecen en código antiguo.
Instalación canónica del núcleo moderno:
install.packages(c("sf", "terra", "tmap", "mapview", "leaflet"))
# Verificar que sf detecta las librerías de sistema correctamente
sf::sf_extSoftVersion()
# GEOS, GDAL, PROJ — comprueba que las tres salen con versión y no NAsf
sf (Simple Features for R) es el paquete de referencia para datos vectoriales en R. Implementa el estándar OGC Simple Features (puntos, líneas, polígonos y sus colecciones) sobre data.frame, de modo que cada fila es una geometría con sus atributos asociados, directamente compatible con el verbo tidyverse y con el resto de la pila de R moderno.
Desarrollado por Edzer Pebesma. Es el sustituto canónico de sp + rgdal + rgeos, todos ya archivados. Cualquier proyecto vectorial nuevo debe empezar aquí.
Cuándo usarlo
- Lectura y escritura de cualquier formato vectorial (Shapefile, GeoPackage, GeoJSON, KML, GML, FlatGeobuf, PostGIS…). El driver lo decide GDAL automáticamente por extensión.
- Operaciones geométricas: intersección, unión, diferencia, buffer, distancia, contención,
validchecks. - Reproyección entre CRS (
st_transform()). - Encadenamiento con
dplyryggplot2(geom_sf()) sin conversiones intermedias.
Cuándo NO usarlo
- Raster: no es su dominio. Usa
terra. - Datasets vectoriales muy grandes (decenas de millones de geometrías): considera
sfcon backendarrow/geoparquet, o motores externos (PostGIS, DuckDB con la extensiónspatial, GeoParquet sobre Apache Arrow).sfcarga todo a memoria. - Operaciones esféricas exigentes sobre grandes mallas globales,
s2(quesfusa por defecto en lon/lat) cubre la mayoría de casos, pero para geodesia avanzadalwgeomo motores especializados son preferibles.
Conceptos clave
- Un objeto
sfes undata.framecon una columna especial de clasesfc(geometry column). Los métodosdplyr(filter,mutate,select, joins) funcionan transparentemente. - CRS siempre explícito:
st_crs(x)devuelve el sistema.st_set_crs(x, 4326)asigna uno sin reproyectar.st_transform(x, 25830)reproyecta los datos. - Geometría esférica por defecto en lon/lat. Desde sf 1.0, las operaciones sobre
EPSG:4326usan el motors2de Google (geodésico). Para forzar el comportamiento planar plano:sf_use_s2(FALSE). st_is_valid()/st_make_valid(), limpieza obligatoria antes de operaciones booleanas. Polígonos con auto-intersecciones rompen GEOS.- Las geometrías se serializan en formato WKB internamente. Para inspección legible,
st_as_text()produce WKT.
Patrón mínimo
library(sf)
library(dplyr)
# Lectura — el driver se infiere por extensión
provincias <- st_read("provincias_es.gpkg", layer = "provincias")
# CRS y reproyección a UTM 30N (península) para análisis métricos
st_crs(provincias)$epsg
provincias_utm <- st_transform(provincias, 25830)
# Operaciones encadenadas con dplyr
andalucia <- provincias_utm |>
filter(ccaa == "Andalucia") |>
mutate(area_km2 = as.numeric(st_area(geom)) / 1e6)
# Buffer + intersección
buffer_costa <- andalucia |>
filter(provincia %in% c("Cadiz", "Malaga", "Almeria")) |>
st_union() |>
st_buffer(dist = 5000) # 5 km en metros (porque estamos en UTM)
# Escritura
st_write(andalucia, "andalucia.gpkg", delete_layer = TRUE)Trampas habituales
- Buffer en grados. Si llamas a
st_buffer(x, 5000)sobre datos enEPSG:4326, el “5000” se interpreta como grados, no metros, produce un polígono planetario absurdo. Reproyecta primero a un CRS métrico apropiado (UTM, ETRS89-LAEA, proyección nacional). - GeoJSON con CRS distinto a 4326. El estándar GeoJSON (RFC 7946) obliga a usar
EPSG:4326. Algunos productores escriben otros CRS en GeoJSON. Consumidores estrictos los rechazan. Para coordenadas proyectadas, usa GeoPackage o FlatGeobuf. - Shapefile y nombres de columna largos. ESRI Shapefile trunca nombres a 10 caracteres. Si el GPKG es opción, prefiérelo siempre, sin límites, soporta múltiples capas y conserva tipos.
sf_use_s2(FALSE)como solución mágica. Suele aparecer en respuestas de Stack Overflow para “arreglar” errores de validez en lon/lat. Funciona, pero al desactivars2cualquier operación geodésica deja de ser geodésica. Mejor:st_make_valid()antes de operar.- Joins espaciales lentos.
st_joinsobre objetos grandes sin índice espacial puede ser O(n·m).sfconstruye índice automáticamente para predicados, pero si tienes problemas de rendimiento, considerasf::st_join(..., largest = TRUE)o un backend tipo DuckDB-spatial.
Enlaces
Relacionados en esta página
terra, contraparte raster (y vector secundario, menos idiomático).sp, predecesor archivado en CRAN.
terra
terra es el paquete moderno de análisis raster en R, y la reescritura oficial de raster por su mismo autor, Robert J. Hijmans. Más rápido, con menor huella de memoria gracias a procesamiento por bloques en C++, y con una API más coherente. Maneja además datos vectoriales como ciudadanos de segunda clase (clase SpatVector), útil cuando se mezclan ambos modelos en un pipeline.
Es el paquete a utilizar para imágenes satélite, modelos digitales de elevación (DEM), climatología, modelos de distribución de especies, índices espectrales (NDVI y derivados) y cualquier álgebra raster a escala.
Cuándo usarlo
- Análisis raster: álgebra de mapas, focal / zonal stats, reclasificación, reproyección.
- Manejo de stacks multi-banda (Sentinel-2, Landsat, MODIS) con
SpatRastermulti-layer. - Conversión entre modelos:
rasterize(),vectorize()/as.polygons(). - Cuando el dataset no cabe en memoria:
terraprocesa por bloques transparentemente.
Cuándo NO usarlo
- Vectores como flujo principal: la API
SpatVectores funcional pero menos idiomática quesf. Para análisis primariamente vectorial, usasfy convierte convect()solo al cruzar con raster. - Datacubes multidimensionales (tiempo, banda, sensor) con metadatos ricos:
stars(también de Edzer Pebesma) está mejor diseñado para datos científicos con dimensiones arbitrarias (xarray-like).terraprivilegia el caso GIS clásico. - Procesamiento en la nube sin descarga: para acceso COG (Cloud-Optimized GeoTIFF) o STAC sobre S3 / GCS,
terrafunciona perogdalcubes+rstacofrecen flujos más fluidos.
Conceptos clave
SpatRasterrepresenta rasters de una o más capas. El dato no necesita caber en memoria (se referencia el archivo en disco).SpatVectores la clase vectorial nativa. Convertir asfconst_as_sf()cuando interese unirse a la API vectorial moderna.- Álgebra de mapas vectorizada:
r1 + r2,log(r),ifel(r > 0, 1, NA)operan layer-by-layer. focal()para vecindarios (móvil),zonal()para estadísticas por zona vectorial,extract()para muestrear raster en geometrías.- Reproyección:
project()(noprojectRaster()como en el viejoraster). Eligemethod = "bilinear"para variables continuas,"near"para categóricas. - Memoria:
terra::terraOptions(memfrac = 0.6, tempdir = "...")controla el ratio de memoria y la carpeta de temporales, crítico en máquinas con poco RAM.
Patrón mínimo
library(terra)
# Lectura — un GeoTIFF de Sentinel-2 (4 bandas: B2 B3 B4 B8)
s2 <- rast("S2_T30STG_20240615.tif")
nlyr(s2)
crs(s2, describe = TRUE)
# NDVI (B8 NIR, B4 Red) — álgebra vectorizada
ndvi <- (s2[["B8"]] - s2[["B4"]]) / (s2[["B8"]] + s2[["B4"]])
names(ndvi) <- "NDVI"
# Estadísticas zonales sobre polígonos
parcelas <- vect("parcelas.gpkg")
ndvi_por_parcela <- extract(ndvi, parcelas, fun = mean, na.rm = TRUE)
# Reproyectar a UTM con remuestreo bilineal
ndvi_utm <- project(ndvi, "EPSG:25830", method = "bilinear")
# Escritura con compresión LZW y BigTIFF si procede
writeRaster(ndvi_utm, "ndvi.tif",
gdal = c("COMPRESS=LZW", "BIGTIFF=IF_SAFER"),
overwrite = TRUE)Trampas habituales
- Confundir
terrayraster. Son paquetes distintos con APIs incompatibles. Ejemplo:raster::projectRaster()no existe enterra, se llamaproject(). Si copias código de tutoriales antiguos, vas a chocar. extract()lenta. Para muchas geometrías pequeñas sobre raster grande,extract()puede ser lenta. Opciones:exact_extract(paqueteexactextractr) suele ser dramáticamente más rápida y produce ponderaciones por píxel.- Comparación de rasters con extensión / resolución distintas.
compareGeom()falla antes de fallar peor. Antes der1 + r2, asegúrate de que comparten CRS, extensión y resolución, o usaresample(). - NoData implícito. Algunos GeoTIFF heredados no tienen
NoDatadeclarado y guardan-9999como valor centinela literal. Eso contamina cualquier estadística. Compruebaminmax(r)y, si procede,r[r == -9999] <- NAoNAflag(r) <- -9999. writeRastersin compresión produce ficheros enormes. Activa siempreCOMPRESS=LZW(sin pérdida) oCOMPRESS=DEFLATEconPREDICTOR=2para enteros.
Enlaces
Relacionados en esta página
sf, contraparte vectorial. Conversiones en ambas direcciones son directas.raster, predecesor en mantenimiento mínimo, evitar para código nuevo.
tmap
tmap es el paquete de cartografía temática en R: una gramática inspirada en ggplot2 pero orientada a mapas estáticos publicables (capas con tm_polygons, tm_dots, tm_raster) con una clase de leyenda, escala, norte, breaks discretos / continuos y paletas pensadas para datos espaciales.
Su valor diferencial frente a ggplot2 + geom_sf es la calidad por defecto: clasificaciones de Jenks o cuantiles bien implementadas, gestión limpia de leyendas continuas, soporte raster nativo, y un switch trivial entre modo estático (tmap_mode("plot")) e interactivo (tmap_mode("view"), que delega en leaflet).
Cuándo usarlo
- Mapas para informes y publicaciones: cuando necesitas leyendas correctas, control fino sobre
breaks, north arrow, scale bar, múltiples facets por variable y composiciones lado a lado. - Coropletas y mapas de calor con clasificaciones estadísticas (Jenks, cuantiles, intervalos iguales, equal-area).
- Prototipado rápido con
tm_shape() + tm_dots()y luego refinamiento.
Cuándo NO usarlo
- Exploración rápida durante el análisis: el setup mínimo de
tmapes más verboso quemapview::mapview(x). Para “ver qué tengo”,mapviewes más eficiente. - Mapas web altamente customizados (popups con HTML rico, control de capas avanzado, tiles personalizados): mejor
leafletdirectamente. - Cuando
ggplot2 + geom_sfya basta. Si el resto de tu informe está enggplot2, mantener coherencia visual puede pesar más que las ventajas cartográficas detmap.
Conceptos clave
- Pipeline:
tm_shape(datos) + tm_polygons("variable", style = "jenks", n = 5, palette = "viridis"). stylecontrola la clasificación:"pretty","jenks","quantile","equal","cont"(continuo).tm_facets(by = "variable")produce small multiples alineados.tmap_mode("view")convierte cualquier mapa estático enleafletinteractivo sin cambiar el código de capas.tmap_mode("plot")vuelve a estático.tmap_save()exporta a PNG, PDF, SVG o HTML. Respeta DPI y dimensiones para publicación.- tmap 4 (publicado en 2024) reescribió la API interna sobre
ggplot2. Compatible hacia atrás en lo esencial, pero algunos parámetros menores han cambiado, revisa el changelog si tu código es de tmap 3.
Patrón mínimo
library(tmap)
library(sf)
provincias <- st_read("provincias_es.gpkg") |>
dplyr::mutate(densidad = poblacion / as.numeric(st_area(geom)) * 1e6)
# Estático: coropleta con clasificación Jenks
tmap_mode("plot")
mapa <- tm_shape(provincias) +
tm_polygons("densidad",
style = "jenks", n = 5,
palette = "YlOrRd",
title = "hab / km²") +
tm_compass(position = c("right", "top"), type = "8star", size = 2) +
tm_scale_bar(position = c("left", "bottom")) +
tm_layout(legend.outside = TRUE, frame = FALSE)
mapa
# Mismo objeto, exportado a PNG 300 DPI
tmap_save(mapa, "densidad_es.png", width = 18, height = 12, units = "cm", dpi = 300)
# Mismo objeto, ahora interactivo (leaflet por debajo)
tmap_mode("view")
mapaTrampas habituales
style = "pretty"por defecto suele producir clasificaciones poco discriminantes. Para datos con sesgo (densidades de población, ingresos),"jenks"o"quantile"casi siempre comunican mejor.- Coropleta sobre datos absolutos (no normalizados por área o población): el mapa termina mostrando “dónde hay más gente”, no la variable de interés. Normaliza antes de mapear.
- Leyendas cortadas cuando el frame es ajustado:
tm_layout(legend.outside = TRUE)resuelve el 80 % de los casos. - Cambio de API entre tmap 3 y 4. Si actualizas el paquete y el código rompe, suele ser un argumento renombrado o un valor por defecto cambiado. La migración suele ser directa.
Enlaces
Relacionados en esta página
mapview, alternativa para exploración rápida. Menos control cartográfico.leaflet, backend interactivo que tmap utiliza en modo"view".
mapview
mapview resuelve un problema cotidiano y específico: inspeccionar visualmente un objeto espacial en una sola línea. mapview(x) abre un mapa interactivo de Leaflet con OpenStreetMap de fondo, popups automáticos por atributos y soporte para vectores (sf) y rasters (terra, raster).
Es la herramienta que sustituye al plot() rápido durante la fase exploratoria: cuando lees una capa nueva y quieres comprobar dónde está, qué resolución tiene, qué columnas trae y si la geometría es coherente.
Cuándo usarlo
- Inspección durante el análisis: ver una capa recién leída, comparar dos capas superpuestas, comprobar resultado de un join espacial.
- Visualización rápida de rasters pequeños / medianos con paletas razonables por defecto.
- Compartir un snapshot puntual:
mapshot()exporta a PNG o HTML.
Cuándo NO usarlo
- Mapas para informe o publicación:
mapviewno está pensado para cartografía rigurosa. Para eso,tmapoleafletdirectamente. - Datasets grandes: por encima de unos miles de geometrías, el renderizado en Leaflet se vuelve pesado. Considera
leafgl(WebGL) o filtrar antes. - Customización fuerte de popups / capas: la API es deliberadamente minimal. Si necesitas control fino,
leafletdirecto es más productivo que pelearse con argumentos demapview.
Conceptos clave
mapview(x)es la única función que vas a usar el 90 % del tiempo.- Combinar capas con
+:mapview(provincias) + mapview(parcelas, color = "red"). zcol = "variable"colorea por una columna. Respeta el tipo (continuo / categórico).mapshot(m, file = "mapa.png")omapshot2()(basado enwebshot2+ Chromium) para exportar, útil cuando solo se necesita una imagen.
Patrón mínimo
library(mapview)
library(sf)
provincias <- st_read("provincias_es.gpkg")
parcelas <- st_read("parcelas.gpkg")
# Inspección de una sola capa
mapview(provincias)
# Coloreado por atributo
mapview(provincias, zcol = "densidad")
# Superposición de capas
mapview(provincias, alpha.regions = 0.3) +
mapview(parcelas, color = "red", col.regions = "red")
# Snapshot a fichero
m <- mapview(provincias, zcol = "densidad")
mapshot(m, file = "snapshot_provincias.png")Trampas habituales
- Renderizado lento con muchas features. Por encima de ~10.000 polígonos el navegador empieza a sufrir. Filtra (
dplyr::filter), simplifica (sf::st_simplify(dTolerance = ...)) o muévete aleafgl. - CRS no lon/lat reproyectados en silencio.
mapviewreproyecta aEPSG:4326automáticamente para Leaflet. Si tu CRS está mal declarado, las geometrías aterrizan en mitad del Atlántico. mapshot2requiere Chromium. Si está sin instalar,mapshot2falla.mapshotoriginal usaphantomjs(también externo). Considerawebshot2::install_chromium()en CI.
Enlaces
Relacionados en esta página
leaflet, backend quemapviewenvuelve. Para control fino, ir directo.tmap, alternativa cuando el output final no es exploratorio.
leaflet
leaflet es el binding de R para la librería JavaScript Leaflet, el estándar de facto para mapas web interactivos ligeros. Permite construir mapas completamente customizables, tiles, marcadores, polígonos, popups HTML, controles de capas, eventos, y embeberlos en informes Quarto / R Markdown, dashboards shiny o páginas estáticas.
Es la base sobre la que se construyen mapview y el modo "view" de tmap. Cuando esos paquetes se quedan cortos en customización, bajar a leaflet directamente es lo natural.
Cuándo usarlo
- Mapas web interactivos para dashboards, informes HTML y aplicaciones
shiny. - Customización de tiles (Carto, Stamen, Mapbox, IGN), popups HTML ricos, control de capas con
addLayersControl(). - Integración con
crosstalkoshinypara enlazar el mapa a filtros y otros widgets.
Cuándo NO usarlo
- Mapas estáticos para informe impreso o PDF: usa
tmapoggplot2 + geom_sf. Exportar leaflet a PNG es posible (mapview::mapshot) pero subóptimo. - Datasets muy grandes: Leaflet vanilla pinta cada feature en DOM. Para decenas de miles de puntos,
leafgl(mismo paquete extendido con WebGL) odeck.glson mejor elección. - Visualizaciones 3D / globo: Leaflet es estrictamente 2D Web Mercator. Para 3D, considera
cesium,mapdeckorayshader(este último orientado a render, no web).
Conceptos clave
- Pipeline en piping explícito:
leaflet() |> addTiles() |> addPolygons(data = x) |> addCircleMarkers(...). - Reproyección automática a
EPSG:3857(Web Mercator). Internamente Leaflet espera lon/lat (EPSG:4326) como entrada y proyecta. Cualquier otro CRS de entrada debe ir víasf::st_transform(x, 4326)primero. - Tiles:
addTiles()usa OSM por defecto.addProviderTiles(providers$CartoDB.Positron)para fondos limpios. Otros providers enleaflet::providers. popup = ~paste0("<b>", nombre, "</b><br>", densidad)admite HTML, útil para tarjetas informativas.- Para Shiny:
leafletOutput+renderLeaflet+leafletProxy()para actualizaciones reactivas sin redibujar el mapa entero.
Patrón mínimo
library(leaflet)
library(sf)
provincias <- st_read("provincias_es.gpkg") |>
st_transform(4326) # leaflet exige lon/lat
pal <- colorNumeric("YlOrRd", domain = provincias$densidad)
leaflet(provincias) |>
addProviderTiles(providers$CartoDB.Positron) |>
addPolygons(
fillColor = ~pal(densidad),
fillOpacity = 0.7,
color = "white",
weight = 1,
popup = ~paste0("<b>", provincia, "</b><br>",
round(densidad), " hab/km²"),
highlightOptions = highlightOptions(weight = 3, color = "#666",
bringToFront = TRUE)
) |>
addLegend(pal = pal, values = ~densidad,
title = "Densidad", position = "bottomright")Trampas habituales
- CRS distinto a
EPSG:4326. Es la fuente número uno de mapas en blanco o features en mitad del océano.st_transform(x, 4326)siempre antes deaddPolygons/addCircles. addTiles()solo, sinsetView(), en un mapa vacío produce un globo terráqueo centrado en (0,0). Si las capas se añaden después conleafletProxy, fijasetView(lng, lat, zoom)al principio.- Popups con encoding incorrecto. Caracteres no ASCII (acentos, eñes) requieren que la fuente de los datos esté en UTF-8. En Windows,
read.csvcon encoding por defecto suele dar problemas,readr::read_csv()odata.table::fread(..., encoding = "UTF-8")los evitan. - Tiles de pago sin token. Mapbox, MapTiler y otros providers requieren API key. Cuando aparece un mapa gris, suele ser un 401 silencioso en la petición de tiles. Comprueba en la consola del navegador.
Enlaces
Relacionados en esta página
sp (legacy)
sp fue durante más de una década la clase canónica de datos espaciales en R: SpatialPoints, SpatialPolygonsDataFrame, SpatialLines y compañía. Su ecosistema dependiente (rgdal para I/O, rgeos para geometría, maptools para utilidades) era el stack obligatorio de cualquier análisis espacial hasta ~2018.
rgdal, rgeos y maptools están archivados en CRAN desde octubre de 2023. sp sigue mantenido a mínimos para garantizar compatibilidad hacia atrás, pero su autor (Edzer Pebesma) recomienda explícitamente migrar a sf para todo proyecto nuevo.
Cuándo usarlo
Esencialmente nunca para código nuevo. Solo aparece en tres contextos:
- Mantenimiento de código heredado que se construyó sobre
spy aún no se ha migrado. - Paquetes downstream que todavía esperan objetos
Spatial*como entrada (cada vez menos, pero existen). - Lectura de tutoriales y libros pre-2020 que enseñan el flujo antiguo.
Cuándo NO usarlo
Para cualquier proyecto nuevo. Usa sf.
Migración
La conversión entre sp y sf es directa:
library(sf)
# sf -> sp
x_sp <- as(x_sf, "Spatial")
# sp -> sf
x_sf <- st_as_sf(x_sp)Casi todas las operaciones de sp + rgeos + rgdal tienen equivalente directo en sf:
sp / rgdal / rgeos |
sf |
|---|---|
readOGR() |
st_read() |
writeOGR() |
st_write() |
spTransform() |
st_transform() |
gIntersection() |
st_intersection() |
gBuffer() |
st_buffer() |
over() |
st_join() / st_intersects() |
proj4string() |
st_crs() |
Trampas habituales
- CRS en formato proj4 obsoleto. El antiguo
proj4stringcon strings tipo"+proj=longlat +datum=WGS84"está siendo retirado por GDAL/PROJ 6+. Migra a EPSG ("EPSG:4326") o WKT2. rgdal::readOGRen código viejo. Si te encuentras un script que falla porthere is no package called 'rgdal', el fix es reemplazarreadOGR(dsn, layer)porst_read(dsn, layer)y, si el resto del código sigue esperando un objetoSpatial*, envolverlo enas(..., "Spatial").
Enlaces
Relacionados en esta página
sf, sucesor moderno y destino obligado de cualquier migración.
raster (legacy)
raster fue el paquete dominante para datos raster en R durante más de una década, también de Robert J. Hijmans. Su sucesor oficial es terra, escrito por el mismo autor con una arquitectura moderna (C++ en lugar de R puro para los cuellos de botella) y una API más limpia.
raster sigue disponible en CRAN y mantenido en modo mínimo, pero el desarrollo activo está en terra. Cualquier proyecto nuevo debe empezar en terra.
Cuándo usarlo
- Mantenimiento de código heredado sobre
raster::RasterLayer/RasterStack/RasterBrick. - Paquetes científicos antiguos que todavía requieren objetos
Raster*como entrada, algunos modelos de distribución de especies, por ejemplo, tardaron años en migrar.
Cuándo NO usarlo
Para código nuevo. Usa terra.
Migración
La conversión es directa en ambas direcciones:
library(terra)
library(raster)
# raster -> terra
r_terra <- rast(r_raster)
# terra -> raster
r_raster <- raster(r_terra) # un layer
b_raster <- brick(r_terra) # multi-layerEquivalencias frecuentes:
raster |
terra |
|---|---|
raster() |
rast() |
brick() / stack() |
rast() (multi-layer) |
projectRaster() |
project() |
extract() |
extract() (compatible) |
calc() |
app() |
overlay() |
lapp() |
crs(r) <- CRS("...") |
crs(r) <- "EPSG:..." |
Trampas habituales
calc()yoverlay()no existen enterra. Si tu código viejo los usa, el equivalente esapp()(funciones por celda sobre una sola capa) olapp()(sobre varias capas alineadas).- CRS proj4 obsoleto. Igual que en
sp: migra a EPSG / WKT2. - Mezclar
rasteryterraen el mismo script funciona porque los objetos son convertibles, pero es propenso a errores. Decide una API y conviértelo todo al inicio del script.