Setup
import wget
import rasterio
from rasterio import mask
import matplotlib.pyplot as plt
import numpy as np
import geopandas as gpd
import pandas as pd
import folium
Download raster images
url2020 = "ftp://ftp.worldpop.org.uk/GIS/Population/Global_2000_2020/2020/UGA/uga_ppp_2020.tif"
#wget.download(url2020)
url2019 = "ftp://ftp.worldpop.org.uk/GIS/Population/Global_2000_2020/2019/UGA/uga_ppp_2019.tif"
#wget.download(url2019)
Open and read raster images
raster_uga = rasterio.open("uga_ppp_2020.tif")
pop_uga_data = raster_uga.read(1)
Sum pixel values
pop_uga_data_pixel_sum = pop_uga_data[pop_uga_data > 0].sum()
Plot raster image
plt.imshow(pop_uga_data);
# Plot using log scale(+1)
plt.imshow(np.log1p(pop_uga_data));
A nicer plot
def plot_raster(rast_data, title='', figsize=(10,10)):
"""
Plots population count in log scale(+1)
"""
plt.figure(figsize = figsize)
im1 = plt.imshow(np.log1p(rast_data),) # vmin=0, vmax=2.1)
plt.title("{}".format(title), fontdict = {'fontsize': 20})
plt.axis('off')
plt.colorbar(im1, fraction=0.03)
title = 'Population Distribution (2020) in Uganda (Log Scaled) \n Est count: {}'.format(pop_uga_data_pixel_sum)
plot_raster(pop_uga_data, title)
Aggregate for a sub-polygon
uga_gdf = gpd.read_file("gadm36_UGA_2.zip")
uga_gdf
uga_gdf['GID_1'].unique()
uga_gdf.plot()
plt.title('Uganda: level 2 regions');
uga_gdf.crs
gtraster, bound = mask.mask(raster_uga, uga_gdf[uga_gdf.NAME_1 == 'Apac'].geometry, crop=True)
gtraster[0][gtraster[0]>0].sum()
For each year, aggregate for all sub-polygons
uga_gdf = gpd.GeoDataFrame.from_file("gadm36_UGA_2.shp")
for year in range(2019, 2021):
# Open and read raster file
raster_uga = 'uga_ppp_{}.tif'.format(year)
pop_raster_uga = rasterio.open(raster_uga)
pop_uga_data = pop_raster_uga.read(1)
# Loop through each polygon contained in the shapefile and use it as the mask to extract values
_results = []
for i in uga_gdf['GID_1']:
roi = uga_gdf[uga_gdf.GID_1 == i]
# using the mask.mask module from Rasterio to specify the ROI
gtraster, bound = mask.mask(pop_raster_uga, roi["geometry"], crop=True)
# values greater than 0 represent the estimated population count for that pixel
_results.append(gtraster[0][gtraster[0]>0].sum())
# Save the estimated counts for each year in a new column
uga_gdf[str(year)] = _results
uga_gdf['growth_rate'] = uga_gdf[['2019', '2020']].pct_change(axis=1)['2020']*100
uga_gdf
Choropleth map
fig, ax = plt.subplots(figsize=(12,8))
uga_gdf.plot(column="growth_rate", scheme='NaturalBreaks', k=5, cmap='coolwarm', legend=True, ax=ax, legend_kwds={'bbox_to_anchor':(0.25, 0.92)})
plt.title('Spatial distribution of population growth: Five natural breaks')
plt.tight_layout()
ax.axis("off")
#plt.savefig('myMap.png',dpi=150, bbox_inches='tight')
plt.show()