Reputation: 8136
I want to add the extracted values from raster stack to the data.frame of the Spatial object using terra
package
f <- system.file("ex/logo.tif", package="terra")
r <- rast(f)
#Plot the raster
plot(r, 1:3)
#Create a vector file
points <- cbind.data.frame(Longitude = 40, Latitude = 40, val = 0.5)
points <- rbind.data.frame(points-1, points, points+1)
points_vect <- vect(points, geom=c("Longitude", "Latitude"), crs=crs(r, proj=T),
type = "points")
#Extract the values from raster using the point vector
terra::extract(r, points_vect, xy = T, method = "simple")
#> ID red green blue x y
#> 1 1 255 255 253 30 30
#> 2 2 149 159 186 40 40
#> 3 3 146 156 207 50 50
As you can see from the output val
filed is not there which we used to get using sp
aurgument of raster
package. Now how can I have the val
field added in the extracted data.frame?
Upvotes: 0
Views: 2016
Reputation: 47036
As long as you are extracting with points (or use a summarizing function with lines and polygons), you can use bind=TRUE
Example data
f <- system.file("ex/logo.tif", package="terra")
r <- rast(f)
points <- data.frame(Longitude = 40, Latitude = 40, val = 0.5)
points <- rbind(points-1, points, points+1)
points_vect <- vect(points, geom=c("Longitude", "Latitude"), crs=crs(r))
points_vect
# class : SpatVector
# geometry : points
# dimensions : 3, 3 (geometries, attributes)
# extent : 39, 41, 39, 41 (xmin, xmax, ymin, ymax)
# coord. ref. : +proj=merc +lon_0=0 +k=1 +x_0=0 +y_0=0 +datum=WGS84 +units=m +no_defs
# names : Longitude Latitude val
# type : <num> <num> <num>
# values : 39 39 -0.5
# 40 40 0.5
# 41 41 1.5
Solution
v <- terra::extract(r, points_vect, xy = T, method = "simple", bind=T)
v
# class : SpatVector
# geometry : points
# dimensions : 3, 6 (geometries, attributes)
# extent : 39, 41, 39, 41 (xmin, xmax, ymin, ymax)
# coord. ref. : Cartesian (Meter)
# names : val red green blue x y
# type : <num> <int> <int> <int> <num> <num>
# values : -0.5 247 254 255 39.5 38.5
# 0.5 149 159 186 40.5 39.5
# 1.5 132 141 182 41.5 40.5
Alternative approaches:
First extract
v <- terra::extract(r, points_vect, xy = T, method = "simple")
v
# ID red green blue x y
#1 1 247 254 255 39 39
#2 2 149 159 186 40 40
#3 3 132 141 182 41 41
Then use cbind
to combine the extracted values v
with points_vect
x <- cbind(points_vect, v)
x
# class : SpatVector
# geometry : points
# dimensions : 3, 9 (geometries, attributes)
# extent : 39, 41, 39, 41 (xmin, xmax, ymin, ymax)
# coord. ref. : +proj=merc +lon_0=0 +k=1 +x_0=0 +y_0=0 +datum=WGS84 +units=m +no_defs
# names : Longitude Latitude val ID red green blue x y
# type : <num> <num> <num> <num> <num> <num> <num> <num> <num>
# values : 39 39 -0.5 1 247 254 255 39 39
# 40 40 0.5 2 149 159 186 40 40
# 41 41 1.5 3 132 141 182 41 41
Or replace the original values like this
values(points_vect) <- v
points_vect
# class : SpatVector
# geometry : points
# dimensions : 3, 6 (geometries, attributes)
# extent : 39, 41, 39, 41 (xmin, xmax, ymin, ymax)
# coord. ref. : +proj=merc +lon_0=0 +k=1 +x_0=0 +y_0=0 +datum=WGS84 +units=m +no_defs
# names : ID red green blue x y
# type : <num> <num> <num> <num> <num> <num>
# values : 1 247 254 255 39 39
# 2 149 159 186 40 40
# 3 132 141 182 41 41
You can also create a new SpatVector
like this
newpts <- vect(v, c("x", "y"), crs=crs(r))
newpts
# class : SpatVector
# geometry : points
# dimensions : 3, 6 (geometries, attributes)
# extent : 39, 41, 39, 41 (xmin, xmax, ymin, ymax)
# coord. ref. : +proj=merc +lon_0=0 +k=1 +x_0=0 +y_0=0 +datum=WGS84 +units=m +no_defs
# names : ID red green blue x y
# type : <num> <num> <num> <num> <num> <num>
# values : 1 247 254 255 39 39
# 2 149 159 186 40 40
# 3 132 141 182 41 41
Upvotes: 6
Reputation: 2584
Since terra::extract
extract values with the order of elements in the vector, you can simply cbind
the field to the extracted dataframe or add the needed field as a column with $
ex <- terra::extract(r, points_vect, xy = T, method = "simple")
ex$val <- points_vect$val
EDIT
The fastest option for extracting raster values is exact_extract
in the exactextract
package, which has an argument (called append_cols
) to attach the vector field to the data.frame resulting from the extraction.
exact_extract
accepts only polygon data (as sf
object) as input, but you can get around this by applying a buffer to your points and summarise the values. See the example below
library(exactextractr)
library(sf)
point_sf <- st_as_sf(points_vect)
point_sf <- st_buffer(point_sf,1)
ex <- exactextractr::exact_extract(stack(r),point_sf,fun="mean",append_cols="val")
#you can summarise values using a vector of functions
ex <- exactextractr::exact_extract(stack(r),point_sf,fun=c("mean","median","max","min"),append_cols="val")
Upvotes: 3