In [1]:
import geopandas as gpd
import pandas as pd
import math
from tqdm import trange, tqdm
import numpy as np
import matplotlib.pyplot as plt
from sklearn.metrics import silhouette_score
from sklearn.cluster import KMeans
import matplotlib.patheffects as pe
import seaborn as sns
from scipy.cluster import hierarchy
import esda
import libpysal
Proposed framework: Leveraging temporal changes of spatial accessibility measurements for better policy implications¶
Supply¶
EV charging stations located in Seoul, South Korea
Source: Korea Environment Corporation EV Monitor. https://www.ev.or.kr/evmonitor.
Demand¶
Floating population over 24 hours
Source: Seoul Open Data Plaza - Floating Population. https://data.seoul.go.kr/dataVisual/seoul/seoulLivingPopulation.do.
Mobility¶
Estimated travel time from demand locations to supply location
Source: iNavi Maps API. https://docs.toast.com/en/Application%20Service/Maps/en/Overview/
In [2]:
# Loading DEMAND data
demand = gpd.read_file('./data/demand.shp')
demand = demand.set_index('GRID_ID')
demand['dummy'] = 'dummy'
study_area = demand.dissolve(by='dummy')
# Plotting the distribution of floating population over 24 hours.
fig = plt.figure(figsize=(18, 10))
colorBrewer = ['#eff3ff', '#bdd7e7', '#6baed6', '#2171b5']
plt_cls = [10, 100, 1000, 10000, 100000]
for qTime in range(24):
ax = fig.add_subplot(4, 6, qTime + 1)
ax.get_xaxis().set_visible(False)
ax.get_yaxis().set_visible(False)
for idx, cls in enumerate(plt_cls):
if idx == 0:
continue
temp_hxg = demand.loc[(plt_cls[idx-1] < demand[f'Flo_WD_{qTime}']) & (demand[f'Flo_WD_{qTime}'] <= cls)]
if temp_hxg.shape[0] > 0:
temp_hxg.plot(ax=ax, color=colorBrewer[idx-1])
ax.text(0.95, 0.95, str(qTime), fontsize=20, ha='right', va='top', transform=ax.transAxes)
fig.text(0.5, 0.9, 'Distribution of floating population over 24 hours', ha='center', size=24)
Out[2]:
In [3]:
# Loading SUPPLY data
supply = gpd.read_file('./data/supply.shp')
supply = supply.set_index('StationID')
# Plot the changes in the degree of supply over 24 hours
fig = plt.figure(figsize=(18, 10))
for qTime in range(24):
ax = fig.add_subplot(4, 6, qTime + 1)
ax.get_xaxis().set_visible(False)
ax.get_yaxis().set_visible(False)
study_area.plot(alpha=0.1, ax=ax)
supply.loc[supply[f'Hour_WD_{qTime}'] == 'Y'].plot(ax=ax, markersize=supply['Weight'], color='black')
ax.text(0.95, 0.95, str(qTime), fontsize=20, ha='right', va='top', transform=ax.transAxes)
ax.text(0.95, 0.12,
'Open Stations: ' + str(supply.loc[supply[f'Hour_WD_{qTime}'] == 'Y'].shape[0]),
fontsize=12, ha='right', va='top', transform=ax.transAxes)
fig.text(0.5, 0.9, 'Locations and weights of operating EV charging stations over 24 hours', ha='center', size=24)
Out[3]:
In [4]:
# Loading MOBILITY data
mobility = pd.read_csv('./data/Mobility.csv')
mobility = mobility.loc[mobility['GRID_ID'].isin(demand.index)]
# Plot the changes in the degree of mobility over 24 hours
fig = plt.figure(figsize=(20, 10))
ax = fig.add_subplot(1, 1, 1)
mobility.loc[:, [f'Tra_WD_{qTime}' for qTime in range(0, 24)]].boxplot(
showfliers=False, fontsize=14, color=dict(boxes='black', whiskers='black', medians='red', caps='black'))
ax.set_xticklabels([f'{qTime}' for qTime in range(24)])
ax.set_xlabel('Hours', fontsize=20)
ax.set_ylabel('Estimated Travel Time (Minutes)', fontsize=20)
Out[4]:
1.2. Preliminary analysis: determining threshold travel time¶
In [5]:
accessible_hxg = pd.DataFrame(index=[f't_{qTime}' for qTime in range(24)], columns=[f'thres_{i}' for i in range(0, 31, 1)])
accessible_hxg = accessible_hxg.fillna(0.0)
for qTime in trange(24):
trvl_time_matrix = pd.DataFrame(index=supply.index, columns=demand.index)
for grid in demand.index:
temp_mobility = mobility.loc[mobility['GRID_ID'] == grid][['Station_ID', f'Tra_WD_{qTime}']]
temp_mobility = temp_mobility.set_index('Station_ID')
trvl_time_matrix[grid] = temp_mobility
cum_pop = []
for threshold_trvl_time in range(0, 31, 1):
accessible_grid = []
for grid in demand.index:
temp_matrix = trvl_time_matrix[grid].loc[trvl_time_matrix[grid] < threshold_trvl_time]
if temp_matrix.shape[0] > 0:
accessible_hxg.at[f't_{qTime}', f'thres_{threshold_trvl_time}'] += 1
accessible_hxg.at[f't_{qTime}', f'thres_{threshold_trvl_time}'] = accessible_hxg.at[f't_{qTime}', f'thres_{threshold_trvl_time}'] / 703
accessible_hxg
Out[5]:
In [6]:
fig, ax = plt.subplots(figsize=(10,8))
accessible_hxg[[f'thres_{i}' for i in range(26)]].transpose().plot(color='black', linewidth=0.5, ax=ax)
for x in range(5, 25, 5):
plt.axvline(x=x, linewidth='1', linestyle='dashed', color='grey')
for y in range(0, 11, 1):
plt.axhline(y=y/10, linewidth='1', linestyle='dashed', color='grey')
plt.axvline(x=15, linewidth='5', linestyle='dashed', color='red')
ax.set_ylim(0, 1.001)
ax.set_xlim(0, 25)
ax.set_yticks([0, 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1.0])
ax.set_yticklabels([f'{i}%' for i in range(0, 110, 10)], fontsize=15)
ax.set_xticks([0,5, 10, 15, 20, 25])
ax.set_xticklabels([f'{i}' for i in range(0, 30, 5)], fontsize=15, rotation=0)
ax.set_ylabel(f'Hexagons accessible \n to EV charging stations', fontsize=25)
ax.set_xlabel('Threshold travel time (Minutes)', fontsize=25)
ax.get_legend().remove()
plt.tight_layout()
plt.show()
In [7]:
# Travel time to the nearest EV charging station from each hexagon (represent minimum travel to obtain service)
nearest_trvl_time = pd.DataFrame(index=mobility['GRID_ID'].unique(), columns=[f'WD_{i}' for i in range(24)])
for grid in tqdm(mobility['GRID_ID'].unique()):
for i in range(24):
nearest_trvl_time.at[grid, f'WD_{i}'] = mobility.loc[(mobility['GRID_ID'] == grid) & (mobility[f'Tra_WD_{i}'] != 0)][f'Tra_WD_{i}'].min()
# Number of accessible EV charging station per each threshold travel time
num_of_accessible_supply_5 = pd.DataFrame(index=mobility['GRID_ID'].unique(), columns=[f'WD_{i}' for i in range(24)])
num_of_accessible_supply_10 = pd.DataFrame(index=mobility['GRID_ID'].unique(), columns=[f'WD_{i}' for i in range(24)])
num_of_accessible_supply_15 = pd.DataFrame(index=mobility['GRID_ID'].unique(), columns=[f'WD_{i}' for i in range(24)])
num_of_accessible_supply_20 = pd.DataFrame(index=mobility['GRID_ID'].unique(), columns=[f'WD_{i}' for i in range(24)])
num_of_accessible_supply_25 = pd.DataFrame(index=mobility['GRID_ID'].unique(), columns=[f'WD_{i}' for i in range(24)])
# Average Percentage of accessible EV charging station per each threshold travel time
percent_accessible_supply = pd.DataFrame(index=[f't_{minutes}' for minutes in [5, 10, 15, 20, 25]], columns=[f'WD_{i}' for i in range(24)])
for i in trange(24):
for grid in mobility['GRID_ID'].unique():
num_of_accessible_supply_5.at[grid, f'WD_{i}'] = mobility.loc[(mobility['GRID_ID'] == grid) &
(mobility[f'Tra_WD_{i}'] <= 5)].shape[0]
num_of_accessible_supply_10.at[grid, f'WD_{i}'] = mobility.loc[(mobility['GRID_ID'] == grid) &
(mobility[f'Tra_WD_{i}'] <= 10)].shape[0]
num_of_accessible_supply_15.at[grid, f'WD_{i}'] = mobility.loc[(mobility['GRID_ID'] == grid) &
(mobility[f'Tra_WD_{i}'] <= 15)].shape[0]
num_of_accessible_supply_20.at[grid, f'WD_{i}'] = mobility.loc[(mobility['GRID_ID'] == grid) &
(mobility[f'Tra_WD_{i}'] <= 20)].shape[0]
num_of_accessible_supply_25.at[grid, f'WD_{i}'] = mobility.loc[(mobility['GRID_ID'] == grid) &
(mobility[f'Tra_WD_{i}'] <= 25)].shape[0]
percent_accessible_supply.at['t_5', f'WD_{i}'] = num_of_accessible_supply_5[f'WD_{i}'].mean() / supply.shape[0]
percent_accessible_supply.at['t_10', f'WD_{i}'] = num_of_accessible_supply_10[f'WD_{i}'].mean() / supply.shape[0]
percent_accessible_supply.at['t_15', f'WD_{i}'] = num_of_accessible_supply_15[f'WD_{i}'].mean() / supply.shape[0]
percent_accessible_supply.at['t_20', f'WD_{i}'] = num_of_accessible_supply_20[f'WD_{i}'].mean() / supply.shape[0]
percent_accessible_supply.at['t_25', f'WD_{i}'] = num_of_accessible_supply_25[f'WD_{i}'].mean() / supply.shape[0]
percent_accessible_supply
Out[7]:
In [8]:
fig, ax = plt.subplots(figsize=(12, 8))
line_style = ['--', '--', '-', '--', '--']
line_colors = ['black', 'black', 'red', 'black', 'black']
for idx, minutes in enumerate([5, 10, 15, 20, 25]):
ax.plot(range(0, 24), percent_accessible_supply.loc[f't_{minutes}'].values, line_style[idx], linewidth=1.5, color=line_colors[idx], label=f'{minutes} Minutes')
ax.plot(range(0, 24), percent_accessible_supply.loc[f't_15'].values, '-', linewidth=5, color='red')
for y in range(0, 10, 1):
plt.axhline(y=y/10, linewidth='1', linestyle='dashed', color='grey')
for x in range(1, 23, 1):
plt.axvline(x=x, linewidth='1', linestyle='dashed', color='grey')
plt.legend(loc='upper center', fontsize='xx-large')
ax.set_xlim(-0.1, 23.1)
ax.set_xticks(list(range(0, 24)))
ax.set_xticklabels([f'{qTime}' for qTime in range(24)], fontsize=15)
ax.set_xlabel('Hours', fontsize=25)
ax.set_ylim(0, 1.001)
ax.set_yticks([0, 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1.0])
ax.set_yticklabels([f'{i}%' for i in range(0, 110, 10)], fontsize=15)
ax.set_ylabel(f'Average EV charging stations \n accessible from hexagons', fontsize=25)
ax.get_legend().remove()
plt.tight_layout()
plt.show()