基于Terra包，按非唯一NAME_1提取多边形栅格中值并添加属性

给定以下R语言（terra包）的栅格与多边形数据集：

f <- system.file("ex/elev.tif", package="terra") 
r <- rast(f)
fs <- system.file("ex/lux.shp", package="terra")
v <- vect(fs)

矢量数据v的NAME_1属性并非唯一，其属性信息如下：

v
# class       : SpatVector 
# geometry    : polygons 
# dimensions  : 12, 6  (geometries, attributes)
# extent      : 5.74414, 6.528252, 49.44781, 50.18162  (xmin, xmax, ymin, ymax)
# source      : lux.shp
# coord. ref. : lon/lat WGS 84 (EPSG:4326) 
# names       :  ID_1   NAME_1  ID_2   NAME_2  AREA   POP
# type        : <num>    <chr> <num>    <chr> <num> <int>
# values      :     1 Diekirch     1 Clervaux   312 18081
#                   1 Diekirch     2 Diekirch   218 32543
#                   1 Diekirch     3  Redange   259 18664

实现步骤与代码 #

要提取相同NAME_1值对应的栅格中值并添加到v的属性表，可按以下步骤操作：

1. 按NAME_1合并多边形
   将所有具有相同NAME_1的多边形合并为单个面要素，仅保留NAME_1属性用于后续匹配：

v_agg <- aggregate(v, by = "NAME_1", fun = NULL)

2. 提取合并多边形的栅格中值
   使用extract函数提取每个合并多边形对应的栅格单元格中值：

# 提取中值，忽略NA值
zonal_stats <- extract(r, v_agg, fun = median, na.rm = TRUE)
# 整理结果，保留NAME_1和中值列
zonal_stats <- cbind(v_agg$NAME_1, zonal_stats[, 2, drop = FALSE])
colnames(zonal_stats) <- c("NAME_1", "elev_median")

3. 关联中值到原矢量属性表
   通过merge函数将统计得到的中值按NAME_1匹配回原矢量数据的属性表：

v <- merge(v, zonal_stats, by = "NAME_1")

4. 验证结果
   查看更新后的矢量属性表，确认elev_median列已添加：

head(v)

内容的提问来源于stack exchange，提问作者Tpellirn