The code works with one raster and one shapefile. I tried to expand the code to use stack() for the two rasterfiles (present and future climate, paired boxplots), because extract(r,v) is available for one raster and one shapefile. Do you know of any other possibilities? Isuppose it has to be within library(terra)
# load libraries
library(terra)
library(ggplot2)
library(raster)
# Load shp and raster from example files of terra
(v <- vect(system.file("ex/lux.shp", package="terra")))
(r1 <- rast(system.file("ex/elev.tif", package="terra")))
(r2 <- rast(system.file("ex/elev.tif", package="terra")))
r1r2<-stack(r1,r2)
Error in .local(x, ...) :
unused argument (new("SpatRaster", ptr = new("Rcpp_SpatRaster", .xData = <environment>)))
# extract, include names from vector attributes
e <- extract(r,v)
e$NAME_2 <- v$NAME_2[e$ID]
str(e)
# plot
ggplot(e, aes(x = NAME_2, y = elevation)) +
geom_boxplot() +
theme_minimal() +
theme(axis.text.x = element_text(angle = 45))
Your question really seems to be "how do I combine two SpatRasters into a single object. In your context, the answer is use c
library(terra)
v <- vect(system.file("ex/lux.shp", package="terra"))
r1 <- rast(system.file("ex/elev.tif", package="terra"))
r2 <- rast(system.file("ex/elev.tif", package="terra")) * 2
x <- c(r1, r2)
names(x) <- c("A", "B")
#and extract values
e <- extract(x, v)
e$name <- v$NAME_2[e$ID]
To make a boxplot for each polygon and layer, see ?boxplot
(or use another plotting approach such as ggplot). For example
f <- reshape(e, varying = names(x), v.names="value", times=names(x), direction="long")
par(mar=c(3,6,1,1), cex.axis=.5)
boxplot(value ~ time:name, data = f, horizontal=TRUE, cex=.25,
boxwex = 0.5, col = c("orange", "yellow"), las=1,
sep = ":", lex.order = F, yaxs = "i", ylab="")