区域统计 R

Question

我在 R 中使用 raster 包并在 shapefile 包含多个 polygons (200) 的光栅上执行 extract。每个多边形都有一列，其中指定了多边形的 class。当我执行 extract 时，我得到一个 dataframe，其中栅格的每个值都是 "assigned" 到多边形：哪个多边形是哪个像素所在的位置。但是，我的多边形包含的信息（即类是什么）消失了。 extract 之后的数据框类似于：

test <- extract(myRaster, myTrainingPolygon, df = TRUE)

ID     band1     band2    band3
1      0.101     0.827    ...
...    ...       ...      ...
200    0.876     0.821    ...

我需要的是

ID     band1     band2    band3    class
1      0.101     0.827    ...      class1
...    ...       ...      ...      ...
200    0.876     ...      ...      class3

我如何获取其中的信息 - 或者首先 - 在执行提取时不要丢失它们？！

Answer 1

?extract 表示 sp 参数控制是否应将提取的值添加到空间对象的数据框中。指定 sp=TRUE 应该可以解决问题。例如：

library(maptools)
library(raster)
data(wrld_simpl)
r <- raster(extent(-180, 180, -90, 90), res=10)
r[] <- runif(ncell(r))
wrld_simpl_new <- extract(r, wrld_simpl, fun=mean, sp=TRUE)

head(wrld_simpl_new)

##     FIPS ISO2 ISO3 UN                NAME   AREA  POP2005 REGION SUBREGION     LON     LAT     layer
## ATG   AC   AG  ATG 28 Antigua and Barbuda     44    83039     19        29 -61.783  17.078 0.9142067
## DZA   AG   DZ  DZA 12             Algeria 238174 32854159      2        15   2.632  28.163 0.1774097
## AZE   AJ   AZ  AZE 31          Azerbaijan   8260  8352021    142       145  47.395  40.430 0.3098710
## ALB   AL   AL  ALB  8             Albania   2740  3153731    150        39  20.068  41.143 0.3746480
## ARM   AM   AM  ARM 51             Armenia   2820  3017661    142       145  44.563  40.534 0.3494729
## AGO   AO   AO  AGO 24              Angola 124670 16095214      2        17  17.544 -12.296 0.2873931

添加的列的名称为 "layer"，因为这是 RasterLayer 的（默认）名称。

Answer 2

始终包含和使用示例数据

library(raster)
r <- raster(ncol=36, nrow=18)
r[] <- 1:ncell(r)
cds1 <- rbind(c(-180,-20), c(-160,5), c(-60, 0), c(-160,-60), c(-180,-20))
cds2 <- rbind(c(80,0), c(100,60), c(120,0), c(120,-55), c(80,0))
p <- spPolygons(cds1, cds2)
p <- SpatialPolygonsDataFrame(p, data.frame(class=c('A', 'B'), stringsAsFactors=FALSE))

现在我们可以做：

v <- extract(r, p)
str(v)

#List of 2
# $ : num [1:38] 326 327 328 329 330 331 332 333 334 335 ...
# $ : num [1:25] 172 173 208 209 244 245 279 280 281 282 ...

创建一个 data.frame，其中包含多边形顺序 ID 和提取的值

d <- data.frame(id=rep(1:length(v), sapply(v, length)), value=unlist(v))

（这种使用 unlist 的特殊方法仅适用于单层）。对于多层做

#d <- data.frame(id=rep(1:length(v), sapply(v, length)), do.call(rbind, v))

使用顺序 ID 和所需的其他多边形属性创建 data.frame

pd <- cbind(id=1:length(p), data.frame(p))

合并两者

m <- merge(pd, d)

区域统计 R

Zonal Statistics R

geography

r

polygon

raster

r-raster