5.2. Tutorial 4¶

The following is a tutorial about the functions of the Semi-Automatic Classification Plugin (SCP). It is assumed that you have a basic knowledge of QGIS.

Tutorial 4: Postprocessing tools for assessing land cover change

5.2.1. Tutorial 4: Postprocessing tools for assessing land cover change ¶

This tutorial aims to analyze land cover change using SCP Postprocessing tools. Basically, we are going to assess land cover change from two raster classifications, and relate the changes to a land use vector file. An overview of several postprocessing tools is also provided.

The following is the video tutorial, and the following text illustrates the phases in detail.

http://www.youtube.com/watch?v=0IUosyr4pRw

The tools can be applied to any land cover classification, but we are going to use Copernicus data, which are freely available (as established by the EU Regulation No 1159/2013 of 12 July 2013) and cover the European countries. Of course, this tutorial is designed for demonstration purposes and it is not endorsed by the European Union. The original Copernicus data (produced with funding by the European Union) are downloaded from https://land.copernicus.eu/ and remain the sole property of the European Union.

Following, a brief description of the data we are going to use.

The Copernicus High Resolution Layers are raster classifications with 20m spatial resolution. Several land cover classes are available, but in this tutorial we are going to use the Imperviousness Density for 2012 and 2015. These data classify the degree of imperviousness (0-100% of impermeable cover of soil), which is the artificially sealed area. The Imperviousness Density was produced using automatic derivation based on calibrated Normalized Difference Vegetation Index. You can find the detailed product specifications here.

The Copernicus Corine Land Cover is a land use/land cover vector produced by standard methodology of photo-interpretation of satellite images. The vector is classified in 44 classes divided in 3 hierarchical levels with minimum mapping unit of 25 hectares. In this tutorial we are considering only the first level of Corine Land Cover 2012, divided in these classes:

artificial surfaces
agricultural areas
forests and semi-natural areas
wetlands
water bodies

5.2.1.1. Refine the classifications with direct editing ¶

You can download the data for this tutorial from this archive , or use your own data (two classification rasters and a land use vector).

For this tutorial, the original Copernicus data were modified by clipping the rasters to a small area over Florence (Italy).

Start QGIS and load the two rasters IMD_2012.tif and IMD_2015.tif that are Copernicus Imperviousness Density for 2012 and 2015 respectively. As you can see, the rasters have values from 0 to 100, representing the degree of imperviousness.

It is useful to refine the classification by photo-interpretation, especially for data produced by semi-automatic processing.

We can use high resolution images or other services such as OpenStreetMap. For example you can follow this tutorial Download the Data to download satellite images, or you can download a subset of a Landsat 8 image, already converted to reflectance, from this link (about 27 MB, data available from the U.S. Geological Survey), unzip the downloaded file, and load the bands in QGIS.

First, we need to define a Band set containing a classification raster (this is required for drawing ROIs).

Open the tab Band set clicking the button in the Меню SCP or the Панель SCP. Click the button to refresh the layer list, and select the IMD_2012 raster (just this raster is sufficient); then click to add selected raster to the Band set 1.

Band set definition

Optionally, we can create a band set for the satellite image to display a color composite; open the tab Band set and select all the Landsat bands in the list; click to add a new band set, then click to add selected rasters to the Band set 2.

Band set of the Landsat image

In QGIS zoom to an area where we want to correct the classification. In this case we are going to manually remove a few pixels pretending they are classification errors.

We need to manually create a ROI, but first check that the Band set 1 is active. Now click the button in the Робоча панель. Left click on the map to define the ROI vertices and right click to define the last vertex closing the polygon. An orange semi-transparent polygon is displayed over the image, which is a temporary polygon (in this case we don’t need to define the Training input).

Manual ROI polygon

Now open the tool Edit raster opening the Меню SCP and the submenu Postprocessing . Select the Input raster, for instance IMD_2012. According to the legend of Imperviousness Density, in Use constant value enter 100 (we want to correct impervious pixels; in case we would like to correct not impervious pixels we would enter the value 0). The other options are fine. Therefore, click RUN to edit the raster.

Attention: the input raster is directly edited; it is recommended to create a backup copy of the input raster before using this tool in order to prevent data loss.

The raster modified

Of course we could repeat these steps to edit any area of the raster.

TIP : Sometimes changes are not immediately visibile because the raster is not refreshed; try to zoom out and zoom in to refresh the view.

5.2.1.2. Classification report ¶

It could be interesting to know the area of each land cover class. In order to get the area statistics, open the Меню SCP and click the tab Classification report under the submenu Postprocessing .

Click the button to refresh the layer list, and select the IMD_2012 raster in Select the classification ; next click RUN to start the calculation; the output report is saved in a text file and displayed in the tab Output.

Report tool

We can repeat the same steps for the IMD_2015 raster. Over the 86% of the area is not impervious.

The classification report

5.2.1.3. Reclassification ¶

Before calculating land cover change it is convenient to reclassify the imperviousness degree into two classes: built-up and not built-up. A possible threshold for the distinction between built-up and not built-up is 30% (for further information read this document ). We can reclassify the raster using the SCP tool, obtaining the simple classification 1 = built-up and 0 = not built-up.

Open the tool Reclassification. In Select the classification select the raster IMD_2012. Click the button twice to add two rows to the table. We need to enter the expressions illustrated in the following table.

Reclassification table

Old value	New value
raster < 30	0
raster >= 30	1

Uncheck the options Use code from Signature list and click RUN to start the reclassification. A new raster will be created (e.g. BU_2012).

Reclassification tool

Now select the IMD_2015 (the reclassification table is the same as before) and click RUN to reclassify the 2015 raster (e.g. BU_2015). Now the two reclassified rasters are loaded in the map and we can assing an appropriate symbology.

Reclassified rasters

5.2.1.4. Remove isolated pixels ¶

We are going to compute the land cover change, but first we may want to remove isolated pixels in order to improve the analysis. In fact, single pixels may not represent real changes between the two classifications, because of geometrical shifts or isolated classification errors. Of course, this step is not always required, and it should be avoided if the purpose of the analysis is to find also the smallest changes.

We are going to use classification sieve for removing single pixels. Open the tool Classification sieve.

In Size threshold leave 2; all patches smaller the the selected number of pixels (i.e. single pixels) will be replaced by the value of the largest neighbour patch. Of course we could increase this value if we want to remove larger patches.

In Select the classification select the raster BU_2012. The option 4 in Pixel connection determines how pixels are considered connected, that is in a 3x3 window diagonal pixels are not considered connected. If we select the option 8 also diagonal pixels are considered connected.

Tool classification sieve

Now click RUN to create the new raster BU_2012_sieve.

Of course, we should repeat these steps also for raster BU_2015 to create the new raster BU_2015_sieve.

The rasters after removing isolated pixels

5.2.1.5. Assess land cover change ¶

Now we can use the tool to assess land cover change between the two classifications 2012 and 2015. Open the tool Land cover change.

This tool is quite straightforward. Click the button to refresh the layer list. In Select the reference classification select the BU_2012_sieve raster, that is the first classification. In Select the new classification select the BU_2015_sieve raster, that is the latest classification.

Uncheck the option Report unchanged pixels, because we want to report only the pixels where the classification changed between 2012 and 2015. Now click RUN to create the new land cover change raster (e.g. change). Also, a text file is created (i.e. a file .csv separated by tab) containing the land cover change statistics.

Land cover change tool

The values of the land cover change raster represent a combination between reference and new classification, as described in the text file. In this case, only the value 1 is present that is the condition where BU_2012_sieve = 0 and BU_2015_sieve = 1.

Land cover change raster

From the report we ca read that 520 pixels changed from 0 to 1, while no pixel changed from 1 to 0 between years 2012 and 2015.

Land cover change report

5.2.1.6. Analyze the context of land cover changes ¶

Now, it could be interesting to compare land cover change to other data such as land use, in order to analyze the context of new built-up areas. We are going to cross the land cover changes to the vector of Corine Land Cover; this way we can differentiate the new built-up areas according to Corine Land Cover classification system.

The original Corine Land Cover data were modified by clipping to a small area over Florence (Italy) and adding a field Class_1 filled with the first level of classification.

Load in QGIS the Copernicus Corine Land Cover shapefile CLC_2012.shp previously downloaded. You can see the symbology of the first level Corine Land Cover classes that are:

artificial surfaces
agricultural areas
forests and semi-natural areas
wetlands
water bodies

A subset of Corine Land Cover

Open the tool Cross classification. Click the button to refresh the layer list. In Select the classification select the change raster, that is our land cover change. Check Use NoData value and set the value 0, in order to exclude unchanged pixels (having value 0 in the change raster) from the analysis.

In Select the reference vector or raster select the vector CLC_2012 and in Vector field select the field Class_1, containing the code of first level classes.

Cross classification tool

Now click RUN to create a new raster of comparison (e.g. change_CLC). The output will report the area of each combination between change code and CLC_2012 code.

From the cross matrix we can evaluate the area in \(m^2\) of built-up changes occurrend in the 5 classes of Corine Land Cover classification.

Cross matrix

		CLC_2012
		1	2	3	4	5
Change	1	157600	48400	2000	0	0

The tool Cross classification can be very useful also for other analyses that involve the comparison with other data, such as population or flood risk, but this could be the subject of other tutorials.

5.2.1.7. Assess the spectral signature of changes ¶

An optional step could be the assessment of the spectral signature of changes. We can download satellite images (see Tutorial 1: Your First Land Cover Classification for the details) and calculate spectral signatures for monitoring the changes through time.

We are going to use the Landsat 8 image downloaded at the beginning of this tutorial for calculating the spectral signature of changes. First, we need to create a Training input to store the spectral signatures calculated from the classes. In the Панель SCP select the tab Входові навчальні дані and click the button to create the Training input (define a name such as signatures.scp). The path of the file is displayed and a vector is added to QGIS layers with the same name as the Training input (in order to prevent data loss, you should not edit this layer using QGIS functions).

Now open the tool Class signature opening the Меню SCP and the submenu Postprocessing . In Select the classification select the raster change_CLC, thus we can distinguish the spectral signatures of changes. In Select input band set enter the number of the band set containing the Landsat 8 bands (i.e. 2). Now click RUN to start the calculation.

The tool for extracting spectral signatures from classes

After a while the spectral signatures are loaded in the Training input.

If the changes involved vegetation, we could calculate spectral signatures for images acquired in different seasons and assess the phenological variations of vegetation through spectral signatures. Also, these spectral signatures could be used as training input for further land cover classifications.

In order to display the signature plot, in the Перелік сигнатур ROI highlight two or more spectral signatures (with click in the table), then click the button . The Графік спектральних сигнатур is displayed in a new window.

The calculated spectral signatures

5.2.1.8. Export the changes to vector format ¶

This is an optional step that may be useful for further analyses and integration with other data. We are going to convert the change raster to vector.

Open the tool Classification to vector. In Select the classification select the change_CLC raster and uncheck the Use code from Signature list. Now click RUN to create a new vector (e.g. change_vector).

Raster to vector tool

In the attribute table of this change_vector you can see the field C_ID that represents the code of the change raster as described in Assess land cover change. Of course we could delete the polygons with code 0 (unchanged area), displaying only changes with code 1.

The vector of changes over the Landsat image