There are many ways to visualize geospatial information using Python, but this time I will use ** Basemap Matplotlib Toolkit ** to visualize location information.

Click here for articles introduced in the past regarding geospatial information visualization, for example.

Note on how to use folium (visualization of location information)

Visualize and understand the regional mesh of Japan on a map
Introduction of "scikit-mobility", a library that allows you to easily analyze human flow data with Python (1)

1. About Basemap

Simply put, Basemap allows you to add matplotlib plotting capabilities while drawing various map projections, map tiles, coastlines, rivers, and more. It is mainly used by geoscientists.

The matplotlib basemap toolkit is a library for plotting 2D data on maps in Python. Basemap does not do any plotting on it’s own, but provides the facilities to transform coordinates to one of 25 different map projections (using the PROJ.4 C library). Matplotlib is then used to plot contours, images, vectors, lines or points in the transformed coordinates. Shoreline, river and political boundary datasets (from Generic Mapping Tools) are provided, along with methods for plotting them. The GEOS library is used internally to clip the coastline and polticial boundary features to the desired map projection region. Basemap is geared toward the needs of earth scientists, particularly oceanographers and meteorologists. (Quoted from https://matplotlib.org/basemap/users/intro.html)

1-1. Basics of Basemap

First, use conda to install baseline.

$ conda install -c anaconda basemap
$ conda install -c conda-forge basemap-data-hires

from mpl_toolkits.basemap import Basemap
import matplotlib.pyplot as plt

fig = plt.figure()
#Specify the range to draw the map.
m = Basemap(llcrnrlat=30, urcrnrlat=50, llcrnrlon=125, urcrnrlon=150) #Create Basemap instance
#Subtract latitude / longitude every 10 degrees. The second option sets where to put the label up, down, left and right.
m.drawparallels(np.arange(-90, 90, 10), labels=[ True,False, True, False]) #Latitude line
m.drawmeridians(np.arange(0, 360, 10), labels=[True, False, False, True]) #Longitude line
m.drawcoastlines() #Border line
plt.show()

1-2. Detailed map settings

You can flexibly change the map resolution and projection method. (There are many detailed settings ...)

class mpl_toolkits.basemap.Basemap(llcrnrlon=None, llcrnrlat=None, urcrnrlon=None, urcrnrlat=None, llcrnrx=None, llcrnry=None, urcrnrx=None, urcrnry=None, width=None, height=None, projection='cyl', resolution='c', area_thresh=None, rsphere=6370997.0, ellps=None, lat_ts=None, lat_1=None, lat_2=None, lat_0=None, lon_0=None, lon_1=None, lon_2=None, o_lon_p=None, o_lat_p=None, k_0=None, no_rot=False, suppress_ticks=True, satellite_height=35786000, boundinglat=None, fix_aspect=True, anchor='C', celestial=False, round=False, epsg=None, ax=None)

projection\resolution	Low resolution	Medium resolution	High resolution
Equirectangular projection
Mercator projection
Lambert projection

By the way, the default value when nothing is set is projection ='cyl' (equirectangular projection) resolution ='c' (coarse resolution)

2. Location information visualization method

2-1. Data source

"Pseudo-human flow data" based on SNS analysis data by the Center for Spatial Information Science, University of Tokyo

This time, I would like to use the pseudo-human flow data that is released for free as a sample.

2-2. Data structure

Read the Kansai area pseudo-human flow data on September 16, 2013 and check the contents.

import pandas as pd
df = pd.read_csv('./Kansai/2013_09_16.csv')
df

'''
First column: User ID
Second column: Gender estimate (1):Male, 2:Female, 0.3: Unknown, NA: Unestimated)
3rd column: Date / time (24 hours every 5 minutes)
Fourth column: Latitude
5th column: Longitude
6th column: Resident category(Major classification)* Character string type
7th column: Resident category(Subcategory)* Character string type
8th column: State(Stay or move)* Character string type
9th column: Resident category ID(Corresponds to the 6th and 7th lines)
'''

This time, we will target the moving user (STAY_MOVE =='MOVE'). Also, while dividing the timestamp into hours and minutes, give a rank for each user and each hour and format it into a form that is easy to handle.

from dfply import *

df = df >> filter_by(X.STAY_MOVE=='MOVE') >> select(columns_to(X.lon, inclusive=True))
df = df >> separate(X.timestamp, ['col1','hour','col2','minute','col3'], sep=[10,13,14,16],convert=True) >> select(~X.col1, ~X.col2, ~X.col3)
df = df >> group_by(X.uid,X.hour) >> mutate(rk=row_number(X.minute))
df

2-3. Basemap initial value setting and location information visualization

Now, I'm going to visualize the location information, but first I will set the map that I will consider as the base this time. There are also methods that can draw the background map of ArcGIS and the borders of prefectures and cities, towns and villages, so I will apply them as well.

For prefectural borders and municipal borders, please download the Shapefile from __ here __. The file structure is as follows.

-gadm
 -gadm36_JPN_1.cpg 
 -gadm36_JPN_1.shp
 -gadm36_JPN_2.dbf
 -gadm36_JPN_2.shx
 -gadm36_JPN_1.dbf
 -gadm36_JPN_1.shx
 -gadm36_JPN_2.prj
 -gadm36_JPN_1.prj
 -gadm36_JPN_2.cpg
 -gadm36_JPN_2.shp

** Initial map setting **

def basemap(): 
    fig = plt.figure(dpi=300)
    m = Basemap(projection="cyl", resolution="i", llcrnrlat=33.5,urcrnrlat=36, llcrnrlon=134, urcrnrlon=137)
    m.drawparallels(np.arange(-90, 90, 0.5), labels=[True, False, True, False],linewidth=0.0, fontsize=8)
    m.drawmeridians(np.arange(0, 360, 0.5), labels=[True, False, False, True],linewidth=0.0, rotation=45, fontsize=8)
    m.drawcoastlines(color='k')
    m.readshapefile('./gadm/gadm36_JPN_1', 'prefectural_bound1', color='k', linewidth=.8) #Prefectural border
    m.readshapefile('./gadm/gadm36_JPN_2', 'prefectural_bound2', color='lightgray', linewidth=.5) #Municipal border
    m.arcgisimage(service='World_Street_Map', verbose=True, xpixels=1000, dpi=300)

** If you do so far, you can plot like matplotlib. ** **

scatter(x, y, *args, **kwargs)

https://basemaptutorial.readthedocs.io/en/latest/plotting_data.html#scatter

tmp1=df[(df['gender']==1) & (df['hour']==9) & (df['rk']==1)]
tmp2=df[(df['gender']==2) & (df['hour']==9) & (df['rk']==1)]

basemap()
plt.scatter(tmp1['lon'],tmp1['lat'],color='b',s=0.5) #Men in blue
plt.scatter(tmp2['lon'],tmp2['lat'],color='r',s=0.5) #Woman in red

Here, only the first log at 9 o'clock is plotted for each user.

hexbin(x, y, **kwargs)

https://basemaptutorial.readthedocs.io/en/latest/plotting_data.html#hexbin

tmp=df[(df['hour']==9) & (df['rk']==1)]

basemap()
plt.hexbin(tmp['lon'],tmp['lat'], gridsize=50, cmap='rainbow', mincnt=1, bins='log',linewidths=0.3, edgecolors='k')

You can also change the size of the hexagrid with gridsize.

3. Bonus

Visualize the frequency distribution by superimposing the precipitation distribution

For example, it is possible to multiply weather data that is likely to affect the movement of people. By the way, September 16, 2013 was the day when the typhoon passed.


From the Japan Meteorological Agency	1km mesh radar echo

Visualize the moving speed and direction of a person

Vector display is also possible by calculating the speed and azimuth from the position information of the start point and end point. (Aside from whether there is a use)

That's it! Until the end Thank you for reading!

[PYTHON] Visualize location information with Basemap