Take the 2-minute tour ×
Geographic Information Systems Stack Exchange is a question and answer site for cartographers, geographers and GIS professionals. It's 100% free, no registration required.

do you have some experience how to automatically get mode value from a raster ( value that appears most often in a dataset)? I am looking for something like "get raster properties" in Spatial analysis in ArcGIS, unfortunately there is no possibility to directly calculate mode. Subsequently, I would like to use this calculation of mode in Model Builder. Based on "mode value" I would like to reclassify raster image into binary file (0= values < mode, 1 = values > mode).

share|improve this question

3 Answers 3

This may not be the most efficient nor the most accurate way, but this is probably the easiest method with model builder, assuming that your raster is categorical:

1) use zonal statistics (majority) with a zone layer with a single zone and the extent of your raster layer.

2) use raster calculator (or Con) to build your binary layer ( Con( "raster" < "majority", 0,1) )

EDIT : Note that you will need to reclassify your raster to integer if it is originally in float (this necessary binning is why I said that it may not be the most accurate method). From what I know, you will not be able to compute the mode of a continuous distribution with the tools available in model builder, but you can exttract the histogram (using zonal histogram), then apply @Jeffrey Evans method. In any case, if you have a very large number of pixels, the majority is good enough for most purposes.

share|improve this answer
    
+1. Is there still no way to read and process a grid's VAT in Model Builder? Incidentally, the conversion from floating point to integer for this purpose is problematic; I suspect that concern is what Dr. Evans may have been trying to address in his answer. I would maintain this question almost requires the assumption that the input is a categorical (integer) grid. –  whuber May 12 at 18:46
    
the conversion is problematic, but you can multiply by a large number before the conversion in order to increase the precision. There is no VAT is the type is not integer, and from what I know you can't use the VAT directly in Model builder (you can in arcpy, using searchCursor) –  radouxju May 12 at 18:57
    
The problematic part is not technical, it's conceptual. The mode depends on the bin size, so what to do? One way is to view the data as a random sample of a larger population (which is usually not true of raster data, which is one problem). A "mode" is considered a mode of the population and is estimated by fitting a smooth density function to the data. That can be done in many ways that produce radically different modes (another problem), so some method to identify an optimal half-width (a sophisticated analog of the bin width) is needed (a third problem). –  whuber May 12 at 19:29
    
If the data is frequency based already then the solution is quite easy and defined elsewhere in the post. For conversion of continuous data, as stated by @whuber, the mode could be significantly influenced by the bandwidth. A common method for optimal bin size (n) is the "Freedman-Diaconis rule" calculated as; n=(max-min)/h where; h=2∗IQR∗n−1/3. However, this will not necessarily return the exact mode of the continuous distribution. I expanded my answer to demonstrate this. –  Jeffrey Evans May 12 at 20:24
    
@Jeffrey The difficulty is that this particular rule is designed to make effective histograms--not to estimate a mode accurately--and most other rules and procedures to estimate bandwidths are intended to achieve unrelated aims such as minimizing bias in estimated moments. But a bigger concern is that very few floating point raster datasets can legitimately be viewed as random samples of a population. But I think this is all moot: the "raster image" referenced in the question almost surely is categorical and can readily be processed as such. –  whuber May 12 at 20:29
up vote 1 down vote accepted

@Jeffrey Evans thank for your inspiration in R, I found out a different approach and really easy:

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

ndvi2011<-raster("n_ndvi2011.img") # read raster .img
plot(ndvi2011)                     # display it
hist(ndvi2011)                     # find a histogram

r<-getValues(ndvi2011) # to get pixel value from raster 
                       # it changes it to a vector, so it can find out modal value

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

# reclassify raster

# create matrix (exemple of needed matrix)
     #    min   mode    0
     #    mode  max 1

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

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

voila! :)

share|improve this answer

This is a vectorization problem and would be very computationally expensive for a vector the size of a raster. To do this in ArcGIS you would have to step out into NumPy. Because of the inherent size of a raster (1 million values is a quite small raster these days), without applying a subsampling approach, this does not seem like a very tractable problem. I encourage you to explore an alternative thresholding approach.

If the data is unimodal, but skewed, the median is comparable to the mode. However, when the distribution exhibits multiple modes, that cannot be smoothed out, then the median becomes biased, but is still more representative of the central tendency than the mean.

I derive the number of modes and value of the (dominate) mode using a spline or a Gaussian kernel density function in R. I do not care for programming in Python but, perhaps somebody could adapt this code to address your problem.

# This example can be run in base R.     
    # FIND NUMBER OF MODES 
    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)
    }

   # VALUE OF MODE
    dmode <- function(x) {
      den <- density(x, kernel=c("gaussian"))
        ( den$x[den$y==max(den$y)] )   
    }  

   # Example with random vector
    x <- runif(1000,0,10000)
       dmode(x)   # value of mode
       nmodes(x)  # Number of modes

    # Plot results with median value 
     plot(density(x))
       abline(v=dmode(x), col="red")
       abline(v=median(x), col="black")
         legend("bottomleft",legend=c("mode","median"),
            lty=c(1,1),col=c("red","black"), bg="white")

If your data is integer it is frequency based and you can just use the bin with the maximum frequency (count). Now if we look at binning the continuous distribution to make this a tractable frequency based problem we can approximate the mode with the correct binning strategy. The "Sturges" algorithm does not always produce desirable results. However, A somewhat stable method for optimal bin size (number of bins) is the "Freedman-Diaconis rule" calculated as; n=(max-min)/h where; h=2∗IQR∗n−1/3.

Here we can observe the sensitivity of a binning approach.

    # Freedman-Diaconis binning
    ( fd.hist <- hist(x,breaks="FD") )

    # Value of binned frequency-based mode
    ( fd.freq.mode <- fd.hist$breaks[which.max(fd.hist$counts)] ) 

    # value of continuous distribution mode
    ( x.mode <- dmode(x) )

    # Residual error of binned mode
    x.mode - fd.freq.mode

# Plot results with mode, median and frequency mode 
 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("mode","median","freq"),
       lty=c(1,1,1),col=c("red","black","blue"), bg="white")
share|improve this answer
    
I do not understand how this is a vectorization problem. Outside of Model Builder, the mode can be found almost instantaneously by sorting the dataset's value attribute table and identifying which value(s) have the greatest frequency. Your code seems to understand "mode" in a different sense as a kind of local maximum of a density of something--but it's not clear how that density would be related to the original data. Are you perhaps assuming the input would be a floating point raster? –  whuber May 12 at 18:44
1  
Yes, my assumption is that the raster is float. I should have thought of the case of an integer and proposed a frequency based solution using the vat's count. The post did not specify and in thinking about a distribution and "mode" I think in terms of a continuous vector of values. –  Jeffrey Evans May 12 at 20:06

Your Answer

 
discard

By posting your answer, you agree to the privacy policy and terms of service.

Not the answer you're looking for? Browse other questions tagged or ask your own question.