In [ ]:
import os

import pandas as pd
import geopandas as gpd
import numpy as np

from shapely.geometry import Point

import pysal as ps
from pysal.esda.getisord import G
from pysal.esda.getisord import G_Local

from matplotlib import pyplot as plt
from pylab import figure, scatter, show
import seaborn as sns
In [ ]:
shpPath = "data/shapefiles/pivoted_averaged_turnstile_data.shp"

stationData = gpd.read_file(shpPath)
nyc_streets = gpd.read_file("visualization/basemap/nyc_streets.shp")

listOfTimes = list(stationData.columns)[1:-1]

ESDA

In [ ]:
output_path = 'visualization/esda'

for i in range(len(listOfTimes)):

    vmin, vmax = -.1, .1

    base = nyc_streets.plot(color='black', linewidth=0.25, alpha=.1)
    fig = stationData.plot(ax=base, column=listOfTimes[i], legend=True, figsize=(10,10), markersize=1, vmin=vmin, vmax=vmax, norm=plt.Normalize(vmin=vmin, vmax=vmax))
    
    fig.axis('off')
    fig.figsize=(10,10)
    fig.set_title('Station Intensity', fontdict={'fontsize': '25', 'fontweight' : '3'})
    fig.annotate(listOfTimes[i], xy=(0.1, .15), xycoords='figure fraction', horizontalalignment='left', verticalalignment='top', fontsize=15)
    
    filepath = os.path.join(output_path, listOfTimes[i] +'_intensity_timemap.png')
    chart = fig.get_figure()
    chart.savefig(filepath, dpi=75)

Weights & Spatial Lags

In [4]:
centroids = np.array(list(zip(stationData.geometry.apply(lambda p: p.x), stationData.geometry.apply(lambda p: p.y))))
threshold = ps.min_threshold_dist_from_shapefile(shpPath)
gravityModel = ps.weights.DistanceBand.from_shapefile(shpPath, threshold, alpha = -2.0, binary = False)

fig = figure(figsize=(15,15))
plt.axis('off')
plt.plot(centroids[:,0], centroids[:,1],'.')
for k,neighbors in gravityModel.neighbors.items():
    origin = centroids[k]
    for neighbor in neighbors:
        segment = centroids[[k,neighbor]]
        plt.plot(segment[:,0], segment[:,1], '-')
plt.title('Minimum Distance Threshold Neighbor Graph')
plt.savefig('visualization/weights/station_distancethreshold_neighbors.png');
plt.savefig('visualization/weights/station_distancethreshold_neighbors.svg', format='svg', dpi=1200)
plt.show()
In [ ]:
threshold = ps.min_threshold_dist_from_shapefile(shpPath)
gravityModel = ps.weights.DistanceBand.from_shapefile(shpPath, threshold, alpha = -2.0, binary = False, idVariable='STATION')

weightsFile = ps.open('data/weights/gravity_weights.gal', 'w')
weightsFile.write(gravityModel)
weightsFile.close;
In [6]:
focus = stationData[stationData.STATION == 'Jay St - MetroTech']
focus = focus.to_crs(epsg=4326)
buffer = focus.buffer(threshold)

neighbors = pd.DataFrame()
for neighbor in gravityModel['Jay St - MetroTech']:
    neighbors = neighbors.append(stationData[stationData.STATION == neighbor])
neighbors = gpd.GeoDataFrame(neighbors.geometry)

fig, ax = plt.subplots(figsize=(20, 20))
plt.axis('off')

ax.set_aspect('equal')
ax.set_xlim(-74.02,-73.95)
ax.set_ylim(40.66,40.72)

nyc_streets.plot(ax=ax, color='black', linewidth=0.25)
focus.plot(ax=ax, marker='o', color='red', markersize=10)
neighbors.plot(ax=ax, marker='o', color='blue', markersize=7.5)
buffer.plot(ax=ax, edgecolor='yellow', alpha=.25)

plt.savefig('visualization/weights/jayst_neighbors.png');
plt.show();