提取栅格值(从堆栈)到for循环中的点 [英] Extract raster values (from Stack) to points in for loop
问题描述
我有一个栅格堆栈和100分.对于每个栅格,我要提取值,并使用三个不同的比例尺/缓冲区进行提取.
I have a raster stack and 100 points. For each raster I want to extract the value and do so using three different scales/buffers.
首先,这是将三个栅格合并为一个堆栈
First, here are three rasters combined into a stack
library(raster)
# Make rasters and combine into stack
set.seed(123)
r1 = raster(ncol=1000, nrow=1000, xmn=0, xmx=1000, ymn=0, ymx=1000)
values(r1) = round(runif(ncell(r1),1,100))
r2 = raster(ncol=1000, nrow=1000, xmn=0, xmx=1000, ymn=0, ymx=1000)
values(r2) = round(seq(1:ncell(r1)))
r3 = raster(ncol=1000, nrow=1000, xmn=0, xmx=1000, ymn=0, ymx=1000)
values(r3) = round(runif(ncell(r1),1,5))
RasterStack <- stack(r1, r2, r3)
然后我生成100个点作为 SpatialPoints
对象
I then generate 100 points as a SpatialPoints
object
#make points
Points <- SpatialPoints(data.frame(xPoints = sample(1:1000, 100),
yPoints = sample(1:1000, 100)))
接下来,我定义要循环通过的三个缓冲区
Next, I define the three buffers that I want to loop through
Scales <- c(60, 500)
为了更好地描述期望的结果,我将首先仅使用单个栅格,而不使用RasterStack.下面的代码定义了一个矩阵(输出),该矩阵在循环中填充,每一列都是在两个不同的 Scales
上提取的 r1
的值.然后将这些列标记在循环外部.
To better describe the desired outcome, I will first use only a single raster, not the RasterStack. The code below defines a matrix (output) which is populated in the loop with each column being the extracted values of r1
at the two different Scales
. The columns are then labeled outside of the loop.
output <- matrix(ncol = length(Scales), nrow = length(Points))
for( i in 1:length(Scales)) {
output[, i] <- extract(r1, Points, method='simple', buffer=Scales[i], fun=mean)
}
colnames(output) <- paste("r1", Scales, sep = "_" )
> head(output)
r1_60 r1_500
[1,] 50.67339 50.42280
[2,] 50.42401 50.42335
[3,] 49.96709 50.44288
[4,] 50.65492 50.52634
[5,] 50.60678 50.43535
[6,] 50.52477 50.48277
我想要相同的输出,但是我不想对单个栅格(例如上面的r1)进行调用,而是要对 RasterStack
中的每个栅格执行此操作.最终结果将是一个矩阵(或data.frame),每个栅格都有两列(r1:r3).如示例中所示,标记将与相应的比例相对应,以便将列标记为 r1_60,r1_500,r2_60,...,r3_500.
I want this same output, but rather than calling a single raster (e.g. r1 above), I want to do this for each raster in the RasterStack
. The final result would be a matrix (or data.frame) that has two columns for each raster (r1:r3). As in the example, labeling would correspond to the respective scale so that the columns were labeled r1_60, r1_500, r2_60, ... , r3_500.
我认为嵌套的 for
循环可以在通过 RasterStack
和通过 Scales
循环的地方工作,但我怀疑可能会有更好的方法
I think a nested for
loop would work where I loop through the RasterStack
and through the Scales
but suspect there might be a better way.
对于真实数据,我有20个栅格,这些栅格分别是1293年和1541年的30,000个位置.我也有5种不同的音阶,因此嵌套的 for
循环将花费很长时间.
For the real data I have 20 rasters that are 1541 by 1293 and around 30,000 locations. I also have 5 different scales so a nested for
loop will take a very long time to run.
添加采用不同的方法,我可以使用以下代码创建数据帧列表,每个数据帧对应于使用给定缓冲区提取的每一层的值.
Addition Taking a different approach, I can use the following code to create a list of data frames, each of which corresponds to the extracted values of each layer using a given buffer.
output <- list()
for(i in 1:length(Scales)){
output[[i]] <- extract(RasterStack, Points, method='simple', buffer = Scales[i], fun = mean)
names(output)[[i]] <- paste("Buffer", Scales[i], sep = "_")
}
从此输出中,如何制作单个6 x 100数据帧,其中每一列都将被标记为"layer_buffer number".例如,layer.1_60,layer.2_60,...,layer.2_500,layer.3_500.
From this output, how can I make a single 6 by 100 data frame where each column would be labeled as the "layer_buffer number". For example, layer.1_60, layer.2_60, ... , layer.2_500, layer.3_500.
我还可以发布一个首选问题.
I can also post a new question of preferred.
推荐答案
为方便起见,我发布了最适合我的解决方案.鉴于栅格数据包错误,我没有使用0缓冲区将值提取到点.
For the sake of closure, I am posting the solution that worked best for me. In light of the raster package bug, I did not extract values to points using the 0 buffer.
Scales <- c(60, 500)
然后,使用前10点,
Points <- Points[1:10]
我使用以下代码为每个缓冲区级别创建了一个列表.
I created a list for each buffer level using the following code.
output <- list()
for(i in 1:length(Scales)){
output[[i]] <- extract(RasterStack, Points, method='simple', buffer = Scales[i], fun = mean)
names(output)[[i]] <- paste("Buffer", Scales[i], sep = "_")
}
Then, following the post linked here, I used the following code to combine the list of data frames into a single data frame.
do.call(cbind,lapply(names(output),function(x){
res <- output[[x]]
colnames(res) <- paste(colnames(res),x,sep="_")
res
}))
返回的df的 head
在下面.
layer.1_Buffer_60 layer.2_Buffer_60 layer.3_Buffer_60 layer.1_Buffer_500
[1,] 50.67339 408657.5 3.013623 50.42280
[2,] 50.42401 449786.5 2.990888 50.42335
[3,] 49.96709 968829.9 2.995279 50.44288
[4,] 50.65492 119448.9 3.009086 50.52634
[5,] 50.60678 141819.5 2.998585 50.43535
[6,] 50.52477 394303.5 2.984253 50.48277
layer.2_Buffer_500 layer.3_Buffer_500
[1,] 435485.7 2.999983
[2,] 460519.9 2.999632
[3,] 775273.5 3.002715
[4,] 273116.8 3.000364
[5,] 289803.0 2.999054
[6,] 426887.0 3.000055
这篇关于提取栅格值(从堆栈)到for循环中的点的文章就介绍到这了,希望我们推荐的答案对大家有所帮助,也希望大家多多支持IT屋!