Import Packages¶
In [1]:
import pandas as pd
import geopandas as gpd
import matplotlib.pyplot as plt
import urllib
import zipfile
import glob
import os
from shapely.geometry import Point
from geopandas import GeoDataFrame
import random as rd
%matplotlib inline
Import and Clean Data¶
For this mini project, we will need 2 layers: Public transit stations and COVID-19 vaccination clinic locations
In [2]:
try:
os.mkdir('downloads')
except:
print('folder exists')
In [3]:
def download_unzip_shp(link, foldername):
# download shapefile as a zipped folder
urllib.request.urlretrieve(link, foldername + ".zip")
# unzip to another folder
with zipfile.ZipFile(foldername + ".zip","r") as zip_ref:
zip_ref.extractall(foldername)
shp_name = glob.glob(foldername + '/*.shp')
return shp_name
In [4]:
# bus stations:
## check this link: https://data.cityofnewyork.us/Transportation/Bus-Stop-Shelters/qafz-7myz
buslink = "https://data.cityofnewyork.us/api/geospatial/qafz-7myz?method=export&format=Shapefile"
busname = "downloads/bus_shp"
bus_shp = download_unzip_shp(buslink, busname)
In [5]:
bus = gpd.read_file(bus_shp[0])
# convert from GCS to PCS for further actions
bus = bus.to_crs(epsg=6434)
bus.plot()
plt.show()
In [6]:
# uncomment to check the shapefile coordinate system
#bus.crs
In [7]:
# metro stations:
## check this link: https://data.cityofnewyork.us/Transportation/Subway-Stations/arq3-7z49
metrolink = "https://data.cityofnewyork.us/api/geospatial/arq3-7z49?method=export&format=Shapefile"
metroname = "downloads/metro_shp"
metro_shp = download_unzip_shp(metrolink, metroname)
In [8]:
metro = gpd.read_file(metro_shp[0])
# convert from GCS to PCS for further actions
metro = metro.to_crs(epsg=6434)
metro.plot()
plt.show()
In [9]:
# join them together:
bshort = bus[['shelter_id','geometry']]
mshort = metro[['name','geometry']]
transit = gpd.GeoDataFrame( pd.concat([bshort, mshort], ignore_index=True) )
transit.plot()
plt.show()
In [10]:
transit['ID'] = transit.index
In [11]:
# uncomment to check the shapefile coordinate system
#metro.crs
In [12]:
# if you would like to download a .csv version: try codes below
# # bus stations:
# urllib.request.urlretrieve("https://data.cityofnewyork.us/api/geospatial/qafz-7myz?method=export&format=Shapefile", "bus.csv")
# # metro stations:
# urllib.request.urlretrieve("https://data.cityofnewyork.us/api/views/kk4q-3rt2/rows.csv?accessType=DOWNLOAD", "metro.csv")
In [13]:
# census tracts:
## check this link: https://data.cityofnewyork.us/City-Government/2010-Census-Tracts/fxpq-c8ku
ctlink = "https://data.cityofnewyork.us/api/geospatial/fxpq-c8ku?method=export&format=Shapefile"
ctname = "downloads/2010 Census Tracts_shp"
ct_shp = download_unzip_shp(ctlink, ctname)
In [14]:
ct = gpd.read_file(ct_shp[0])
# convert from GCS to PCS for further actions
ct = ct.to_crs(epsg=6434)
ct.plot()
plt.show()
In [15]:
# uncomment to check the shapefile coordinate system
#ct.crs
In [16]:
# vaccine station
# the processing notebook is in the floder
# several geocoding errors are fixed manually
vs = pd.read_csv('Vac_station_final.csv')
vs.head(2)
Out[16]:
In [17]:
## generate geo-spatial dataframe using csv files
geometry = [Point(xy) for xy in zip(vs.longitude, vs.latitude)]
vs = GeoDataFrame(vs, geometry=geometry)
vs.set_crs(epsg=4326, inplace=True)
vs = vs.to_crs(epsg=6434)
#vs.head(1)
In [18]:
vs.plot()
plt.show()
In [ ]:
Analysis¶
In [19]:
plt.rcParams["figure.figsize"] = (10,10) # change the plots to a larger size
First, select a random station¶
In [20]:
# select 1 random station:
station_code = rd.randrange(3896)
station = transit[transit['ID'] == station_code]
# plot both the select station and the census tracts
ax = ct.plot()
station.plot(ax=ax, color='brown')
plt.show()
Then, generate 1 mile radius buffer and clip census tract to 3 mile bounding box¶
In [21]:
radius = 1 * 5280 ## convert 1 mile radius to 5280 ft
# generate buffer around the selected station
buffer = station.geometry.buffer(radius)
bounding = station.buffer(radius * 3, cap_style=3) ## 3 is square bounding
ct_clip = gpd.clip(ct,bounding)
ax = ct_clip.plot()
buffer.plot(ax=ax, color='yellow')
station.plot(ax=ax, color='brown')
plt.show()
Next, clip vaccination station layer with bounding box, and futher select vaccination stations in 1-mile buffer¶
In [22]:
vac_clip = gpd.clip(vs,buffer)
other_vac_clip = gpd.clip(vs,bounding)
ax = ct_clip.plot()
buffer.plot(ax=ax, color='yellow')
station.plot(ax=ax, color='brown')
other_vac_clip.plot(ax=ax, marker = "o", markersize=20, color = 'violet')
vac_clip.plot(ax=ax, marker = "*", markersize=260, color = 'hotpink')
plt.show()
In [23]:
vac_clip.head()
Out[23]:
In [ ]:
In [24]:
# if you would like to generate shapefile/csv output, uncomment lines below
# os.mkdir('outputs')
# vac_clip.to_file(driver='ESRI Shapefile', filename = "outputs/vac_clip.shp")
# vac_clip.to_csv('outputs/vac_clip.csv')
In [ ]: