Extract values from layers in terra package-CodePudding

I have a terra::rast with n layers. I would like to extract values for certain layers only with the package terra (I prefer answers with this package rather than raster). It seems I do not get how to use the argument layer in the function terra::extract. Here is a minimum (almost) working example:

rs1 = terra::rast(nrows = 3, ncols = 4)
terra::values(rs1) = 4:15
rs2 = terra::rast(nrows = 3, ncols = 4)
terra::values(rs2) = 25:36
rs3 = terra::rast(nrows = 3, ncols = 4)
terra::values(rs3) = 97:108
rs4 = terra::rast(nrows = 3, ncols = 4)
terra::values(rs4) = 51:62

rs = c(rs1, rs2, rs3, rs4)
names(rs) = paste0("layer", 1:nlyr(rs))

coords = matrix(data = c(0, 5, 10, -54, 0, 12), byrow = FALSE, ncol = 2)
colnames(coords) = c("x", "y")

terra::extract(x = rs, y = coords)[, c("layer1", "layer4")] # Works but not elegant in terms of memory
terra::extract(x = rs, y = coords, layer = c("layer1", "layer4")) # Does not work
terra::extract(x = rs, y = coords, layer = c(1, 4)) # Does not work

This code provides the following error:

Error in .local(x, y, ...) : length(layer) == nrow(y) is not TRUE

I do not understand why the number of layers should equal the number of rows in y. I would imagine it should be at max the number of layers in x. How to extract the values at the given coordinates for the layers 1 and 4 exclusively (or whatever other layers, this is only an example)?

CodePudding user response：

I guess you are just looking for the subset() function of terra. Please find the reprex below.

Reprex

library(terra)

extract(x = subset(rs, c(1,4)), y = coords)
#>   layer1 layer4
#> 1     14     61
#> 2     10     57
#> 3     10     57

^{Created on 2021-11-02 by the reprex package (v0.3.0)}