Training Material for UN Open GIS OpenData

Introduction

Scope

The following educational material has been drafted within the framework of the OSGeo UN Committee Educational Challenge - Open Geospatial Data and software for UN sustainable development goals. The overarching goal is to show that at this time, the combination of open (geo)data globally available and the significant developments of the free and open source solutions for geospatial is sufficient to initiate geospatial analysis, at worldwide level, at small and intermediate scales, to better understand our ecosystem. In that respect, we have employed OSGeo software solutions to process global open geospatial datasets to answer one selected indicator for a sustainable development goal. The selected indicator is 9.1.1 Proportion of the rural population who live within 2 km of an all-season road (C0901010) which supports the target of developing quality, reliable, sustainable and resilient infrastructure, including regional and transborder infrastructure, to support economic development and human well-being, with a focus on affordable and equitable access for all. The indicator has been chosen after a close analysis of all SGDs and the corresponding indicators as to comply with the following:

1. to have a spatial dimension;
2. to not be an indicator that is already addressed through another initiative, such as the GEO Wetlands Initiative, WHO Interactive Air Pollution Maps, GEO AquaWatch, ESA CoastColour etc.;
3. if possible, to not be yet the subject of a published methodology.

Audience

We have prepared this educational material for researchers, educators and professionals in local, regional, national or international agencies with minimal to intermediate geospatial information knowledge. We assume our audience has already basic knowledge of geospatial data structures, formats and that they have already used a GIS software, as to have basic skills and understanding of how to work with geospatial and tabular data. In that respect, we have limited the interactions with the command line, however we have inserted references to it.

Data and software used

Datasets and software used For the calculation of the SDG indicator, we have only used QGIS 3.4. We have also taken into consideration that changes that might occur from one version to another and thus focused on the functions used more than on a step-by-step guide. The datasets used for our exercise are the following:

Acquired knowledge

After going through the entire educational material, one will be able to:

• Have a broader view on what are the types of geospatial data open at global scale, as well as what are their limitations.
• Have a more deeper understanding of working with geospatial data using a dedicated software
• Consistent knowledge of QGIS software fundamentals
• Learn how to create cartographic representations of the obtained results

Educational Material

Open geospatial information and its role in answering UN Sustainable Development Goals

Population dataset description
Global Administrative Units Dataset description
For road related data, we have decided to use OpenStreetMap data as it is the only homogeneously designed globally available dataset. Without doubt, the amount and the quality of the available data for various regions around the world can vary consistently. However, given the clear and consistent definition of each map element and tag, this exercise should be reproducible in any other part of the world.
Yet, given our area of interest, the Tabora county from Tanzania, we must take into consideration specific developments for Africa, more precisely, the Highway Tag Africa - Topology of Road Network in African countries, and furthermore, the East Africa Tagging Guidelines.
However, with consideration to the global replicability of our educational material, we will also insert specifications on a more general scale. Of course, it must be acknowledged that the workflow presented here could require other adjustments with respect to the specificity of the road dataset used in calculation.

Preparing the geospatial data

For the scope of this exercise we have chosen the Tabora county of Tanzania. As we strive to create an educational material that can be applied no matter the region of interest, a decision was made to use the available datasets, on a global level. The following table presents the datasets used:

