Network and Interconnection in Python Maps

A hands-on tutorial using Python packages including geopandas, networkx, and matplotlib Basemap

Mar 9, 2022

9 min read

Data and statistics related to demography and development indicators are at the heart of creating robust policies for development. Moreover, visual and spatial representations of these indicators provide crucial insights that are pertinent to provide rigour and transparency in policy or investment decisions. With this knowledge in mind, I was interested in visualising the spatial distribution of population and population density across my home country, Nepal, and how the geographical regions have been connected via the domestic aviation sector.

There are several packages available not only for data analysis, but also for constructing maps, geospatial analysis, and geocoding in Python. Some of the main advantages of using Python for geospatial analysis include its user-friendly syntax, potential support from a large community (e.g. via GitHub, Stack Overflow, etc.), and open access- one does not need to pay for expensive licenses.

Python packages such as geopandas, networkx, and matplotlib Basemap are useful for visualising network and interconnection such as transportation routes, pipelines, and power system interconnectors. In the first section of this post, I am going to walk you through how I constructed maps showing the spatial distribution of population and population density of Nepal. In the second section, I am going to describe how I plotted the domestic airports in Nepal and the flight routes between the capital city and other locations. I am going to use packages such as rasterio, geopy, geopandas, and matplotlib for this purpose. Let’s get started.

Spatial distribution of population and population density of Nepal

I found the open-access data of the population and population density of Nepal for 2020 on the website of WorldPop. The WorldPop initiative provides the open access archive of spatial demographic datasets for Central and South America, Africa, and Asia to support development, disaster response, and health applications (WorldPop, 2022). First, I downloaded the required files from the URL link directly to my local folder using wget.download(URL, folder). The data was available in the TIF (Tagged Image File) format. TIF is a high-quality graphics format used for storing raster graphics (two-dimensional image as a rectangular matrix or grid of square matrix), which is popular among graphic artists (FileInfo, 2022).

I used the rasterio package to open the TIF files for both population as shown below:

#Population file
file1 = "data/npl_ppp_2020.tif"

#Open file
raster_file1 = rasterio.open(file1)

#Read raster data
pop_data = raster_file1.read(1)

The data in the TIF file is read as numpy array. I created a simple function, which took the natural logarithm of the raster data and then displayed it as an image.

def plot_raster(raster_data, title = "", figsize = (8, 6)):
 """
 Plot population count in log scale (+1)
 """
 plt.figure(figsize = figsize)
 im1 = plt.imshow(np.log1p(raster_data), cmap = "viridis",)
 plt.title(f"{title}", fontdict = {"fontsize": 14})
 plt.axis("off")
 plt.colorbar(im1, fraction = 0.03)

As a result, I get a visually appealing map of population distribution in Nepal using a natural logarithm scale as shown below. With the brightness level of pixels, it is clear that the highest population of Nepal is in the capital city, Kathmandu, which is located in the central region.

Geographically Nepal is divided into three regions: the Himalayas in the North bordering China, Terai (plains) in the South bordering India, and Hills in between. It is clear from the map that the population decreases as the elevation increases from South to North in Nepal. In the map, we can also see three distinct white tree-like structures. These white lines are due to NaN values and represent the three main rivers of Nepal along with their tributaries: Koshi river in the East, Gandaki river in the Center, and Karnali river in the West.

👁 Log-scaled spatial distribution of population in Nepal in 2020. Image by Author.

Log-scaled spatial distribution of population in Nepal in 2020. Image by Author.

The population density map of Nepal paints a slightly different picture as depicted below. While the population of the Hilly region is sparse as compared to that of the Southern Terai region, the population density is still moderate in the Hills. There are some regions in the South with a population density close to zero (shown by dark patches), which represent the national parks of Nepal. The population density is still highest in Kathmandu, followed by places in the Southern region bordering India. On the contrary, the population density of the Himalayan region in the North is very low due to the extreme weather and difficult terrains.

👁 Log-scaled spatial distribution of population density in Nepal in 2020. Image by Author.

Log-scaled spatial distribution of population density in Nepal in 2020. Image by Author.

Visualising domestic aviation in Nepal

Being a mountainous country comes with its own set of pros and cons. While Nepal is blessed with natural beauty and biodiversity, the mountainous terrain makes several places in the Hilly and Himalayan region inaccessible through road transport. Therefore, having domestic airports is quite important for connectivity, freight and tourism especially in remote places in the context of Nepal.

I got the political and administrative boundaries shapefile of Nepal from Open Data Nepal. I read the shapefile using geopandas. And I plotted the interactive map of Nepal using the explore() feature in the latest version of geopandas as shown below.

👁 Exploring the political and administrative boundaries shapefile of Nepal using geopandas. Image by Author.

Exploring the political and administrative boundaries shapefile of Nepal using geopandas. Image by Author.

Geocoding

From the website of the Civil Aviation Authority of Nepal, I got the list of 34 operational airports in Nepal To plot these airports on the map of Nepal, I needed to get the longitude and latitude of each airport. Geocoding is the process of transforming addresses into geographic coordinates (longitude and latitude). The geopy library in Python is useful for both geocoding and reverse geocoding as shown in the code snippet below. Using geopy, I got the coordinates of 31 operational airports in Nepal.

