[GIS] How does Zonal Statistics work exactly

layersqgisrastersagazonal statistics

I'm a bit confused on how the zonal statistics tools work in QGIS (either the 'Zonal Statistics' tool under the Raster tools from the QGIS geoalgorithms toolbox or the 'Raster statistics from polygons' tool under the SAGA Vector<->Raster tools).

Basically I have a raster map layer and some small polygons in a separate layer. I want to assign each polygon the value of the raster below. In case a polygon overlaps multiple raster cells, I just want to have the maximum.

Yet the output doesn't make any sense and I can't get behind the reasoning of those results (which are obviously not the maximum of the underlain raster cell):

Any ideas what I'm doing wrong or what I didn't consider yet? The results with the python API in qgis.analysis.QgsZonalStatistics() (according to this question here) are much better, but it does only offer count, mean and sum…but for some reason no max or min.

Best Answer

Update as April 2018 and QGIS 2.18

The current version of Zonal Statistics Plugin can provide:

Count: to count the number of pixels
Sum: to sum the pixel values
Mean: to get the mean of pixel values
Median: to get the median of pixel values
StDev: to get the standard deviation of pixel values
Min: to get the minimum of pixel values
Max: to get the maximum of pixel values
Range: to get the range (max - min) of pixel values
Minority: to get the less represented pixel value
Majority: to get the most represented pixel value
Variety: to count the number of distinct pixel values

Original answer

In short, this approach with QGIS, as far as I'm aware is not possible to view the raster pixel maximum as a result. Only the count, mean and sum. I've provided some resources that covers this.

Here is the doc page for the plugin Zonal Statistics Plugin

With the icon_zonal_statistics Zonal statistics plugin, you can analyze the results of a thematic classification. It allows you to calculate several values of the pixels of a raster layer with the help of a polygonal vector layer (see figure_zonal_statistics). You can calculate the sum, the mean value and the total count of the pixels that are within a polygon. The plugin generates output columns in the vector layer with a user-defined prefix.

Also, this QGIS Tutorial Sampling Raster Data using Points or Polygons provides a sample data to achieve zonal statistics return. I've included screenshot they provide as a example of their workflow and results.

Enable the Zonal Statistics Plugins. This is a core plugin so it is already installed. See Using Plugins to know to how enable core plugin

Once the processing finishes, select the counties layer. Use the Identify tool and click on any county polygon. You will see three new attributes added to the layer: ZS_count, ZS_mean and ZS_sum. These attributes contain the count of raster pixels, mean of raster pixel values and sum of raster pixel values respectively. Since we are interested in average temperature, the ZS_mean field will be the one to use.

Related Solutions

[GIS] How to calculate raster statistics for polygons

I was struggling to do exactly the same thing, but for various reasons I'm committed to using QGIS. I tried using v.rast.stats using the GRASS plugin and also via the Sextante plugin. The latter approach failed, because it seems to attach the stats to a temporary vector layer which it then deletes. The GRASS plugin worked, but it doesn't deal with overlapping polygons.

After doing some digging around (in the source of the promising-sounding ZonalStats plugin), I found that QGIS actually has zonal statistics methods built into the API, and these also have Python bindings. So as long as you are only looking for count, sum and mean statistics for your polygon features, the Python Console (Plugins > Python Console) is currently the easiest way to attach the stats to the polygon attribute table.

Select your raster layer in the TOC, and type the following into the console (it grabs the source file name of your raster layer)

>>> rasterfile = qgis.utils.iface.mapCanvas().currentLayer().source()
Select your vector layer, and execute the following command in the console (it grabs the vector layer itself)

>>> vectorlayer = qgis.utils.iface.mapCanvas().currentLayer()
Execute the following three commands in the console (they pass the vector layer and raster file to QGIS' built-in zonal stats calculator)

>>> import qgis.analysis

>>> zonalstats = qgis.analysis.QgsZonalStatistics(vectorlayer,rasterfile)

>>> zonalstats.calculateStatistics(None)

The results will be appended as extra fields in the polygon layer.

Zonal Statistics

Note that if you want to append a prefix to the new fields (rather than just having them as sum, mean and count), then you would replace the first command in step 3 with the following (assuming you wanted to use 'zonal_' as the prefix:

>>> zonalstats = qgis.analysis.QgsZonalStatistics(vectorlayer,rasterfile,"zonal_")

[GIS] Zonal Statistics for multiple rasters

I have something fairly close, although it only currently works for average and rejects max and min.

Save all your raster layers in one folder with the same extension (e.g. .tif)
Load shapefile (polygons) into QGIS

Open Python (Plugins > Python console) and run the following

import glob, os, qgis.analysis
vectorlayer = qgis.utils.iface.mapCanvas().currentLayer()
rasterfolder = 'C:/Users/…/Folder/'
os.chdir(rasterfolder)
for lyr in glob.glob("*.tif"):
    qgis.analysis.QgsZonalStatistics(vectorlayer, lyr, attributePrefix=lyr, rasterBand=1).calculateStatistics(None)

Using the vector layer loaded in QGIS was an easy way to append the results to the attribute table. The code to loop through rasters in a folder was from here.

If there is a way to specify the statistic using QgsZonalStatistics.Mean that would be helpful to know.

Best Answer

Related Solutions

[GIS] How to calculate raster statistics for polygons

[GIS] Zonal Statistics for multiple rasters

Related Question