I have a list of xy coordinates (originally from a raster) that I want to associate with their political unit (e.g., point a is in New York, point b is in Vermont). Associating a point with a separate raster or polygon is possible in terra using extract, as shown below. However, I can't figure out how to retain the xy coordinates when extracting using vector data. When using raster data, setting xy=TRUE works (in both the documentation and this answer to a similar question). When using vector data, this gives an error: unused argument (xy = TRUE)
Is there any way to do this using vector data? Or else is there a more optimized method to use when linking to vector data?
library(terra)
r <- rast(system.file("ex/elev.tif", package="terra")) #raster data
v <- vect(system.file("ex/lux.shp", package="terra")) #vector data
points <- spatSample(r, 5, as.points=TRUE) #generate points in the region
crds(points)
x y
[1,] 5.929167 49.67083
[2,] 6.279167 49.94583
[3,] 6.487500 49.84583
[4,] 6.004167 49.49583
[5,] 6.379167 49.62083
# associate the points with the vector data - but no xy coordinates
terra::extract(v, points)
id.y ID_1 NAME_1 ID_2 NAME_2 AREA POP
1 1 3 Luxembourg 8 Capellen 185 48187
2 2 NA <NA> NA <NA> NA NA
3 3 NA <NA> NA <NA> NA NA
4 4 3 Luxembourg 9 Esch-sur-Alzette 251 176820
5 5 2 Grevenmacher 12 Grevenmacher 210 29828
# with vector data, xy=TRUE gives an error
terra::extract(v, points, xy=TRUE)
Error in .local(x, y, ...) : unused argument (xy = TRUE)
# with raster data, xy=TRUE works
terra::extract(r, points, xy=TRUE)
ID elevation x y
1 1 323 5.929167 49.67083
2 2 NA 6.279167 49.94583
3 3 NA 6.487500 49.84583
4 4 358 6.004167 49.49583
5 5 283 6.379167 49.62083
You can combine the coordinates with the extracted values like this:
vp <- terra::extract(v, points)
data.frame(crds(points), vp)
# x y id.y ID_1 NAME_1 ID_2 NAME_2 AREA POP
#1 6.054167 49.69583 1 3 Luxembourg 11 Mersch 233 32112
#2 6.445833 49.85417 2 NA <NA> NA <NA> NA NA
#3 5.795833 49.80417 3 1 Diekirch 3 Redange 259 18664
#4 6.070833 50.02083 4 1 Diekirch 1 Clervaux 312 18081
#5 5.979167 49.82917 5 1 Diekirch 2 Diekirch 218 32543
You note that
# with raster data, xy=TRUE works
terra::extract(r, points, xy=TRUE)
That is correct, but the coordinates returned are for the centers of the cells of r, not the coordinates of points
If your points came from a raster, perhaps rasterize is what you ought to be using instead.