👁 Geocoding and reverse geocoding using geopy library in Python. Image by Author.

Geocoding and reverse geocoding using geopy library in Python. Image by Author.

Plotting operational airports in the map of Nepal

Having nepal as the geopandas dataframe and df containing Airports, Latitude, and Longitude columns, I plotted the airports in the map of Nepal using the following code:

fig, ax = plt.subplots(figsize = (8, 6))
plt.rcParams["font.size"] = 14

nepal.plot(color = "silver", ax = ax)

for lat, lon in zip(df["Latitude"],df["Longitude"]):
 ax.plot(lon, lat, marker = "o", color = "red")

plt.title("Operational airports in Nepal", fontweight = "bold")

👁 Map of Nepal with 31 operational airports. Image by Author.

Map of Nepal with 31 operational airports. Image by Author.

The map above shows the distribution of 31 operational airports across Nepal.

Convex Hull

I was curious to know if I want to cover the largest possible area of Nepal by going on a round trip, what would be the flight route I had to opt for in a hypothetical case. To get the answer to this question, I used the ConvexHull method. In geometry, a convex hull or convex envelope is the smallest convex set that contains it. I created an array of points containing latitude and longitude of each airport of Nepal and passed it inside scipy.spatial.ConvexHull() method.

In the 2-D case, the unique elements of the flattened simplices attribute of the ConvexHull object holds the pairs of indices of the points that make up the line segments of the convex hull. In this way, I was able to get the indices of the airports that would complete the convex hull.

from scipy.spatial import ConvexHull, convex_hull_plot_2d
points = np.array([[lon, lat] for lat, lon in zip(airport_latitudes, airport_longitudes)])
hull = ConvexHull(points)

hull_indices = np.unique(hull.simplices.flat)
hull_points = points[hull_indices, :]

fig, ax = plt.subplots()
nepal.plot(ax = ax, color = "silver")

ax.plot(points[:,0], points[:,1], 'o')
ax.plot(hull_points[:,0], hull_points[:,1], "ro", alpha = 0.25)

for simplex in hull.simplices:
 ax.plot(points[simplex, 0], points[simplex, 1], 'r - ')

plt.title("Convex Hull of the operational airports in Nepal")

👁 Convex hull of the operational airports in Nepal. Image by Author.

Convex hull of the operational airports in Nepal. Image by Author.

Domestic flights from Kathmandu to other locations

As discussed in the first section of this post, Kathmandu, the capital city of Nepal, has both the highest population and population density in Nepal. The only operational international airport of Nepal (at the time of writing) is located in Kathmandu. A separate unit of the same airport also serves as the domestic airport in Kathmandu. There are 12 possible domestic flights from Kathmandu to other airports in Nepal (Flight Connections, 2022). The following dataframe flights_from_ktm shows the name, longitude, and latitude of 12 airports having flights from Tribhuvan International Airport in Kathmandu.

👁 Dataframe showing the list of 12 airports along with their latitude and longitude that has a connection with Tribhuvan International Airport in Kathmandu, Nepal.

Dataframe showing the list of 12 airports along with their latitude and longitude that has a connection with Tribhuvan International Airport in Kathmandu, Nepal.

I used the following code to plot the flight routes from Kathmandu to the other 12 airports in Nepal.

👁 Domestic flight routes with connection to Kathmandu, Nepal. Image by Author.

Domestic flight routes with connection to Kathmandu, Nepal. Image by Author.

The plot above shows the flight connection between the airport in Kathmandu with 12 other airports in Nepal. The airport in Kathmandu is represented by a blue marker, and the others by red.

I came to learn recently that flights do not take a straight line as it normally appears on the flat map. Airplanes follow "great circle" routes to account for the curvature of the earth (Wonderopolis, 2017 and Forbes, 2020). To represent these curvatures in the flight routes, I used a simple function called draw_curve() to plot curved paths between two points instead of using a straight line.

#Function to draw a curved line between two points
def draw_curve(p1, p2):
 """Return curves lines between two points p1 and p2."""
 #np.cosh(x) is equivalent to 0.5*(np.exp(x) - np.exp(-x))
 a = (p2[1] - p1[1])/ (np.cosh(p2[0]) - np.cosh(p1[0]))
 b = p1[1] - a * np.cosh(p1[0])
 x = np.linspace(p1[0], p2[0], 100)
 y = a * np.cosh(x) + b
 return x, y

👁 Domestic flight routes with connection to Kathmandu, Nepal are shown as curved lines. Image by Author.

Domestic flight routes with connection to Kathmandu, Nepal are shown as curved lines. Image by Author.

Consequently, I get curved lines between the Tribhuvan International Airport in Kathmandu and the other 12 airports, that the former has a connection.

Conclusion

In this post, I described how I represented the spatial distribution of population and population density data in the map of Nepal in the first section. In the subsequent section, I explained the geocoding technique using the geopy library to get the coordinates of the domestic airports in Nepal, and the ConvexHull() method of scipy.spatial. Furthermore, I demonstrated a way to plot the distribution of operational airports on the map of Nepal, and the flight routes between Kathmandu and other places of Nepal using the geopandas package.

The data and notebook for this post are available in this GitHub repository. Thank you for reading!