MappingR - Visualización de datos espaciales en R
Beatriz M. Meriño - bvergaram@gmail.com
04 de Junio de 2024
CONFIGURACIONES IMPORTANTES
Clonar Git-Hub
Crear nuevo proyecto
Version control
Git
BASES DE DATOS LIBRES
Límites Administrativos. GISCO
GISCO es un repositorio abierto de datos geoespaciales que incluye varios conjuntos de datos como países, líneas costeras, etiquetas o niveles NUTS. Los conjuntos de datos se proporcionan normalmente en varios niveles de resolución (60M/20M/10M/03M/01M) y en 3 proyecciones diferentes (4326/3035/3857).
Tenga en cuenta que la librería no proporciona metadatos sobre los archivos descargados, la información está disponible en la página web API webpage.
Página completa con ejemplos y viñetas en https://ropengov.github.io/giscoR/
Buscar información para extraer límites administrativos:
Extraer límites administrativos de Chile.
Límites Administrativos usando geodata.
- geodata es una librería para descargar datos
geográficos. Facilita el acceso a datos de clima, elevación, suelo,
cultivos, presencia de especies y límites administrativos, y es un
sucesor de la función
getData()del paquete raster.
Buscar información sobre la función gadm
Descargar limites administrativos para Chile
#chl <- gadm(country="Chile", level=1, path = ("data/layers"))
chl <- readRDS("data/layers/gadm/gadm41_CHL_1_pk.rds")
plot(chl)
Datos ambientales usando geodata
Worldclim. Temperatura Mínima usando geodata
WorldClim es una base de datos de alta resolución espacial de datos meteorológicos y climáticos mundiales. Puede descargar datos meteorológicos y climáticos de cuadrícula para históricos (casi actual) y futuros.
Hay datos climáticos mensuales de temperatura mínima, media y máxima precipitación, radiación solar, velocidad del viento, presión de vapor de agua y precipitaciones totales.
chile.tmin <- worldclim_country(country = "Chile",
var = "tmin",
path = ("data/climatic"))
terra::plot(mean(chile.tmin), plg = list(title ="Temperatura Minima Anual (C)"))
Worlclim. Bioclim usando geodata
- También hay 19 variables «bioclimáticas». Las variables bioclimáticas derivan de los valores mensuales de temperatura y precipitaciones para generar variables con mayor significado biológico. A menudo se utilizan en modelos de distribución de especies y en técnicas de modelado ecológico relacionadas.
chile.bios <- worldclim_country(country = "Chile",
var = "bio",
path = ("data/climatic"))
terra::plot(chile.bios)
Soilgrid usando geodata
- SoilGrids es un sistema de cartografía digital que utiliza información global sobre el perfil del suelo y datos de covarianza para modelar la distribución espacial de las propiedades del suelo en todo el mundo.
#soil <- soil_world(var = c("nitrogen", "ocd", "phh2o",
# "sand", "silt"),
# depth = 60, path = "data/climatic")
#chile.soil <- crop(soil,chile)
#writeRaster(chile.soil, filename = "data/climatic/chile.soil.tif")
chile.soil <- rast("data/climatic/chile.soil.tif")
terra::plot(chile.soil)
Elevación usando geodata
chile.elevacion <- elevation_30s(country="Chile",
path=("data/elevation"),
mask = TRUE)
terra::plot(chile.elevacion)
Elevación usando elevatr
- El paquete elevatr proporciona acceso a datos de elevación de Datos abiertos de AWS Terrain Tiles y el Open Topography Global datasets API para modelos de elevación raster. Para los datos de elevación de puntos, pueden utilizarse USGS Elevation Point Elevation Point Query Service o las elevaciones puntuales pueden derivarse de los AWS Tiles.
#chile_elevacion <- get_elev_raster(locations = chile,
# z = 9,
# clip = "locations",
# neg_to_na = TRUE,
# override_size_check = TRUE)
#writeRaster(chile_elevacion, "data/elevation/CHL_elv_.tif")
chile_elevacion <- rast("data/elevation/CHL_elv_6.tif")
terra::plot(chile_elevacion)
Mosaico Satelital usando ggmap
- ggmap es un paquete de R que facilita la recuperación de mapas raster de servicios cartográficos en línea como GoogleMaps, Stadia Maps, y OpenStreetMap y representarlos gráficamente con el programa ggplot2.
Stadia Maps
Stadia Maps ofrece recursos cartográficos de varios estilos, incluidos los actualizados tiles from Stamen Design. Se requiere una clave API.
revisar: registrarse y hay un nivel gratuito para uso no comercial. Una vez que tengas tu clave API, invoca el registro función:
register_stadiamaps(«YOUR-API-KEY», write = F).
#Reemplaza con tu clave API
#??register_stadiamaps
#https://client.stadiamaps.com/signup/
#Eliminar el #
# register_stadiamaps(key = "TU-API-AQUI")bbox <- c(left = -80, bottom = -60, right = -65, top = -15)
map <- get_stadiamap(bbox, zoom = 6, maptype = "stamen_terrain_background",
where = "data/layers") ## ℹ © Stadia Maps © Stamen Design © OpenMapTiles © OpenStreetMap contributors.## ℹ 48 tiles needed, this may take a while (try a smaller zoom?)
PLOT 1: Anotación del mapa
Este tutorial es una versión modificada del proyecto y la sección de datos espaciales del Coding Club. Puedes modificar el código con sus propios datos para replicar las funciones.
Cargar Datos de Ocurrencias
Cargar Mapas de referencia
## Warning: plotting the first 9 out of 121 attributes; use max.plot = 121 to plot
## all
## Warning: plotting the first 10 out of 168 attributes; use max.plot = 168 to
## plot all
Plotear Mapa con Puntos de Ocurrencias
raw.map <- plot.data %>%
ggplot() +
geom_sf(data = regions, colour = "white", fill = "gray40") +
geom_point(aes(x = Longitude,
y = Latitude,
color = Species.name),
alpha = 0.5) +
labs(color = "Species occurrence") +
xlim(-77, -66) +
theme_bw() +
scale_size_area()
print(raw.map)
Ajustar y añadir referencias geográficas
adj.map <- plot.data %>%
ggplot() +
geom_sf(data = regions,
color = "white",
fill = "gray40") +
geom_point(aes(x = Longitude,
y = Latitude,
color = Species.name),
alpha = 0.5) +
labs(color = "Species") +
xlim(-74, -70) +
ylim(-35.5, -28.5) +
theme_test() +
guides(size = "none") +
annotation_scale(location = "bl") +
annotation_north_arrow(pad_y = unit(0.7, "cm"),
pad_x = unit(0.6, "cm"),
height = unit(1, "cm"),
width = unit(1, "cm"),
which_north = "T")
print(adj.map)
Añadir Referencia Continental
inset <- ggplot(data = SurAmerica) +
geom_sf(color = "gray40",
fill = "gray40") +
annotate("rect",
xmin = -72,
xmax = -69,
ymin = -25.5,
ymax = -28.5,
color = "red",
fill = NA) +
xlim(-80, -36) +
labs(x = NULL,
y = NULL) +
theme_test() +
theme(axis.text = element_blank(),
axis.ticks = element_blank(),
axis.ticks.length = unit(0, "pt"),
axis.title = element_blank(),
plot.margin = margin(0, 0, 0, 0, "cm"))
print(inset)
Insertar Mapa de Referencia en el Mapa Ajustado
final.map <- ggdraw() +
draw_plot(adj.map) +
draw_plot(inset, x = 0.234,
y = 0.665,
width = 0.3,
height = 0.3)
print(final.map)
Exportar Plots
Mantener las dimensiones de un plot es crucial. Siempre lea las recomentaciones para autores de la revista y evite que el tamaño del plot supere una hoja A4 impresa. Si necesita corregir detalles usando un editor vectorial como Inkscape o Illustrator, el formato vectorial .pdf o .svg es una buena práctica. La revista siempre señala el formato con el cual debe entregar sus salidas gráficas.
Guardar en formato PDF
PLOT 2: Plot DEM
Los modelos de elevación digital (DEM) pueden ser utilizados para crear mapas complejos. Basados en un modelo de elevación, se estiman valores como pendientes y sombras para crear un efecto de relieve sin tener que renderizar un gráfico en 3D. Considera siempre factores como el área que quieres gráficar, la resolución y el tamaño de tu archivo ráster y la paleta de colores. Este tutorial es una dataptación del tutorial de tidyterra y como crear efectos de sombreado.
Cargar Mapas de referencia
- Usaremos Chile para descargar datos de elevación y “Sur América” como referencia continental.
regiones <- ne_states(country = "Chile",
returnclass = "sf")
SurAmerica <- ne_countries(continent = "South America",
returnclass = "sf")Elevación usando geodata
chile.elevacion <- elevation_30s(country="Chile",
path=("data/elevation"),
mask = TRUE)
names(chile.elevacion) <- "elv"
terra::plot(chile.elevacion)
#ajustar en caso de que el dem tenga valores negativos
chile.elevacion <- chile.elevacion %>%
mutate(elv = pmax(0, elv))
terra::plot(chile.elevacion)
Cortar por área de interés
- Ajustaremos el área de interes a Chile central para poder visualizar el efecto del sombreado en la pendiente. Recordar que en áreas más grandes el tiempo de proceso puede demorar más y se necesitaran mayores recursos computacionales.
#limites del raster
#
elv_limits <- minmax(elv) %>% as.vector()
# redondear los limites de los quiebres a 500
elv_limits <- c(floor(elv_limits[1] / 500),
ceiling(elv_limits[2] / 500)) * 500
# min value = 0.
elv_limits <- pmax(elv_limits, 0)
# Comparar
minmax(elv) %>% as.vector()## [1] 0 5497## [1] 0 5500#plot dem
grad <- hypso.colors(10)
dem.plot<- autoplot(elv) +
scale_fill_gradientn(colours = grad, na.value = NA)## Scale for fill is already present.
## Adding another scale for fill, which will replace the existing scale.
Crear efecto de sombreado
- Calcular las características del terreno a partir de los
datos de elevación. Los valores de elevación deben estar en las mismas
unidades que las unidades del mapa (normalmente metros) para datos
ráster proyectados (planos). Básicamente lo que queremos crear es una
capa que se aproxime a la «textura» potencial de la superficie en
función de la elevación y la posición del sol. Esto es sencillo con las
funciones
terra::terrain()yterra::shade().
slope <- terrain(elv, "slope", unit = "radians")
aspect <- terrain(elv, "aspect", unit = "radians")
hill <- shade(slope, aspect, 0, 90)
names(hill) <- "shades"
# Sombreado
pal_greys <- hcl.colors(1000, "Grays")
hillshade <- ggplot() +
geom_spatraster(data = hill) +
scale_fill_gradientn(colors = pal_greys, na.value = NA)
print(hillshade)
Ajustar la paleta de colores
pal_greys <- hcl.colors(1000, "Grays")
index <- hill %>%
mutate(index_col = rescale(shades,
to = c(1, length(pal_greys)))) %>%
mutate(index_col = round(index_col)) %>%
pull(index_col)
vector_cols <- pal_greys[index]
hill_plot <- ggplot() +
geom_spatraster(data = hill,
fill = vector_cols,
maxcell = Inf,
alpha = 1)
print(hill_plot)
Base plot con tidyterra
- tidyterra proporciona gradientes adicionales cuyos colores se colocan de forma desigual con el objetivo de proporcionar una mejor comprensión de los mapas. Por ejemplo, por convención los colores verdes se asocian a elevaciones bajas mientras que los naranjas, marrones y blancos se asocian a elevaciones altas en algunas de las paletas de elevación más comunes (también conocidas como tintes hipsométricos). Paletas de colores para representar la hipsometría y la batimetría con tidyterra. Referencias acá.
base_plot <- hill_plot +
geom_spatraster(data = elv, maxcell = Inf) +
scale_fill_hypso_tint_c(
limits = elv_limits,
palette = "dem_screen",
alpha = 0.5,
labels = label_comma())
print(base_plot)
Ajustar base plot
base_text_size <- 10
base.map <- base_plot +
guides(fill = guide_legend(title = " Altitude m.",
direction = "vertical",
keywidth = 1,
keyheight = 1,
label.position = "right",
title.position = "top",
override.aes = list(alpha = 0.5))) +
labs(x="Longitude", y="Latitude") +
theme_minimal() +
theme(plot.background = element_rect("grey97", colour = NA),
plot.margin = margin(20, 20, 20, 20),
plot.caption = element_text(size = base_text_size * 0.5),
plot.title = element_text(face = "bold",
size = base_text_size * 0.9),
axis.text = element_text(size = base_text_size),
legend.position = "right",
legend.title = element_text(size = base_text_size * 0.8),
legend.text = element_text(size = base_text_size * 0.8),
legend.key = element_rect("grey50"),
legend.spacing.x = unit(2, "pt"))
print(base.map)
Plot datos de ocurrencia
Crear una paleta de colores personalizada
## J. chilensis O. degus
## "#FF9933" "#0077B5"occ.data <- occ.data %>%
mutate(Species.name = gsub("_", ". ", Species))
##Subset de palma
palma.occs <- subset(occ.data, Species == "J_chilensis")
head(palma.occs)species.plot <- base.map +
geom_point(data = palma.occs,
aes(x = Longitude, y = Latitude),
alpha = 0.5,
colour= "gray20",
size=2.5) +
geom_point(data = palma.occs,
aes(x = Longitude,
y = Latitude,
colour = Species.name),
alpha = 0.75, size=2) +
scale_color_manual(values = palmas.pal) +
labs(color = "Species occurence") +
annotation_scale(location = "bl") +
annotation_north_arrow(pad_y = unit(0.7, "cm"),
pad_x = unit(0.6, "cm"),
height = unit(1, "cm"),
width = unit(1, "cm"),
which_north = "T")
print(species.plot)
Añadir Referencia Continental
inset <- ggplot(data = SurAmerica) +
geom_sf(color = "gray40",
fill = "gray40") +
annotate("rect",
xmin = -72,
xmax = -69,
ymin = -25.5,
ymax = -28.5,
color = "red",
fill = NA) +
xlim(-80, -36) +
labs(x = NULL,
y = NULL) +
theme_void() +
theme(panel.background = element_rect(fill = "transparent",
color = NA),
plot.background = element_rect(fill = "transparent",
color = NA),
axis.text = element_blank(),
axis.ticks = element_blank(),
axis.ticks.length = unit(0, "pt"),
axis.title = element_blank(),
plot.margin = margin(0, 0, 0, 0, "cm"))
print(inset)
Insertar Mapa de Referencia en el Mapa Ajustado
final.map <- ggdraw() +
draw_plot(species.plot) +
draw_plot(inset, x = 0.225, y = 0.62, width = 0.3, height = 0.3)
print(final.map)
Exportar Plots
Guardar en formato PDF
Guardar en formato JPG
GUARDAR IMAGEN
Mantener un registro de los archivos de análisis es crucial. No queremos volver a esperar un par de días para corregir el plot porque algunos análisis pueden tardar. Crear un documento en formato R markdown es una buena idea para trabajar bajo los estándares de ciencia reproducible porque te permite compartir todo tu código hasta el plot final.