Topic	Name collection/dataset	Abstract	Indicators	Produce/collector	Owner	License	Type of data	Format	Scale/spatial resolution	Edition	CRS	Other URL
Administrative units	Database of Global Administrative Areas	GADM provides maps and spatial data for all countries and their sub-divisions.	administrative units	University of California, Berkeley,Museum of Vertebrate Zoology, and theInternational Rice Research Institute (Global Administrative Areas 2009)	GDAM	The data are freely available for academic use and other non-commercial use. Redistribution, or commercial use is not allowed without prior permission.	vector	Geopackage, shapefile, geodatabase. KMZ, R formats	n/a	April 2018	Geographic WGS84	https://gadm.org/metadata.html
World Population	WorldPop	Alpha version 2010 and 2015 estimates of numbers of people per grid square, with national totals adjusted to match UN population division estimates (http://esa.un.org/wpp/) and remaining unadjusted.	Settelments, Population numbers, birth and pregnancy, age structures, poverty spatial distribution etc.	GeoData Institute, University of Southampton	GeoData Institute, University of Southampton	CC BY 4.0	raster	GeoTIFF	100 m	July 2013	Geographic WGS84	http://www.worldpop.org.uk/data/methods/
World Population	Global Rural-Urban Mapping project (GRUMP), v1	To provide a polygon representation of urban areas with city or agglomeration name and time series population estimates.	urban geometries	Socioeconomic Data and Applications Center (sedac)	Socioeconomic Data and Applications Center (sedac)	CC BY 4.0	vector	shapefile	30 arc-second	2006	Geographic WGS84	n/a
World Population	Global Human Built-up And Settlement Extent (HBASE) Dataset From Landsat, v1 (2010)	To provide high spatial resolution estimates of global urban extent derived from global 30m Landsat satellite data for the target year 2010 and a companion dataset to the Global Man-made Impervious Surface (GMIS) dataset.	urban extent	Socioeconomic Data and Applications Center (sedac)	Socioeconomic Data and Applications Center (sedac)	CC BY 4.0	raster	GeoTiff	30 m	2017	Geographic WGS84, UTM	n/a

Step 1

Add all listed data to your project. Create group layers as you bring in the data, so it is easier when you start processing to navigate through all datasets. Create and save your project, so you can pick up the work from where you left it.
/ [Layer]-->[Add layer] - will open the Data Source Manager that allows you to load the data.: the administrative unites (GDAM dataset), the population numbers (WorldPop - we will use the TZA_popmap15adj_v2b.tif file), the urban extent (we will use the global_urban_extent_polygons_v1.01.shp file).
As mentioned, we will use OpenStreetMap data for the roads geometry and condition. Bringing OSM data will require you install a new plugin - OSM Downloader.

[Plugins]-->[Manage and install plugins]. 

Step 2

The datasets are in various projections, either the Geographic projection EPSG 4326 or the Pseudo_Mercator EPSG 3857. A geographic coordinate system is based on a spheroid and uses angular units (degrees). Thus, when using QGIS calculator, for example, it returns values in decimal degree and not meters. You can see the used units in a projection's description that can be retrieved from epsg.io. As we will work with road geometries, we must reproject all the datasets in a projected coordinate system, which is based on a 2D plane (with the spheroid projected on a 2D plane) and uses linear units, such as meters. For our study, we identify a suitable CRS for our region of interest, the Tabora county in Tanzania. To do that we will use epsg.io. After a quick search, we find WGS 84 / UTM zone 36S-EPSG: 32736 to be the appropriate for our region.
To reproject vector data using QGIS, we have to save the file with the desired projection.
Click on the vector layer you want to reproject and choose [Export]-->[Save features as..]
For raster datasets, we will use gdalwrap that is available as a processing tool in the Processing toolbox.
[Processing]-->[Toolbox] We can search by typing the keyword 'reproject' in the search bar.

• Fig1 reproject.jpg

  Fig. 1 Reproject tool for a raster dataset

Then, we will cut all layers by the boundary of the selected county, Tabora.
Firstly, export from the administrative units level 1, Tabora county. Secondly, clip all layers by its geometry. fig2_clip.jpg|Clip tool identified in the Processing toolbox fig3_clip_rastere.jpg| gdalwarp utility allows cutting a raster by a vector extent </gallery>

Step 3

Step 3 produces the rural areas of the Tabora county. According to Wikipedia, Tanzania is divided into regions (GDAM administrative level 1), districts (GDAM administrative level 2) and wards, (GDAM administrative level 3).
We will calculate the Rural Access Index on wards, thus we will extract the rural regions from the administrative units level 3. The resulting dataset will be vector type. [Vector]-->[Geoprocessing tools]-->[Difference]

Step 4

Step 4 prepares the dataset from which we will extract the population number for all rural areas of each administrative unit. We will use zonal statistics, a tool now available in the Processing Toolbox.

• Fig4 zonalStats.jpg

  Zonal statistics ran on the administrative units - rural extent

Step 5

Bring in the roads!
Step 5 is the most time consuming processing stage and, more over, it may vary when this exercise will be applied to other regions in the world. For the roads geometry and condition we will use the OpenStreetMap data available.
As we will see, Tanzania is very well represented on the OSM map, even if we will encounter various situations mainly regarding quality of network connectivity, an important aspect for our study. That is because in 2015, the Crowd2Map Tanzania was launched and during the following years, there have been significant crowd mapping campaigns.
After step 1 and 2, the road dataset should be imported into QGIS, clipped, reprojected in EPSG 32736 and saved as a geopackage file.
Next, we will do a preliminary cleaning, by eliminating all roads segments that are not suitable for cars. The SDG indicator 9.1.1 Proportion of the rural population who live within 2 km of an all-season road refers to roads that are suitable for any kind of vehicle (average modern automobile), thus we will filtrate by: "highway" = 'cycleway' or "highway" = 'pedestrian' or "highway" = 'path' or "highway" = 'footway' or "highway" = ‘residential’ or "highway" = ‘service’. We have also deleted roads under construction, because that means that the roads can not be used for access.
OpenStreetMap defines map features by tags. Each tag has 2 elements: the key and the value. The key is used to describe a topic, category, or type of feature and the value describes the specific form of the key-specified feature. In our filter above, the tag "highway" = 'cycleway' is formed by the key "highway" and the value 'cycleway', which means that the map feature is a road for bicycles. More details on how OpenStreetMap is build can be found in the OSM wiki.
There is a number of tags considered in OSM road data that gives important indication related to the road condition. These are:

• Key:smoothness with 8 possible values: excellent, good, intermediate, bad, very_bad, horrible, very_horrible, impassable. According to the OSM community, in developing countries, the HDM model is using some values of smoothness to help define road quality. The HDM stands for Highway Development and Management (HDM-4) and it represents the World Bank tool for the management of roads network, particularly in developing countries.
• Key: surface
• Key:surface:grade
• Tag: tracktype with a gradual variation from grade1 - solid to grade5 - soft.

Considering the East Africa Tagging Guidelines, there is an emphasis in using surface=* with the generic separation of surface=paved for sealed roads and surface=unpaved for the others. Analyzing the available data for the area of interest, the Tabora county of Tanzania, we can identify that for a number of road segments, we have road quality particular tags identified. We can see that in the other_tags column. Looking at the structure of the attributes registered in this column, we observe a number of things:

• The tags (key+value) are separated by comma;
• The keys are separated by their values through "=>" ;
• Tags are not necessarily in the same order;
• Not all tags that are relevant for our exercise are registered.

To simplify the processing steps, we will break the other_tags column into one column per tag. We will do this using an external solution, LibreOffice. The file containing the attributes for each segment is the .dbf file (dBase file) which can be opened and manipulated by tabular data solutions as well. It is advisable though that any change on the file to be done with caution, as it can corrupt the file, and thus make it unusable.

-->open .dbf file using LibreOffice --> Text to column by separator comma

Next, we create one column for each road condition key available: surface, smoothness and tracktype with the values specific for each road segment. Afterwards, we populate the three new columns (surface, smoothness and tracktype) with the specific value for each road segment. For this, we will use the Field Calculator of QGIS introducing the following expression:

 CASE 
  WHEN  "other_tags" LIKE '%surface%' 
  THEN  regexp_replace ("other_tags", 'surface"=>"', )
  ELSE surface
 END

where,
Other_tags - column from where we take the attribute value, and will be subsequently replaced by columns Tag1, Tag2 and so on for each of the three columns.
Surface - is the key name we are looking for in the attribute list.

After we have populated the three new columns (surface, smoothness and tracktype), we can extract information that will allow us to better understand the roads dataset we will use to calculate the RAI. We will make use of the statistics functionality of QGIS to show the total length of roads in our area of interest, but also various lengths, as for example the lengths of paved roads, or bad roads lengths and so on.

To calculate the roads length, we open the attribute table of the dataset, we create a new column and we write $length.

• Fig5 length.jpg

  Calculating the roads length.

The next step leading to the RAI calculation is to determine the road condition that will be taken into account. This represents a difficult assessment to make, mainly due to lack of harmonization between OpenStreetMap defined tags and official internationally used methodologies for road condition calculations. For the purpose of our exercise, we will consider the proposed definition of what a “road in good condition” given in the World Bank 2016 document Measuring Rural Access: Using New Technologies (Transport & ICT. WB, 2016.):
"A road in good condition refers to:
* Paved road with IRI less than 6 meters/km and unpaved road with IRI less than 13 meters/km, when IRI data are available * Paved road in excellent, good, or fair condition and unpaved road in excellent or good condition, when IRI data are not available but other road condition data, such as the PCI or visual assessment by class value, are available.
In our case, there is no available IRI data available, thus our data falls in the second category, where we have a qualitative description: Paved road in excellent, good, or fair condition and unpaved road in excellent or good condition. In order to work with available OSM road data, we will chart a compatibility map between OSM tags that describe road condition and the given classification.
Analysing our data with GroupStats, we discovery the following situation:

• Fig6 roadCondition1.jpg

  Caption1

• Fig6 roadCondition2.jpg

  Caption2

GroupStats is a plugin for easy calculation of statistics for features in a vector layer. The plugin has a control panel that allows the user to make calculations (as in our example, length) for specific group features, as can be seen in attached file 1.

• Fig7 groupStats.jpg

  Example of calculation the length for roads, per categories

A look at the complete dataset tags shows that we can find any kind of combination of the tags we are interested in. We can observe that, for our test dataset, we have 100% the HIGHWAY tag completed for all roads segments.

File:File1 highway

A complete look at the dataset using GroupStats

Continuing the analysis, with consideration to the definitions given by the OSM when assigning the specific values, we will proceed in mapping the correspondence. And so, giving the tags descriptions and the definition of a road in “good or fair condition in rural areas”, we have considered the following mappings:

• Fig8 mappingTags.jpg

  Mapping OSM tags to international road condition definitions

Next, we will clean the roads dataset of the categories that we have excluded, regardless of the highway=* value. Using GroupStats, we measure that we will exclude from RAI calculation. The filter GroupStats must be done by "smoothness" = 'very_bad' OR "smoothness" = 'very_horrible' OR "smoothness" = 'impassable' OR "smoothness" = 'horrible' OR "tracktype"='grade4' OR "tracktype" = 'grade5' and it will find 3082 matching features, with a total length of 3359.96 KM.

• Fig9 excludedRoadsbyCondition.jpg

  GroupStats analysis on bad quality OSM road segments

The spatial distribution of the roads eliminated from our analysis is presented in the map below:

• 

  Excluded road condition categories

After a quick glance, we can see that a significant central rural area is served by the roads with a too poor condition to be included in the RAI calculation. Further on, we will work with a road layer that has been cleaned of the roads that are not in good condition. To do that, select in the attribute table using the following expression ["smoothness" = 'very_bad' OR "smoothness" = 'very_horrible' OR "smoothness" = 'impassable' OR "smoothness" = 'horrible' OR "tracktype"='grade4' OR "tracktype" = 'grade5'] and then delete the selected features.
Running GroupStats, we discover the following situation for Tabora county, with respect to road condition accepted in the RAI calculations:

• Fig10 OSMRoads RAI.jpg

  The total length of road segments in good or fair condition

The next step in preparing the road dataset resides in solving the irregularities related to the road network available. With consideration to the RAI calculation, to have a correct topological road network is desirable, but not a condition. As OpenStreetMap is an open collaborative mapping project, there is a high possibility that datasets collected have all sorts of inconsistency. In the following processing steps, we are aiming at identifying and correcting these inconsistencies as to eliminate as much as possible artificial results for RAI. One situation to resolve is when we have segments with shorter lengths that are further away from any other connected road by at least 2 km. Why would this situation yield artificial values for RAI? Because even if rural population has access to such a road in good conditions it offers no real connectivity, it gets them only so far. Before proceeding, let us analyse the current situation: identify the percentage and spatial distribution of segments shorter than 2 km. As we can see in the attribute table, the road dataset for our exercise has 21811 features. After filtering, we identify 19709 segments that are shorter than 2km.