如何使用 terra 包将从栅格堆栈中提取的值添加到 Spatial 对象的 data.frame？
How to add the extracted values from raster stack to the data.frame of the Spatial object using terra package?
我想使用 terra 包
将从栅格堆栈中提取的值添加到空间对象的 data.frame
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
正如您从输出中看到的那样，val 文件不存在我们过去使用 raster 包的 sp 参数获得的文件。现在如何在提取的 data.frame 中添加 val 字段？
由于 terra::extract 按照向量中元素的顺序提取值，您可以简单地 cbind 将字段添加到提取的数据框中，或者使用 $ 添加所需的字段作为列
ex <- terra::extract(r, points_vect, xy = T, method = "simple")
ex$val <- points_vect$val
编辑
提取栅格值最快的选项是 exactextract 包中的 exact_extract，它有一个参数（称为 append_cols）将矢量场附加到 data.frame 提取的结果。
exact_extract 仅接受多边形数据（作为 sf 对象）作为输入，但您可以通过对点应用缓冲区并汇总值来解决此问题。请参阅下面的示例
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")
只要您是用点提取（或使用线和多边形的汇总函数），就可以使用 cbind 和 terra
示例数据
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
使用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
并将提取的值 v 与 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
或者像这样替换原来的值
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
或者像这样创建一个新的SpatVector
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
我想使用 terra 包
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
正如您从输出中看到的那样，val 文件不存在我们过去使用 raster 包的 sp 参数获得的文件。现在如何在提取的 data.frame 中添加 val 字段？
由于 terra::extract 按照向量中元素的顺序提取值，您可以简单地 cbind 将字段添加到提取的数据框中，或者使用 $ 添加所需的字段作为列
ex <- terra::extract(r, points_vect, xy = T, method = "simple")
ex$val <- points_vect$val
编辑
提取栅格值最快的选项是 exactextract 包中的 exact_extract，它有一个参数（称为 append_cols）将矢量场附加到 data.frame 提取的结果。
exact_extract 仅接受多边形数据（作为 sf 对象）作为输入，但您可以通过对点应用缓冲区并汇总值来解决此问题。请参阅下面的示例
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")
只要您是用点提取（或使用线和多边形的汇总函数），就可以使用 cbind 和 terra
示例数据
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
使用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
并将提取的值 v 与 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
或者像这样替换原来的值
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
或者像这样创建一个新的SpatVector
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