Advanced Raster Analysis (QGIS3)

In the previous tutorial Basic Raster Styling and Analysis (QGIS3), you learnt about performing raster algebra with Raster Calculator. This tutorial builds on these techniques and shows you how to use other raster analysis tools from the Processing Toolbox. You will learn how to process with Land Use Land Cover (LULC) rasters in QGIS to extract certain types of landcover classes and map changes.

Overview of the task

We will use the South African National Land Cover dataset to identify and extract informal settlements in the City of Johannesburg, South Africa. We will also use a change assessment dataset to identify urban growth patterns in the city from 2014 to 2018.

Other skills you will learn

  • How to reproject raster data to another projection.

  • How to load an excel file in QGIS.

  • How to adjust the transperancy of a raster layer in QGIS.

Get the data

We will download the following datasets for this tutorial

  1. The South African National Land Cover 2018 dataset: The new South African National Land-Cover 2018 dataset has been generated from 20 meter multi-seasonal Sentinel 2 satellite imagery and contains 73 landcover classes.

  2. The South African National Land Cover 2018 Change Assessments: This dataset compares the change in 20 landcover classes from 2014 to 2018.

  3. COJ Boundary: A boundary shapefile for the City of Johannesburg, South Africa.

The Environmental Geographical Information Systems (E-GIS) provides access to environmental geospatial data for South Africa. We will download the South African National Land-Cover (SANLC) rasters from this portal.

  1. Visit the GIS Data Download page. Click I agree to accept the Conditions of Use and proceed.

../../_images/data1.png
  1. You will need to create a free account to download the dataset. Click I want to create an account and follow the instructions to create an account.

../../_images/data2.png
  1. Once logged in, search for South African National Land Cover (SANLC) 2018 Computer Automated Land Cover (CALC). This dataset is provided in 2 different projections. For this tutorial, we will downnload the ALBERS dataset. Click to it to download the SA_NLC_2018_Albers_CALC_data.zip file.

../../_images/data3.png
  1. Next, search for New South African National Land Cover (SANLC) 2014 and 2018 Change Assessment Datasets Computer Automated Land Cover (CALC) and click on the SA_NLC_2014_2018_CLASS_CHANGE_CALC (DATASET AND REPORT) to download the SA_NLC_2014_2018_CLASS_CHANGE_ALBERS_CALC.tif.vat.zip file.

../../_images/data4.png

The City of Johannesburg publishes spatial datasets as part of the Spatial Development Framework 2040 (SDF) for Johannesburg data. We will download the boundary shapefile from this site.

  1. Click the download link http://bit.ly/joburg-sdf-16.

../../_images/data5.png
  1. Click on the SDF Shapefiles directory.

../../_images/data6.png
  1. Download the SDF Shapefiles.zip file, and unzip it to a folder.

../../_images/data7.png

For your convenience, a clipped version of the required layers is available from the link below:

SA_NLC_2018_ALBERS_CALC.tif

SA_NLC_2014_2018_CLASS_CHANGE_ALBERS_CALC.tif

COJ_Boundary.zip

Data Source: [SANLC] [COJ]

Procedure

  1. Browse to the downloaded and unzipped folder in the browser. Expand it and drag and drop the SA_NLC_2018_ALBERS_CALC.tif in canvas.

../../_images/011.png
  1. Once the layer is loaded, you can notice the CRS will be set as Unknown CRS on the bottom right. Double-click on it to open the Project Properties - CRS dialog box.

../../_images/021.png
  1. At the bottom, you will see a preview of the projection extent. This Unknown CRS is a custom Lambert Equal-Area Projection defined for the country of South Africa. We will later reproject this layer to another projection. Click OK.

../../_images/031.png
  1. Load other two layers SA_NLC_2014_2018_CLASS_CHANGE_ALBERS_CALC.tif and COJ_Boundary. You will see that the raster layer cover the entire country. For our analysis, we are only interested in the area covered by the COJ_Boundary layer. We will now clip the raster layer to this region. Go to Processing ‣ Toolbox ‣ GDAL ‣ Raster extraction ‣ Clip raster by mask layer tool. Double-click to open it.

../../_images/041.png
  1. In the Clip Raster by Mask Layer, select SA_NLC_2018_ALBERS_CALC as the Input layer, then COJ_Boundary as Mask layer. We also have an option to reproject the data to another projection. It is a good practice to keep all your data layers in the same projection. We will reproject the rasters to match the CRS to that of the COJ_Boundary layer. Select EPSG:4326 - WGS 84 as the Target CRS.

../../_images/051.png
  1. The default output data format is GeoTiff. GeoTiff files can get very large if they are not compressed. A good practice is to always apply a loss-less compression when creating new raster layers. Expand Advanced Parameters and choose High Compression as the Profile. Next, click the ... button next to Clipped (mask) and select Save to file… to enter layer name as SA_NLC_2018_Clipped. Click Run.

../../_images/061.png
  1. Once the algorithm finishes, do not close the window. We will apply the same operation to the other raster layer. Switch to the Parameters tab and change the Input layer to SA_NLC_2014_2018_CLASS_CHANGE_ALBERS_CALC. Keep all other options but change the output layer name to SA_NLC_2014_2018_CLASS_CHANGE_Clipped. Click Run.

../../_images/071.png
  1. Both clipped layers will now be loaded in canvas. Select the original layer and click Remove Layer to remove them.

../../_images/081.png
  1. All the three remaining layers are now in the same CRS. We can switch the project CRS to the CRS of the layers now. Right click on any clipped layer and choose Layer CRS ‣ Set Project CRS from Layer.

../../_images/091.png
  1. Now the project CRS will be set to EPSG:4326. Bring the SA_NLC_2018_Clipped layer to top.

