Population Density Map
Population Density Map
M. Bahrami
Wolfram Research, Inc.
Population density maps are among the most valuable datasets, e.g. when investigating the epidemic outbreak of an infectious disease. The Humanitarian Data Exchange (HDX) seems to have world’s most detailed population density datasets. These high-resolution maps, created by Facebook, estimate the number of people living within 30-meter grid tiles. In this short computational essay, I will discuss how to covert these maps into a computable datasets in the Wolfram Language, with desirable size and shape of grid tiles. An an example, I shall focus on Qatar and show how to get desired maps using GeoHistogram, a very powerful built-in function in the Wolfram Language.
The HDX high resolution population density maps can be obtained through their website. The dataset is provided in in two formats, GeoTIFF and CSV, for distribution of various populations. I will download and process the CSV file of the overall population density. The same methods can be applied for any other desired population dataset. Since Qatar’s dataset is small, I shall import it directly. But for larger datasets, I suggest downloading the file and then importing it.
Import zipped data from HDX and get the CSV file:
In[]:=
qatar=Import,"population_qat_2019-07-01.csv";
Return the number of tiles on the original data:
In[]:=
Length[qatar]
Out[]=
307646
Return the headings of the original data:
In[]:=
qatar〚1〛
Out[]=
{Lat,Lon,Population}
Generate a list of two sublists as geo positions and a corresponding population:
In[]:=
{geo,ppl}={qatar〚2;;,;;2〛,qatar〚2;;,-1〛};
Generate a list whose elements are {geoposition->population}:
In[]:=
data=Thread[(GeoPosition/@geo)ppl];
Generate a density histogram:
In[]:=
GeoHistogram[data,PlotStyleEdgeForm[Gray],GeoScaleBar"Metric",PlotLegendsAutomatic,GeoRange]
Out[]=
Specify the diameter of each tile:
In[]:=
Table[GeoHistogram[data,i*,PlotStyleEdgeForm[Gray],GeoScaleBar"Metric",PlotLegendsAutomatic,GeoRange,PlotLabel"Grid Tile ("<>ToString[i]<>"km)"],{i,10,25,5}]
Out[]=
In[]:=
Table[GeoHistogram[data,i*^2,PlotStyleEdgeForm[Gray],GeoScaleBar"Metric",PlotLegendsAutomatic,GeoRange,PlotLabel"Grid Tile ("<>ToString[i]<>")"],{i,100,800,200}]
2
km
Out[]=
One can also set three different tile shapes: “Hexagon”, “Rectangle” and “Triangle”.
In[]:=
Table[GeoHistogram[data,{i,250*^2},PlotStyleEdgeForm[Gray],GeoScaleBar"Metric",PlotLegendsAutomatic,GeoRange,PlotLabel"Tile Shape: "<>i],{i,{"Hexagon","Rectangle","Triangle"}}]
Out[]=
Therefore, given a tile shape and size, we need to extract the geo positions of tiles generated by GeoHistogram and their corresponding values.
In[]:=
gh=GeoHistogram[data,{"Hexagon",250*^2}];
If you look at the full form of the GeoHistogram output, you will notice that the grid tiles are wrapped in Tooltip.
In[]:=
Short[polygons=Cases[gh,Tooltip[Polygon[y_?MatrixQ],x_]Polygon[GeoPosition[y]]x,∞]]
Out[]//Short=
PolygonGeoPosition1.28428×,PolygonGeoPosition174644.,63,PolygonGeoPosition102.292
6
10
However, the actual positions of polygon points are not correct.
In[]:=
GeoGraphicsEdgeForm[Red],polygons[[All,1]],EdgeForm[Black],Polygon
Out[]=
A closer look shows that the polygon points in the GeoHistogram output are also projected using GeoGridPosition. Accordingly, in order to extract polygons with correct geo positions, we need to include the projection using GeoGridPosition too.
In[]:=
polygonCounts=Cases[gh,Tooltip[y:Polygon[_?MatrixQ],x_]:>GeoPosition[GeoGridPosition[y,(GeoProjection/.Options[gh,GeoProjection])]]x,∞];
In[]:=
GeoGraphicsEdgeForm[Red],polygonCounts[[All,1]],EdgeForm[Black],Polygon
Given this result, we now put the data in the format of a list with sublists as {polygon points, center of polygon, population density} and export it as CSV file (might be useful for users of other programming languages).
Summary
Summary
Given the data from HDX, generate grid tiles of population density for a given shape and size of tiles:
Visualize each tile with a color specified by cfunc, defined before:
Acknowledgment
Acknowledgment
I’d like to thank my colleagues, Alan Joyce (Wolfram|Alpha team) and José Martín-García (Algorithms R&D team), for their valuable input in helping with datasets and important suggestions on code.