Cómo obtener automáticamente el modo (modus) de un ráster en ArcGIS

Question

¿Tienes experiencia en cómo obtener automáticamente el valor de la moda de un ráster (valor que aparece con más frecuencia en un conjunto de datos)? Estoy buscando algo similar a "obtener propiedades del ráster" en el análisis espacial en ArcGIS, desafortunadamente no hay posibilidad de calcular la moda directamente. Posteriormente, me gustaría utilizar este cálculo de la moda en el Constructor de Modelos. Con base en el "valor de la moda" me gustaría reclasificar la imagen del ráster en un archivo binario (0= valores < moda, 1 = valores > moda).

Answer 1

Este puede que no sea el método más eficiente ni preciso, pero probablemente sea la forma más sencilla con Model Builder, siempre y cuando tu ráster sea categórico:

1) Utiliza estadísticas zonales (mayoría) con una capa de zona con una sola zona y la extensión de tu capa de ráster.

2) Utiliza la calculadora de ráster (o Con) para construir tu capa binaria ( Con( "raster" < "mayoría", 0,1) )

EDITAR: Ten en cuenta que necesitarás reclasificar tu ráster a un entero si originalmente está en punto flotante (este agrupamiento necesario es por qué dije que puede que no sea el método más preciso). Hasta donde sé, no podrás calcular la moda de una distribución continua con las herramientas disponibles en Model Builder, pero puedes extraer el histograma (usando histograma zonal), luego aplicar el método de @Jeffrey Evans. En cualquier caso, si tienes un gran número de píxeles, la mayoría es suficiente para la mayoría de los propósitos.

Answer 2

@Jeffrey Evans gracias por tu inspiración en R, descubrí un enfoque diferente y realmente fácil:

library("raster")
library("rgdal")

ndvi2011<-raster("n_ndvi2011.img") # leer raster .img
plot(ndvi2011)                     # mostrarlo
hist(ndvi2011)                     # encontrar un histograma

r<-getValues(ndvi2011) # para obtener el valor del píxel del raster
                       # lo convierte en un vector, para poder encontrar el valor modal

# valores necesarios
mode<-modal(r, na.rm=TRUE)
min<-minValue(ndvi2011)
max<-maxValue(ndvi2011)

# reclassificar raster

# crear matriz (ejemplo de matriz necesaria)
     #    min   mode    0
     #    mode  max 1

m<-c(min, mode, 0, mode, max, 1)
m<-matrix(m, ncol=3, byrow=T)

# raster reclasificado
rec_ndvi11<-reclassify(ndvi2011, m, right = T)

¡voilà! :)

Answer 3

Este es un problema de vectorización y sería muy costoso computacionalmente para un vector del tamaño de un ráster. Para hacer esto en ArcGIS, tendrías que utilizar NumPy. Debido al tamaño inherente de un ráster (1 millón de valores es un ráster bastante pequeño en estos días), sin aplicar un enfoque de submuestreo, esto no parece ser un problema muy tratable. Te animo a explorar un enfoque de umbral alternativo.

Si los datos son unimodales, pero sesgados, la mediana es comparable con la moda. Sin embargo, cuando la distribución presenta múltiples modas que no se pueden suavizar, entonces la mediana se vuelve sesgada, pero sigue siendo más representativa de la tendencia central que la media.

Derivo el número de modas y el valor de la moda (dominante) utilizando un spline o una función de densidad kernel gaussiana en R. No me interesa programar en Python, pero quizás alguien podría adaptar este código para abordar tu problema.

# Este ejemplo se puede ejecutar en base R.
    # ENCONTRAR NÚMERO DE MODAS
    nmodes <- function(x) {
       den <- density(x, kernel=c("gaussian"))
       den.s <- smooth.spline(den$x, den$y, all.knots=TRUE, spar=0.8)
         s.0 <- predict(den.s, den.s$x, deriv=0)
           s.1 <- predict(den.s, den.s$x, deriv=1)
         s.derv <- data.frame(s0=s.0$y, s1=s.1$y)
         nmodes <- length(rle(den.sign <- sign(s.derv$s1))$values)/2
           if ((nmodes > 10) == TRUE) { nmodes <- 10 }
         if (is.na(nmodes) == TRUE) { nmodes <- 0 }
       return(nmodes)
    }

   # VALOR DE LA MODA
    dmode <- function(x) {
      den <- density(x, kernel=c("gaussian"))
        ( den$x[den$y==max(den$y)] )
    }

   # Ejemplo con vector aleatorio
    x <- runif(1000,0,10000)
       dmode(x)   # valor de la moda
       nmodes(x)  # Número de modas

    # Mostrar resultados con valor de la mediana
     plot(density(x))
       abline(v=dmode(x), col="red")
       abline(v=median(x), col="black")
         legend("bottomleft",legend=c("modo","mediana"),
            lty=c(1,1),col=c("red","black"), bg="white")

Si tus datos son enteros, son basados en frecuencia y simplemente puedes utilizar el bin con la frecuencia máxima (conteo). Ahora, si analizamos la discretización de la distribución continua para hacer de esto un problema basado en frecuencia manejable, podemos aproximar la moda con la estrategia de discretización correcta. El algoritmo "Sturges" no siempre produce resultados deseables. Sin embargo, un método algo estable para el tamaño óptimo del bin (número de bins) es la "regla de Freedman-Diaconis" calculada como; n=(máx-min)/h donde; h=2∗RIQ∗n−1/3.

Aquí podemos observar la sensibilidad de un enfoque de discretización.

    # Discretización de Freedman-Diaconis
    ( fd.hist <- hist(x,breaks="FD") )

    # Valor de la moda basada en frecuencia discreta
    ( fd.freq.mode <- fd.hist$breaks[which.max(fd.hist$counts)] )

    # valor de la moda de la distribución continua
    ( x.mode <- dmode(x) )

    # Error residual de la moda discreta
    x.mode - fd.freq.mode

# Mostrar resultados con moda, mediana y moda de frecuencia
 plot(density(x))
   abline(v=dmode(x), col="red")
   abline(v=median(x), col="black")
   abline(v=fd.freq.mode, col="blue")
legend("bottomleft",legend=c("modo","mediana","freq"),
       lty=c(1,1,1),col=c("red","black","blue"), bg="white")