Monday, September 6, 2010

Detecting hot spot patterns in point incidents data with SAGA GIS

SAGA GIS has a number of grid analysis commands that can be used to detect patterns (i.e. find hot spots) in point data e.g. accident incidents. An example incident distribution on a map is shown in the figure below.

The point vector data may be storing additional attributes about the incident such as the number of occurrences.

Just by looking at the vector point data, it is impossible to detect a pattern or hot spot incidents in the area. By forming a graphical occurrence distribution surface from the vector point data with a grid analysis tool like SAGA GIS, it would be so much easier to visualize the areas of high incident density. An example of this surface is shown below, overlaid with other map features for reference.



Below are the steps I took to generate the grid surface for visualizing the incident patterns.

Import vector point data

  1. Start up SAGA GIS
  2. Select Modules | File | GDAL/OGR | OGR: Import Vector Data.

    The OGR: Import Vector Data dialog box appears.
  3. Click the Files row. Click the [...] button on the right. Browse and select the point vector data file e.g. Incidents.shp. Click Open.

    The selected file is displayed in the OGR: Import Vector Data dialog box.

  4. Click Okay.

    The point data is imported into SAGA GIS.

Convert the vector point data to a grid layer
  1. Select Modules | Grid | Gridding | Shapes to Grid.

    The Shapes to Grid dialog box appears.
  2. Click the Shapes drop down list and select the previously imported vector points e.g. 01.Incidents.
  3. Click the Attribute drop down list and select the field to get the new grid cell values, e.g. OCCURRENCE.
  4. Click the Preferred Target Grid Type drop down list and select a suitable data type e.g. Integer (2 byte).


  5. Click Okay.

    The User Defined Grid dialog box appears.
  6. In the Cellsize field, type in the new cell size e.g. 30.


  7. Click Okay.

    The vector points are rasterized into a new grid layer.
Workaround the Aggregate module's Sum No Data cell value bug
SAGA GIS 2.0.5's Aggregate module has a problem summing grid cells with negative values, i.e. cells with no data and usually represented with a large negative number e.g. -99999. In order for the SUM function to work, those large negative values have to be set to 0. This can be done by using the Calculus module as described below.
  1. In the Workspace pane, select the newly created point grid layer e.g. 01.Incidents [OCCURRENCE].


  2. In the Object Properties pane on the right, click the Settings tab. Then change the No Data field values to "0; 0".


  3. Click Apply.
  4. Select Modules | Grid | Calculus | Grid Calculator.

    The Grid Calculator dialog box appears.

  5. In the Grid system drop down list, select the newly created grid layer's system, e.g. "30;354x 420y;...".
  6. Select the Grids row. Click the [...] button

    The Grids dialog box appears.

  7. On the list on the left side, select the newly created grid layer, e.g. 01.Incidents[OCCURRENCE]. Click the >> button.

    The selected grid layer is moved to the list on the right.
  8. Click Okay.

    The Grids dialog box is closed.
  9. In the Formula field, type in ifelse(lt(a,0),0,a).

    Note: this formula tells the Calculator to replace all grid cell values below 0 with the value 0.

  10. Click Okay

    A new grid layer is created with no negative values.

Using Aggregate to sum adjacent grid cells
  1. Select Modules | Grid | Construction | Aggregate.

    The Aggregate dialog box appears.
  2. In the Grid system drop down list, choose the previously created grid layer's system, e.g. "30;354x 420y;...".
  3. In the Grid drop down list, select the newly created grid layer without negative values, e.g. ifelse(lt(a,0),0,a).
  4. In the Aggregation Size field, type in an appropriate size (in cells), e.g. 20. In the Method drop down list, choose Sum.


  5. Click Okay.

    The aggregated grid layer is created.

    This is how the aggregated layer looks like in the map view.

Resample the aggregated grid layer to match the original point grid layer.
  1. Select Modules | Grid | Construction | Resampling.

    The Resampling dialog box appears.
  2. In the Grid system drop down list, select the aggregated grid layer's system, e.g. "600;17x21y...".
  3. In the Grid drop down list, choose the aggregated grid layer, e.g. "ifelse(lt(a,0),0,a)".
  4. In the Target Grid drop down list, choose grid.


  5. Click Okay.

    The Choose Grid dialog box appears.
  6. In the Grid system drop down list, choose the original point grid layer's system, e.g. "30;354x420y...".


  7. Click Okay.

    The Down-Scaling dialog box appears.

  8. Click Okay.

    The aggregated grid layer is resampled to match the original point grid layer's dimensions.

No comments: