如果单元格包含点，则将栅格单元格值更改为 NA

Question

我在单个 stack 中包含多个栅格图层，还有一个包含点坐标的 SpatialPoints 对象。如果单元格包含一个点，我正在尝试将栅格图层的值更改为 NA。

我在下面提供了一个可重现的例子。然而，为了说明这一点，我的真实数据由 10 层 'habitat'（高程、树冠覆盖等）组成，它们都包含在一个堆栈中。 rasters 覆盖的区域大约为 34 x 26 平方公里。此外，我还有大约 10,000 个来自动物的 GPS 位置。

因此，对于可重现的示例。

library(raster)
library(sp)

制作三个rasters并组合成一个stack

set.seed(123)
r1 <- raster(nrows=10, ncols=10)
r1 <- setValues(r1, sample(c(1:50), 100, replace = T))

r2 <- raster(nrows=10, ncols=10)
r2 <- setValues(r2, sample(c(40:50), 100, replace = T))

r3 <- raster(nrows=10, ncols=10)
r3 <- setValues(r3, sample(c(50:55), 100, replace = T))

Stack <- stack(r1, r2, r3)
nlayers(Stack)

制作 SpatialPoints

Pts <- SpatialPoints(data.frame(x = sample(-175:175, 50, replace = T),
                                y = sample(-75:75, 50, replace = T)))

绘制其中一个栅格和点

plot(Stack[[2]])
plot(Pts, add = T)

现在，是否可以将至少包含一个点的每个单元格的值更改为 NA？在 stack.

上执行此操作会很棒

Answer 1

我会使用 extract()，它可以被要求 return 每个点所在的栅格像元的像元编号：

ii <- extract(Stack, Pts, cellnumbers=TRUE)[,"cells"]
Stack[ii] <- NA

## Check any one of the layers to see that this worked:
plot(Stack[[2]])
plot(Pts, add=TRUE)

或者，可以使用 rasterize()，但对于非常大的栅格，它可能（？）会更慢：

ii <- !is.na(rasterize(Pts, Stack))
Stack[ii] <- NA

如果单元格包含点，则将栅格单元格值更改为 NA

Change raster cell values to NA if cell contains a point

r

raster

sp