22. January 2014
In some applications you want to calculate the density of points. It sounds very easy and in fact it is using QGIS. Let me show you how to create a point density raster. Especially, let me show you how to do this with the heatmap plugin in QGIS.
You will need the common installation of QGIS 2.0.1 Dufour and an installed heatmap plugin (see how to install a plugin here). We will use some field data from Indonesia with 500.000 data points. You may download the data (12MB) we will use here.
Fire up QGIS and set the CRS of the project according to the used CRS in your data. In our case it is WGS84 Zone 49S which is EPSG:32749 and use the Layer->Add Delimited Text Layer function to add the downloaded test data to your project:
I’ve chosen to create a spatial index as it improves the handling of this more or less bigger dataset. After this you need to determine the CRS of the data itself:
Now your project should look like this:
Creating a point density map is the purpose of the heatmap function/plugin. As you might have this installed prior you should follow Raster->Heatmap->Heatmap. We need to fill the parameters for the heatmap:
I’ve changed the standard parameters a little. The colums Output raster and Output format should be clear. Avoid some special ooutput formats and try to stick with standard GeoTIFF as it is recommended by the author Arun. Radius is the radius to search for points around a location in m or map units (in our case it’s meters as we are using a projected Coordinate Reference System (CRS) with a base unit of m). As we are using the Quartic (biweight) kernel it defines the direct distance to the point itself. I am using 564m as the radius so I can easily state that the values in the raster describe the density for a square kilometer. By using a different kernel this could change as the form of the kernel changes. The Quartic (biweight) seems like a circle kernel.
The rows and columns as well as the cell sizes determine the spatial resolution of the target raster. The finer the raster the longer it takes to compute the raster. Keep this in mind. There are also more parameters to adjust but you may want to press the help button to get some insight in the functionality of the plugin itself.
The resulting map is an image going from black (means small number of points around) to white (large number of points around the given location). We will change this to be more informative. First open the properties of the raster go to the histogram tab and compute the histogram. You will see that the curves reaches till approx. 13 on the x-axis which means that the highest values are about 13 points per km²:
We will use this value in our Singleband pseudocolor representation. Therefore switch over to the Style tab and change the Render Type to Singleband pseudocolor. I am using the color ramp YlOrRd and set the Max value to 13. After this click on Classify and on Ok:
By using the identify results button you can now check the value at each location:
You may ask, why the values are not integer. According to the filter we have used, points near a certain location have bigger influence than points far away a given location. So the points are weighted in their influence. This is useful as we might have an underlying uncertainty in the location of points. But if we have several points near each other it should be more certain that the point locations are reliable.
For visual impression I like to see some indicators where the values change. Therefore I am using the Contour function in the Raster dialog. The contour creation function is part of the GdalTools￼ plugin by Faunalia. Well the name Contour is a little critical. Basically this function takes a raster image and creates lines of identical values. In our case it is the point density. Go to Raster->Extraction->Contour:
The interval is set to one and I also would like to have a new column in the line shape file where the value of the current line is stored. The results are better to interpret afterwards. You can easily see the slight right orientation of the data: