【问题标题】:How to calculate the most abundant pixel value in a layer stack of different classifications of the same area (majority vote)?如何计算同一区域不同分类的层堆栈中最丰富的像素值(多数票)?
【发布时间】:2019-09-28 08:54:42
【问题描述】:

我有多个分类为栅格,每一个都属于一个不同的组,它们具有代表某些类别的离散整数像素值(例如:111 代表稳定的针叶林)。我想找出在所有分类中每个像素最丰富的值,我将其堆叠到栅格堆栈中。

有这样的功能吗?如果没有任何想法如何解决此任务?

library(raster)
r1 <- raster(nrows = 4, ncols = 4, xmn = 0)
r1[] <- sample(1:4, 16, replace = T)

r2 <- raster(nrows = 4, ncols = 4, xmn = 0)
r2[] <- sample(1:4, 16, replace = T)

r3 <- raster(nrows = 4, ncols = 4, xmn = 0)
r3[] <- sample(1:4, 16, replace = T)

r <- stack(r1, r2 ,r3)
r[2,2]

例如 1,2,3 或 4 类在 Pixel [2,2] 中最丰富? 最后,我希望有一个具有最丰富像素值的栅格层和一个具有该值的分数的栅格层(留在上面的示例中,值 1、1 和 4,因此值 1 将以66%)。

【问题讨论】:

    标签: r spatial raster r-raster


    【解决方案1】:

    如果您在示例中使用随机值,则应修复随机种子以使其可重现,仅供参考。

    这应该可以满足您的需要:

    library(raster)
    
    set.seed(42)
    
    r1 <- raster(nrows = 4, ncols = 4, xmn = 0)
    r1[] <- sample(1:4, 16, replace = T)
    
    r2 <- raster(nrows = 4, ncols = 4, xmn = 0)
    r2[] <- sample(1:4, 16, replace = T)
    
    r3 <- raster(nrows = 4, ncols = 4, xmn = 0)
    r3[] <- sample(1:4, 16, replace = T)
    
    r <- stack(r1, r2 ,r3)
    
    f <- function(x){
    
      n <- length(x)
    
      if (length(unique(x)) == n){
    
        return(c(NA,NA))
    
      }else{
    
        v <- as.integer(names(which.max(table(x))))
    
        p <- (sum(x == v) / n) * 100
    
        return(c(v,p))
    
      }
    }
    

    我们创建了一个函数,该函数将向量作为输入,计算多数值并返回它在向量长度中所占的百分比。然后将该函数传递给calc 以应用于堆栈r。所以你会得到一个 2 波段的栅格。

    如果没有重复值,那么如果唯一值的长度等于向量的长度,则返回NA

    不确定当您有两个具有相同多数值的类时会发生什么,可能该函数只采用其中一个您也可以为此添加检查。

    所以如果我们现在查看像素[2,2],我们可以查看结果:

    值为1,分数为66%(与其他类似)

    > r[2,2]
         layer.1 layer.2 layer.3
    [1,]       3       1       1
    
    plot(calc(x = r, fun = f2))
    

    【讨论】:

      【解决方案2】:

      您可以使用raster::modal 方法。它有一个参数ties,您可能需要考虑。

      library(raster)
      set.seed(42)
      r <- raster(nrows = 4, ncols = 4, xmn = 0)
      s <- stack(lapply(1:12, function(i) setValues(r, sample(1:4, 16, replace = TRUE))))
      
      x <- modal(s)
      plot(x)
      

      【讨论】:

        猜你喜欢
        • 2021-10-30
        • 2017-09-16
        • 2022-01-14
        • 1970-01-01
        • 2021-12-19
        • 2017-11-18
        • 1970-01-01
        • 2021-05-26
        • 1970-01-01
        相关资源
        最近更新 更多