[R][r] 光栅操作，从缓冲区中提取值而不重叠

Question

我和我的研究伙伴正在撰写一篇关于巴西蝙蝠的适宜性/种群数量的论文 我们正在努力操纵光栅来实现我们的方法， 模型中将使用大约 200 个栅格（不同物种的不同场景） 我们为此问题选择并裁剪了一个

我们正在尝试做的事情如下：

• 根据像素周围的缓冲区从每个栅格中提取信息，这些像素的值高于阈值（在本例中为 0.8）
• 理想情况下，重叠的缓冲区将仅选择要保留的像素之一。创建一组生成非重叠缓冲区的像素。
• 此外，对于下一步，我们希望根据当前值和线性模型更改栅格值（这部分很好）
• 值更改后，我们希望合并来自同一场景的栅格以寻找物种之间的相似性

我们一直在努力寻找一种不会产生重叠缓冲区或以更规范的方式聚合像素的方法。

我们的栅格具有 30 弧秒（大约 1 平方公里）的分辨率，我们只能聚合像素，但我们需要特定区域，并且不可能只聚合整个像素以生成 6 平方公里等值 另一种可能性是聚合部分像素以访问例如 6 平方公里的区域。有可能吗？

到目前为止，我已经尝试了两种方法（提取和聚合），但都没有返回所需的结果。以下是目前使用的代码：

library(raster)
# library(ggplot2)
library(dplyr)

raster_values <- c(0.7151692, 0.7234125, 0.7242436, 0.8838134, 0.9855102, 0.9921246, 0.9679778, 0.8245632, 0.8965716, NA, NA, NA, 0.721549, 0.6988058, 0.7333487, 0.8138089, 0.829727, 0.8689544, 0.8607966, 0.9794912, 0.9012381, NA, NA, NA, 0.7118917, 0.7103891, 0.7527786, 0.7792872, 0.8320968, 0.8082156, 0.920545, 0.9788723, 0.7859345, NA, NA, NA, 0.6824703, 0.7039984, 0.7589136, 0.7905939, 0.9024587, 0.9848859, 0.969553, 0.8404503, NA, NA, NA, 0.9922802, 0.6883243, 0.731391, 0.7682831, 0.8586601, 0.9850862, 0.9705435, 0.9888299, 0.8164479, 0.666971, NA, NA, 0.757661, 0.7628939, 0.7868114, 0.7910978, 0.7650495, 0.8227689, 0.8148086, 0.8691386, 0.7376462, 0.9176324, NA, NA, 0.998813, 0.7585487, 0.7721, 0.6508481, 0.6098195, 0.7708853, 0.7119401, 0.625409, 0.6886432, 0.8641906, NA, NA, 0.9991203, 0.7227083, 0.6550816, 0.5863016, 0.7050957, 0.6629267, 0.6550342, 0.6217013, 0.8762864, 0.9989462, NA, NA, NA, 0.6901366, 0.7041199, 0.6307223, 0.6305411, 0.7033732, 0.7092581, 0.7340803, 0.7865254, 0.9964261, 0.7940102, NA, NA, 0.6661926, 0.7381653, 0.6544684, 0.6170949, 0.641997, 0.7506128, 0.9248958, 0.9903375, 0.7662657, 0.7847621, NA, NA, 0.6582187, 0.7361452, 0.6488761, 0.6309077, 0.6542051, 0.781707, 0.940975, 0.9350743, 0.7824667, NA, NA, NA, 0.7215696, 0.7388574, 0.6753907, 0.6716958, 0.7162136, 0.7920918, 0.8702987, 0.9929227, 0.9775091, NA, NA, NA)

adeq_raster <- raster(nrows=12, ncols=12, xmn=-49, xmx=-48.5, ymn=-28, ymx=-27.5, crs = "+proj=longlat +datum=WGS84 +no_defs ", vals=raster_values)
adeq_raster_df <- as.data.frame(adeq_raster, xy = T)

# ggplot()+
#   geom_raster(data = adeq_raster_df, aes(x = x, y = y, fill = layer))+
#   scale_fill_viridis_c()+
#   coord_quickmap()


# - Extract method
#this method works fine, but does not avoid overlapping buffer.
adeq_raster_df_points <- filter(adeq_raster_df, layer >= 0.80) %>% 
  dplyr::select(x,y) #seleciona apenas as duas colunas de localização
#transfor the points df into a sp format to be used in the extract function and maintain the coordenates values.
adeq_raster_df_points_sp <- sp::SpatialPoints(coordinates(adeq_raster_df_points))

adeq_raster_extract <- extract(adeq_raster,adeq_raster_df_points_sp, buffer = 5000, fun = mean, na.rm=TRUE, sp = TRUE)
adeq_raster_extract_df <- as.data.frame(adeq_raster_extract)


# - Aggregate method
# - This method have not been so promissing, beucase it restrain the agregation to a certain number of pixels and we need buffers based in distance (3,5 or 10 km).
adeq_raster_5km <- aggregate(adeq_raster, 2, fun=mean, expand=F, na.rm=TRUE)#, filename='output/adeq_raster_5km')
adeq_raster_5km_df <- as.data.frame(adeq_raster_5km, xy = T)

我们欢迎任何有关如何处理缓冲区的建议。

Answer 1

请在代码中提供示例数据——应该没有理由下载数据（请参阅 R 帮助文件和这些页面以获取 100 多个示例）。提出一个重点问题也很好。

我认为，这里是为了获得你想要的缓冲区，使用terra 包。你可以和raster 和朋友做类似的事情。

示例数据：

library(terra)
r <- rast(nrow=20, ncol=20, xmin=0, xmax=1, ymin=0, ymax=1, crs="+proj=utm +zone=1 +datum=WGS84")
set.seed(123)
values(r) <- runif(ncell(r))
r <- ifel(r < .97, NA, r)

制作点和缓冲

p <- as.points(r)
b <- buffer(p, 0.1)

找到不相交的缓冲区

x <- relate(b, relation="intersects")
i <- rowSums(x, na.rm=TRUE) == 1
bb <- b[i,]

或者你想要合并重叠缓冲区？在这种情况下，您可以使用聚合（可能后跟disaggregate 来创建单个多边形）

aa <- aggregate(b)

插图

plot(r, xlim=c(-0.1,1.1), ylim=c(-0.1, 1.1))
lines(b, lty=2)
lines(aa, col="blue", lwd=2)
lines(bb, col="red", lwd=3)

如果您现在想为每个聚合缓冲区采样一个单元格

a <-disaggregate(aa)
set.seed(3) 
e <- extract(r, a, function(i) sample(na.omit(i),1))

==== 与您的新示例数据 ====

你也可以这样想

library(terra)

v <- c(0.72, 0.72, 0.72, 0.88, 0.99, 0.99, 0.97, 0.82, 0.9, NA, NA, NA, 0.72, 0.7, 0.73, 0.81, 0.83, 0.87, 0.86, 0.98,0.9,NA,NA,NA,0.71,0.71,0.75,0.78,0.83,0.81,0.92,0.98,0.79,NA,NA,NA,0.68,0.7,0.76,0.79,0.9,0.98,0.97,0.84,NA,NA,NA,0.99,0.69,0.73,0.77,0.86,0.99,0.97,0.99,0.82,0.67,NA,NA,0.76,0.76,0.79,0.79,0.77,0.82,0.81,0.87,0.74,0.92,NA,NA,1,0.76,0.77,0.65,0.61,0.77,0.71,0.63,0.69,0.86,NA,NA,1,0.72,0.66,0.59,0.71,0.66,0.66,0.62,0.88,1,NA,NA,NA,0.69,0.7,0.63,0.63,0.7,0.71,0.73,0.79,1,0.79,NA,NA,0.67,0.74,0.65,0.62,0.64,0.75,0.92,0.99,0.77,0.78,NA,NA,0.66,0.74,0.65,0.63,0.65,0.78,0.94,0.94,0.78,NA,NA,NA,0.72,0.74,0.68,0.67,0.72,0.79,0.87,0.99,0.98,NA,NA,NA)
r <- rast(nrows=12, ncols=12, xmin=-49, xmax=-48.5, ymin=-28, ymax=-27.5, vals=v)
x <- r > 0.8

x <- ifel(r > 0.8, 1, NA)
p <- disaggregate(as.polygons(x))
b <- buffer(p, .01)
plot(x)
lines(b)