在 R 中的栅格堆栈上找到第二高的值 [英] Find second highest value on a raster stack in R

问题描述

在 R 中，我可以使用 max/min 命令轻松计算地理参考栅格堆栈中每个像元的最大值/最小值.

In R I can easily compute the max/min value of each cell in a georeferenced raster stack using the max/min commands.

set.seed(42)
require(raster)
r1 <- raster(nrows=10, ncols=10)
r2=r3=r4=r1
r1[]= runif(ncell(r1))
r2[]= runif(ncell(r1))+0.2
r3[]= runif(ncell(r1))-0.2
r4[]= runif(ncell(r1))
rs=stack(r1,r2,r3,r4)
plot(rs)
max(rs)
min(rs)

然而，我一直试图找到一种方法来找到堆栈中第二高的值.就我而言，堆栈上的每个栅格都表示特定模型跨空间的性能.我想比较第一个和第二个最佳值，以确定从亚军开始的最佳模型有多好，而不必将我的堆栈转换为矩阵，然后再转换回栅格.有什么想法或建议吗??

However, I have been trying to find a way to find the second highest values across a stack. In my case, each raster on the stack denotes performance of a particular model across space. I would like to compare the first vs second best values to determine how much better is the best model from its runner up without having to convert my stack to a matrix and then back into a raster. Any ideas or suggestions??

推荐答案

您可能想要使用 calc()，根据您的具体情况调整下面的代码.只是为了表明它如宣传的那样工作，我分别绘制了通过在 4 层 RasterStack 对象的每个单元格中找到的最高、第二高、第三和第四高值形成的层.

You'll probably want to use calc(), adapting the code below to your precise situation. Just to show that it works as advertised, I've separately plotted layers formed by taking the highest, second highest, third, and fourth highest values found in each cell of the 4-layer RasterStack object.

zz <- range(cellStats(rs, range))

par(mfcol=c(2,2))
plot(calc(rs, fun=function(X,na.rm) X[order(X,decreasing=T)[1]]), main="1st",zlim=zz)
plot(calc(rs, fun=function(X,na.rm) X[order(X,decreasing=T)[2]]), main="2nd",zlim=zz)
plot(calc(rs, fun=function(X,na.rm) X[order(X,decreasing=T)[3]]), main="3rd",zlim=zz)
plot(calc(rs, fun=function(X,na.rm) X[order(X,decreasing=T)[4]]), main="4th",zlim=zz)

或者，更紧凑、更高效，只需构建一个新的栅格堆栈来保存重新排序的值，然后绘制其图层:

Or, more compactly and efficiently, just construct a new raster stack holding the reordered values and then plot its layers:

zz <- range(cellStats(rs, range))
rs_ord <- calc(rs, fun=function(X,na.rm) X[order(X,decreasing=T)])

par(mfcol=c(2,2))
plot(rs_ord[[1]], main="1st", zlim=zz)
plot(rs_ord[[2]], main="2nd", zlim=zz)
plot(rs_ord[[3]], main="3rd", zlim=zz)
plot(rs_ord[[4]], main="4th", zlim=zz)

这篇关于在 R 中的栅格堆栈上找到第二高的值的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持IT屋！