../../_images/101.png
  1. Click on SA_NLC_2018_Clipped and use the identify tool in the Attributes Toolbar to click on the image and inspect the pixel values. You will see that the pixel values range from 1 to 73. These values represent a distinct land use/land cover class.

../../_images/111.png
  1. The dataset classes are described in the SANLC 2018 Presentation, which can be downloaded from the EGIS Portal. For this exercise, we are interested in the informal settlements represented by class numbers 51 through 54.

../../_images/121.png
  1. Let’s extract pixels belonging to these classes. Go to Processing ‣ Toolbox ‣ Raster analysis ‣ Raster calculator tool. Double-click to open it.

../../_images/13.png
  1. The source image has only 1 band. The @1 suffix indicates the band number. Enter the following expression to select pixels from class 51-54.

"SA_NLC_2018_Clipped@1" >= 51 AND "SA_NLC_2018_Clipped@1" <= 54
../../_images/14.png
  1. Scroll down and click the ... button next to Reference layer(s). Select the SA_NLC_2018_Clipped layer and click OK.

../../_images/15.png
  1. Next, click the ... button next to Output and select Save to File….

../../_images/16.png
  1. Name the output file residential_informal.tif and click Run.

../../_images/17.png
  1. Once the process finishes, a new layer, residential_informal will be added to QGIS. This raster layer has only two-pixel values - 1 where our expression evaluated true and 0 where it was false. The pixels that appear white are the ones belonging to the informal settlement classes. We will style this layer better so we can see the informal settlements clearly. Click the Open the layer styling panel button.

../../_images/18.png
  1. Select the residential_informal layer and change the renderer to be Paletted/Unique values. Click the Add values manually (+) button.

../../_images/19.png
  1. Change the Value to 1 and enter Residential Informal as Label. Select a color of your choice.

../../_images/20.png
  1. We can now see all the informal settlements in the city of Johannesburg. It would be helpful to see them in context with a base map. We have access to a variety of base maps from the QuickMapServices plugin. Once you install the plugin, go to Web ‣ QuickMapServices ‣ OSM ‣ OSM Standard to add the OpenStreetMap layer.

../../_images/21.png
  1. Now you can easily identify and verify whether our analysis correctly identified the informal settlements. You can select the residential_informal layer and switch to the Transparency tab in the Layer styling panel. You can reduce the Global Opacity to see both the extracted pixels and the base-map together.

../../_images/22.png
  1. You have now completed the first part of the tutorial. Now we will use the SA_NLC_2014_2018_CHANGE_Clipped raster layer to identify regions that were urbanized between 2014 and 2018. Turn off all layers except SA_NLC_2014_2018_CHANGE_Clipped, then click the Open the layer styling panel button. Switch to the Transparency tab and enter 0 in Additional no data value. This will set the pixels with value 0 to transparent.

../../_images/23.png
  1. Use the Identify tool in the Attributes Toolbar to click on the image and inspect the pixel values. You will see that the pixel values range from 21-420. Each value indicates a transition from one of the 73 source classes to another class.

../../_images/24.png
  1. Your data download comes with a spreadsheet named lcccodes.xlsx. This sheet has a sheet 03 urban_change_codes that gives more details about each pixel value. We are interested in all pixel values where any 2014 class changed into a 2018 built-up class. In the image below, these are highlighted in blue.

../../_images/25.png
  1. Our goal is to map changes in the built-up class. We will apply a transformation on the SA_NLC_2014_2018_CHANGE_Clipped layer so all the pixel values are mapped from their original values to either of the following values.

1

All pixels which were a built-up class in both 2014 and 2018

2

All pixels which changed from a non built-up class in 2014 to a built-up class in 2018.

0

All remaining pixels

  1. To do this, we need to create a table specifying these rules. As QGIS is able to read spreadsheets directly, it is the most convenient method to create this table. Our spreadsheet should have 3 columns, MIN, MAX, and OUTPUT. Each row should be the range of input raster values that should be assigned an output value. Create a spreadsheet as shown below and save it to your computer as reclass.xlsx. You may also download a ready-to-use copy from this link - reclass.xlsx

../../_images/27.png
  1. Locate the reclass.xlsx file in the browser. Drag-and-drop it to the main window.

../../_images/28.png
  1. A new layer Sheet1 will be added to the Layers panel. Right-click on it and select Open Attribute Table. Verify that the sheet was imported correctly and you have 3 columns named MIN, MAX and OUTPUT. Open the Processing Toolbox ‣ Reclassify by layer tool.

../../_images/29.png
  1. In the Reclassify by layer dialog, select SA_NLC_2014_2018_CHANGE_Clipped as the Raster layer. Select Sheet1 as the Layer containing class breaks. Select MIN, MAX and OUTPUT fields for their respective fields.

../../_images/30.png
  1. Expand the Advanced Parameters section. Change the Range boundaries to min <= value <= max. Click the button for Reclassified raster and enter the output file name as builtup_change.tif. Click Run.

../../_images/31.png
  1. Once the processing finishes, a new layer builtup_change with pixel values 0-2 will be added to the canvas. In the Layer styling panel, choose Paletted/Unique values, then click Classify.

../../_images/32.png
  1. Choose the color of your choice for each category and label the 0, 1 , and 2 pixel values as Non Built-up, Existing Built-up and New Built-up.

../../_images/33.png
  1. Now in the Transparency tab, reduce the Global Opacity, and turn on the OSM Standard layer to see both the builtup_change pixels and the base-map together.

../../_images/34.png

If you want to give feedback or share your experience with this tutorial, please comment below. (requires GitHub account)