Sign in for more options. Each user starts with 1,000 free credits.
KzO9K
In [1]:
url = 'https://data.cityofnewyork.us/resource/5rq2-4hqu.csv'
In [ ]:
In [2]:
import pandas as pd
In [3]:
pd.read_csv(url + '?$select=count(*)')
Out[3]:
In [4]:
t = pd.read_csv(url + '?$limit=1')
t.iloc[0]
Out[4]:
In [5]:
# Import Tree data set
In [6]:
# Also import air Pol dataset
In [7]:
import pandas as pd
def load(
endpoint_url,
selected_columns=None,
buffer_size=1000,
search_term_by_column=None,
**kw,
):
buffer_url = (f'{endpoint_url}?$limit={buffer_size}')
if selected_columns:
select_string = ','.join(selected_columns)
buffer_url += f'&$select={select_string}'
for column, search_term in (search_term_by_column or {}).items():
buffer_url += f'&$where={column}+like+"%25{search_term}%25"'
print(buffer_url)
tables = []
if endpoint_url.endswith('.json'):
f = pd.read_json
else:
f = pd.read_csv
t = f(buffer_url, **kw)
while len(t):
print(len(tables) * buffer_size + len(t))
tables.append(t)
offset = buffer_size * len(tables)
t = f(buffer_url + f'&$offset={offset}', **kw)
return pd.concat(tables, ignore_index=True, sort=False)
In [8]:
# endpoint_url = 'https://data.cityofnewyork.us/resource/nwxe-4ae8.csv'
endpoint_url = url
selected_columns = 'tree_id', 'latitude', 'longitude'
buffer_size = 100000
t = load(endpoint_url, selected_columns, buffer_size)
In [9]:
tree_table = t
In [10]:
tree_table[:3]
Out[10]:
In [11]:
import gzip
from invisibleroads_macros.disk import uncompress
from os.path import exists
from urllib.request import urlretrieve
def download(target_path, source_url):
if not exists(target_path):
urlretrieve(source_url, target_path)
return target_path
def download_zip(target_folder, source_url):
archive_path = download(target_folder + '.zip', source_url)
return uncompress(archive_path, target_folder)
def download_gz(target_path, source_url):
archive_path = download(target_path + '.gz', source_url)
with gzip.open(archive_path, 'rb') as f:
open(target_path, 'wb').write(f.read())
return target_path
In [12]:
# Load air pollution raster
source_url = (
'https://data.cityofnewyork.us/download/'
'q68s-8qxv/application%2Fzip')
source_folder = '/tmp/pollution'
download_zip(source_folder, source_url)
Out[12]:
In [ ]:
In [13]:
# Load air pollution raster
source_url = (
'https://data.cityofnewyork.us/download/'
'q68s-8qxv/application%2Fzip')
source_folder = '/tmp/pollution'
download_zip(source_folder, source_url)
Out[13]:
In [14]:
# Show nitrogen dioxide images
from glob import glob
sorted(glob('/tmp/pollution/*no2*.ovr'))
Out[14]:
In [15]:
# Show nitrogen dioxide images in aa2
from glob import glob
sorted(glob('/tmp/pollution/aa2_no2300m/*'))
Out[15]:
In [16]:
# Install packages
# If this keeps failing, try starting the notebook with more memory (e.g. MEDIUM or LARGE)
import subprocess
assert subprocess.call('pip install rasterio'.split()) == 0
In [17]:
# Load image
import rasterio
pollution_raster = rasterio.open('/tmp/pollution/aa2_no2300m/prj.adf')
In [18]:
# Get coordinate reference (defines how each pixel corresponds to a specific location)
pollution_proj4 = pollution_raster.crs.to_string()
pollution_proj4
Out[18]:
In [19]:
# Get coordinate reference (defines how each pixel corresponds to a specific location)
pollution_proj4 = pollution_raster.crs.to_proj4()
pollution_proj4
Out[19]:
In [20]:
# Load image array
x = pollution_raster.read()
In [21]:
x.shape
Out[21]:
In [22]:
pollution_raster.height, pollution_raster.width
Out[22]:
In [23]:
# Get coordinates for specific pixel
pollution_raster.xy(pollution_raster.height // 2, pollution_raster.width // 2)
Out[23]:
In [24]:
x.max()
Out[24]:
In [25]:
x.min()
Out[25]:
In [26]:
x.shape
Out[26]:
In [27]:
# Enable inline plots
%matplotlib inline
In [28]:
array = x[0]
array.shape
Out[28]:
In [29]:
# from matplotlib.colors import Normalize
# normalize = Normalize(vmin=x.min(), vmax=x.max())
In [30]:
import matplotlib.pyplot as plt
plt.imshow(array, vmin=0);
In [31]:
import geotable
In [32]:
# Convert neighborhood geometries into pollution coordinate reference
source_proj4 = neighborhood_proj4 = geotable.LONGITUDE_LATITUDE_PROJ4
target_proj4 = pollution_proj4
from geotable.projections import get_transform_shapely_geometry
f = get_transform_shapely_geometry(source_proj4, target_proj4)
In [33]:
# get zipcode gem
# for each zipcode gem, get the tree count
# sum the air poll.
#that will be our data set table
In [34]:
# get zipcode gem
In [35]:
url = 'https://data.cityofnewyork.us/download/i8iw-xf4u/application%2Fzip'
In [36]:
# # Get UTM projection
# import geotable
# utm_proj4 = geotable.load_utm_proj4(url)
In [37]:
# # Get zipcode area in square meters
# nyc_zipcode_table = geotable.load(url, target_proj4=utm_proj4)
# nyc_zipcode_table.iloc[0]['geometry_object'].area
In [38]:
from geotable.projections import LONGITUDE_LATITUDE_PROJ4
In [ ]:
In [39]:
nyc_zipcode_table = geotable.load(url, target_proj4=LONGITUDE_LATITUDE_PROJ4)
In [40]:
nyc_zipcode_table[:3]
Out[40]:
In [60]:
utm_proj4 = geotable.load_utm_proj4(url)
utm_geotable = geotable.load(url, target_proj4=utm_proj4)
utm_geotable['Area in Sq Km'] = utm_geotable['geometry_object'].apply(lambda g: g.area / 1000000)
nyc_zipcode_table['Area in Sq Km'] = utm_geotable['Area in Sq Km']
nyc_zipcode_table
Out[60]:
In [61]:
from shapely.geometry import MultiPoint, Point
In [ ]:
In [62]:
tree_sample_count = 100000
In [ ]:
In [63]:
if tree_sample_count < 0:
tree_sample_count = len(tree_table)
tree_sample_count
Out[63]:
In [64]:
import numpy as np
import random
sample_tree_labels = random.sample(list(tree_table.index), tree_sample_count)
sample_tree_labels = np.array(sample_tree_labels)
sample_tree_labels[:3]
Out[64]:
In [65]:
sample_tree_table = tree_table.loc[sample_tree_labels]
sample_tree_table[:3]
Out[65]:
In [66]:
sample_tree_xys = sample_tree_table[['longitude', 'latitude']].values
sample_tree_xys[:3]
Out[66]:
In [67]:
sample_tree_points = [Point(_) for _ in sample_tree_xys]
In [68]:
sample_tree_collection = MultiPoint(sample_tree_points)
In [69]:
def get_tree_count(r):
g = r['geometry_object']
if g.area:
selected_tree_collection = g.intersection(sample_tree_collection)
if selected_tree_collection.is_empty:
return 0
try:
return len(selected_tree_collection)
except TypeError:
return 1
relative_distances, relative_indices = query_kd_tree(g)
return len(relative_indices)
In [70]:
from geopy.distance import distance as get_distance
def get_geodesic_distance_in_meters(g1, g2):
latitude_longitude1 = g1.centroid.coords[:2][::-1]
latitude_longitude2 = g2.centroid.coords[:2][::-1]
return get_distance(
latitude_longitude1,
latitude_longitude2).meters
def get_sum_distance(r):
g = r['geometry_object']
if g.area:
selected_tree_collection = g.intersection(sample_tree_collection)
if selected_tree_collection.is_empty:
return 0
g_centroid = g.centroid
return sum(get_geodesic_distance_in_meters(
g_centroid, _) for _ in selected_tree_collection)
relative_distances, relative_indices = query_kd_tree(g)
return sum(relative_distances)
In [71]:
whos
In [72]:
nyc_zipcode_table[:3]
Out[72]:
In [73]:
source_sample_count = len(nyc_zipcode_table)
target_geotable = nyc_zipcode_table[:source_sample_count].copy()
target_geotable[
'Tree Count'
] = target_geotable.apply(get_tree_count, axis=1)
# target_geotable[
# 'Sum of Tree Distances'
# ] = target_geotable.apply(get_sum_distance, axis=1)
# target_geotable[
# 'Average Tree Diameter'
# ] = target_geotable.apply(get_average_diameter, axis=1)
target_geotable[:5]
Out[73]:
In [74]:
# Convert neighborhood geometries into pollution coordinate reference
from geotable.projections import get_transform_shapely_geometry
source_proj4 = geotable.LONGITUDE_LATITUDE_PROJ4
target_proj4 = pollution_proj4
f = get_transform_shapely_geometry(source_proj4, target_proj4)
# For each neighborhood polygon, get sum of nitrogen dioxide particulate matter
from rasterio.mask import mask
In [75]:
def get_pollutant_sum(old_polygon):
new_polygon = f(old_polygon)
output_image, output_transform = mask(pollution_raster, [new_polygon], crop=True)
output_image[output_image < 0] = 0 # Zero null values so that we do not get -inf
return output_image.sum()
# Drop rows that do not have a geometry
target_geotable.dropna(subset=['geometry_object'], inplace=True)
# Define new column
target_geotable['pollutant_sum'] = target_geotable[
'geometry_object'].apply(get_pollutant_sum)
In [81]:
target_geotable[:3]
Out[81]:
In [82]:
target_geotable['Tree Count Per Sq km'] = target_geotable['Tree Count'] / target_geotable['Area in Sq Km']
In [83]:
# target_path = target_folder + '/augmented-dataset.csv'
target_geotable.to_csv("zipCodeAndTreeCountDataset.csv", index=False)
# print(f'augmented_table_path = {target_path}')
In [84]:
# train the model
In [ ]:
In [ ]:
In [61]:
# for each zipcode gem, get the tree count
# sum the air poll.
#that will be our data set table
In [ ]: