<jupyter_start><jupyter_text>Traffic Ticket Quotas in Seattle Agustin Lew 12/10/2019 HCDE512  Import required libraries<jupyter_code>from IPython.core.interactiveshell import InteractiveShell
InteractiveShell.ast_node_interactivity = "all"
%matplotlib inline
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import matplotlib
import geopy
import geopy.distance
import gmplot<jupyter_output><empty_output><jupyter_text>Load data from the raw_data folder<jupyter_code>%pwd
%cd ..
home_directory = %pwd
%cd src
df_seattle = pd.read_csv(home_directory + '/raw_data/wa_seattle_2019_08_13.csv')
df_seattle.head()<jupyter_output>/home/ablew/anaconda3/lib/python3.7/site-packages/IPython/core/interactiveshell.py:3049: DtypeWarning: Columns (10,11,20,21,22,23,26) have mixed types. Specify dtype option on import or set low_memory=False.
  interactivity=interactivity, compiler=compiler, result=result)<jupyter_text>Information on raw data<jupyter_code>print("The raw dataframe has " + str(len(df_seattle)) + " rows and " + str(len(df_seattle.columns)) + " columns.")
print("The dataset has " + str(df_seattle.officer_id_hash.nunique()) + " unique officers.")<jupyter_output>The dataset has 1652 unique officers.<jupyter_text>Column names and types<jupyter_code>df_seattle.dtypes<jupyter_output><empty_output><jupyter_text>Data cleaning* Removed missing values for the columns 'officer_id_hash', 'date', 'time', 'lat'* Changed date and time columns to pandas datetime type* Included only records newer than 2013-01-01* Included only vehicular traffic stops* Selected columns: date, time, lat, lng, officer_id_hash, type, citation_issued, warning_issued, outcome, raw_type_description<jupyter_code>df_seattle = df_seattle[~df_seattle.officer_id_hash.isna()]
df_seattle = df_seattle[~df_seattle.time.isna()]
df_seattle = df_seattle[~df_seattle.date.isna()]
df_seattle = df_seattle[~df_seattle.lat.isna()]
df_seattle['date'] = pd.to_datetime(df_seattle['date'], errors='coerce')
df_seattle['time'] = pd.to_datetime(df_seattle['time'], format='%H:%M:%S', errors='coerce')
df_seattle = df_seattle[df_seattle.date >= '2013-01-01']
df_seattle = df_seattle[['date', 'time', 'lat', 'lng', 'officer_id_hash', 'type', 'citation_issued', 'warning_issued', 'outcome', 'raw_type_description']]
df_seattle = df_seattle[df_seattle.type == 'vehicular']
df_seattle.head()<jupyter_output><empty_output><jupyter_text>Exports clean dataframe to the src folder<jupyter_code>df_seattle.to_csv(home_directory + '/src/clean_seattle.csv', index=False)<jupyter_output><empty_output><jupyter_text>Precinct coordinates (obtained from Google Maps)* N = North Precinct* S = South Precinct* SW = Southwest Precinct* E = East Precinct* W = West Precinct<jupyter_code>sw_precinct = (47.536655, -122.361851)
s_precinct = (47.538801, -122.293376)
w_precinct = (47.616313, -122.336638)
e_precinct = (47.615069, -122.317191)
n_precinct = (47.703082, -122.334784)<jupyter_output><empty_output><jupyter_text>Creates officers dataframe, one row for each officer, showing the average lat/long coordinates, citation_ratio, warning_ratio, and the distance from their average location to each precinct.  Each officer was assigned the closest precinct as their precinct.<jupyter_code># Creates officers dataframe and creates additional columns
df_officers = df_seattle.groupby(df_seattle.officer_id_hash).sum().citation_issued.reset_index()
df_officers['lat'] = df_seattle.groupby(df_seattle.officer_id_hash).mean().reset_index().lat
df_officers['lng'] = df_seattle.groupby(df_seattle.officer_id_hash).mean().reset_index().lng
df_officers['citation_ratio'] = df_seattle.groupby(df_seattle.officer_id_hash).mean().reset_index().citation_issued
df_officers['warning_ratio'] = df_seattle.groupby(df_seattle.officer_id_hash).mean().reset_index().warning_issued

# Calculates distance from the average lat/long values to each precinct
df_officers['S'] = [geopy.distance.geodesic((df_officers.lat.iloc[i], df_officers.lng.iloc[i]), s_precinct).km for i in range(len(df_officers))]
df_officers['SW'] = [geopy.distance.geodesic((df_officers.lat.iloc[i], df_officers.lng.iloc[i]), sw_precinct).km for i in range(len(df_officers))]
df_officers['W'] = [geopy.distance.geodesic((df_officers.lat.iloc[i], df_officers.lng.iloc[i]), w_precinct).km for i in range(len(df_officers))]
df_officers['E'] = [geopy.distance.geodesic((df_officers.lat.iloc[i], df_officers.lng.iloc[i]), e_precinct).km for i in range(len(df_officers))]
df_officers['N'] = [geopy.distance.geodesic((df_officers.lat.iloc[i], df_officers.lng.iloc[i]), n_precinct).km for i in range(len(df_officers))]

precinct_names=['S', 'SW', 'W', 'E', 'N']
precinct_list = []
for i in range(len(df_officers)):
    S = df_officers.S.values[i]
    SW = df_officers.SW.values[i]
    W = df_officers.W.values[i]
    E = df_officers.E.values[i]
    N = df_officers.N.values[i]
    
    precinct_ix = np.argmin([S,SW,W,E,N])
    precinct_list.append(precinct_names[precinct_ix])
    
df_officers['precinct'] = precinct_list

df_officers.head()
str(round(len(df_officers[df_officers.citation_issued == 0]) / df_officers.officer_id_hash.nunique() * 100, 2)) + "% of officers have 0 citations issued."<jupyter_output><empty_output><jupyter_text>Exports officer dataframe to the src folder<jupyter_code>df_officers.to_csv(home_directory + '/src/officers.csv', index=False)<jupyter_output><empty_output><jupyter_text>Creating time series of different granularities: hourly, daily, weekly, monthly for each officer<jupyter_code>daily_timeseries = df_seattle.groupby([df_seattle.officer_id_hash, df_seattle.date]).count().reset_index()[['officer_id_hash', 'date', 'citation_issued']]
daily_timeseries = daily_timeseries.merge(df_officers[['officer_id_hash', 'precinct']], on='officer_id_hash', how='left')
daily_timeseries.head()
weekly_timeseries = df_seattle.groupby([df_seattle.officer_id_hash, pd.Grouper(key='date', freq='W-MON')]).count().reset_index()[['officer_id_hash', 'date', 'citation_issued']]
weekly_timeseries = weekly_timeseries.merge(df_officers[['officer_id_hash', 'precinct']], on='officer_id_hash', how='left')
weekly_timeseries.head()
monthly_timeseries = df_seattle.groupby([df_seattle.officer_id_hash, pd.Grouper(key='date', freq='M')]).count().reset_index()[['officer_id_hash', 'date', 'citation_issued']]
monthly_timeseries = monthly_timeseries.merge(df_officers[['officer_id_hash', 'precinct']], on='officer_id_hash', how='left')
monthly_timeseries.head(10)<jupyter_output><empty_output><jupyter_text>Exports timeseries dataframes of different granularities to the src folder<jupyter_code>daily_timeseries.to_csv(home_directory + '/src/daily_ts.csv', index=False)
weekly_timeseries.to_csv(home_directory + '/src/weekly_ts.csv', index=False)
monthly_timeseries.to_csv(home_directory + '/src/monthly_ts.csv', index=False)<jupyter_output><empty_output><jupyter_text>From the time series dataframes, we create new dataframes on the precinct level<jupyter_code>daily_precinct = daily_timeseries.groupby(['precinct','date']).sum().reset_index()
daily_precinct.head()
weekly_precinct = weekly_timeseries.groupby(['precinct', pd.Grouper(key='date', freq='W-MON')]).sum().reset_index()
weekly_precinct.head()
monthly_precinct = monthly_timeseries.groupby(['precinct', pd.Grouper(key='date', freq='M')]).sum().reset_index()
monthly_precinct.head()
daily_precinct.to_csv(home_directory + '/src/daily_precinct.csv', index=False)
weekly_precinct.to_csv(home_directory + '/src/weekly_precinct.csv', index=False)
monthly_precinct.to_csv(home_directory + '/src/monthly_precinct.csv', index=False)<jupyter_output><empty_output><jupyter_text>Plotting number of monthly citations per precinct<jupyter_code>plt.rcParams['figure.figsize'] = [12, 7]
plt.plot(monthly_precinct[monthly_precinct.precinct == 'E'].date,
         monthly_precinct[monthly_precinct.precinct == 'E'].citation_issued,
         label='E')
plt.plot(monthly_precinct[monthly_precinct.precinct == 'N'].date, 
         monthly_precinct[monthly_precinct.precinct == 'N'].citation_issued,
         label='N')
plt.plot(monthly_precinct[monthly_precinct.precinct == 'W'].date, 
         monthly_precinct[monthly_precinct.precinct == 'W'].citation_issued,
         label='W')
plt.plot(monthly_precinct[monthly_precinct.precinct == 'S'].date, 
         monthly_precinct[monthly_precinct.precinct == 'S'].citation_issued,
         label='S')
plt.plot(monthly_precinct[monthly_precinct.precinct == 'SW'].date, 
         monthly_precinct[monthly_precinct.precinct == 'SW'].citation_issued,
         label='W')
plt.plot(monthly_precinct[monthly_precinct.precinct == 'SW'].date, 
         monthly_precinct.groupby(monthly_precinct.date).mean().citation_issued,
         label='Mean', color='black')
plt.legend()
plt.grid()
plt.title('Citations per Precinct')
plt.savefig(home_directory + '/results/citations_precinct.png')
plt.show()<jupyter_output>/home/ablew/anaconda3/lib/python3.7/site-packages/pandas/plotting/_converter.py:129: FutureWarning: Using an implicitly registered datetime converter for a matplotlib plotting method. The converter was registered by pandas on import. Future versions of pandas will require you to explicitly register matplotlib converters.

To register the converters:
	>>> from pandas.plotting import register_matplotlib_converters
	>>> register_matplotlib_converters()
  warnings.warn(msg, FutureWarning)<jupyter_text>Periodicity We will check periodicity using autocorrelation plots, resampling to smooth the time series<jupyter_code>daily_precinct2 = daily_precinct
daily_precinct2 = daily_precinct2.set_index(daily_precinct2.date)

pd.plotting.autocorrelation_plot(daily_precinct2[(daily_precinct2.precinct == 'E')].citation_issued.resample("1y").mean(), label='E')
pd.plotting.autocorrelation_plot(daily_precinct2[(daily_precinct2.precinct == 'S')].citation_issued.resample("1y").mean(), label='S')
pd.plotting.autocorrelation_plot(daily_precinct2[(daily_precinct2.precinct == 'W')].citation_issued.resample("1y").mean(), label='W')
pd.plotting.autocorrelation_plot(daily_precinct2[(daily_precinct2.precinct == 'N')].citation_issued.resample("1y").mean(), label='N')
pd.plotting.autocorrelation_plot(daily_precinct2[(daily_precinct2.precinct == 'SW')].citation_issued.resample("1y").mean(), label='SW')
plt.savefig(home_directory + '/results/autocorrelation_day.png')
plt.title('Autocorrelation by Day per Precinct')<jupyter_output><empty_output><jupyter_text>Weak negative autocorrelation every 1.5/2 days.<jupyter_code>weekly_precinct2 = weekly_precinct
weekly_precinct2 = weekly_precinct2.set_index(weekly_precinct2.date)

pd.plotting.autocorrelation_plot(weekly_precinct2[(weekly_precinct2.precinct == 'E')].citation_issued.resample("1y").mean(), label='E')
pd.plotting.autocorrelation_plot(weekly_precinct2[(weekly_precinct2.precinct == 'S')].citation_issued.resample("1y").mean(), label='S')
pd.plotting.autocorrelation_plot(weekly_precinct2[(weekly_precinct2.precinct == 'W')].citation_issued.resample("1y").mean(), label='W')
pd.plotting.autocorrelation_plot(weekly_precinct2[(weekly_precinct2.precinct == 'N')].citation_issued.resample("1y").mean(), label='N')
pd.plotting.autocorrelation_plot(weekly_precinct2[(weekly_precinct2.precinct == 'SW')].citation_issued.resample("1y").mean(), label='SW')
plt.savefig(home_directory + '/results/autocorrelation_week.png')
plt.title('Autocorrelation by Week per Precinct')<jupyter_output><empty_output><jupyter_text>Negative autocorrelation every 3 weeks.<jupyter_code>monthly_precinct2 = monthly_precinct
monthly_precinct2 = monthly_precinct2.set_index(monthly_precinct2.date)

pd.plotting.autocorrelation_plot(monthly_precinct2[(monthly_precinct2.precinct == 'E')].citation_issued, label='E')
pd.plotting.autocorrelation_plot(monthly_precinct2[(monthly_precinct2.precinct == 'S')].citation_issued, label='S')
pd.plotting.autocorrelation_plot(monthly_precinct2[(monthly_precinct2.precinct == 'W')].citation_issued, label='W')
pd.plotting.autocorrelation_plot(monthly_precinct2[(monthly_precinct2.precinct == 'N')].citation_issued, label='N')
pd.plotting.autocorrelation_plot(monthly_precinct2[(monthly_precinct2.precinct == 'SW')].citation_issued, label='SW')
plt.savefig(home_directory + '/results/autocorrelation_month.png')
plt.title('Autocorrelation by Month per Precinct')<jupyter_output><empty_output><jupyter_text>No autocorrelation for months.  Hypothetical Quotas For each granularity, we calculate which officers passed the hypothetical ticket quota (range 1-10) at least once.<jupyter_code>pass_list_day = []
for i in range(10):
    day_test = daily_timeseries
    day_test['pass'] = day_test.citation_issued > i
    day_test2 = day_test.groupby('officer_id_hash').mean().reset_index()
    pass_list_day.append(len(day_test2[day_test2['pass'] != 0]) / day_test2.officer_id_hash.nunique())
pass_list_week = []
for i in range(10):
    week_test = weekly_timeseries
    week_test['pass'] = week_test.citation_issued > i
    week_test2 = week_test.groupby('officer_id_hash').mean().reset_index()
    pass_list_week.append(len(week_test2[week_test2['pass'] != 0]) / week_test2.officer_id_hash.nunique())
pass_list_month = []
for i in range(10):
    month_test = monthly_timeseries
    month_test['pass'] = month_test.citation_issued > i
    month_test2 = month_test.groupby('officer_id_hash').mean().reset_index()
    pass_list_month.append(len(month_test2[month_test2['pass'] != 0]) / month_test2.officer_id_hash.nunique())
plt.plot(range(10), pass_list_week, label = 'week')
plt.plot(range(10), pass_list_month, label = 'month')
plt.plot(range(10), pass_list_day, label = 'day')
plt.grid()
plt.title('Proportion of officers passing quota')
plt.ylabel('Proportion of officers')
plt.xlabel('Hypothetical ticket quota')
plt.legend()
plt.savefig(home_directory + '/results/passing_officers.png')
plt.show()<jupyter_output><empty_output><jupyter_text>As we can see, the proportion drops off rapidly after the 1-ticket quota. If we assume that majority = 60%, the majority of police officers won't pass the traffic ticket quota when the quota is set to 2 tickets per day, per week or per month. Even if the quota was 1 per day, it would follow by linearity that officers would pass 7 per week or 30 per month, but none of those cases occur.   Creates interactive heatmap in an html file in the results folder. Please input Google Maps API Key obtained from https://cloud.google.com/maps-platform/pricing/<jupyter_code>latitudes = df_seattle[(df_seattle.citation_issued == True)].lat
longitudes = df_seattle[(df_seattle.citation_issued == True)].lng

# Creating the location we would like to initialize the focus on. 
# Parameters: Lattitude, Longitude, Zoom
gmap = gmplot.GoogleMapPlotter(47.596971, -122.327652, 12)

# Overlay our datapoints onto the map
gmap.heatmap(latitudes, longitudes)

# Input API Key here
gmap.apikey = 'ENTER-API-KEY-HERE'

# Generate the heatmap into an HTML file
gmap.draw(home_directory + "/results/my_heatmap.html")<jupyter_output><empty_output><jupyter_text>Creating a dataframe grouped by weekday and plotting histogram<jupyter_code>df_weekday = df_seattle[df_seattle.citation_issued == True].groupby(df_seattle.date.dt.weekday_name).count()
df_weekday = df_weekday.reindex(index = ['Monday','Tuesday','Wednesday', 'Thursday', 'Friday', 'Saturday', 'Sunday'])
ax = df_weekday.citation_issued.plot(kind='bar', color='blue')
plt.title('Citations per Weekday')
plt.ylabel('Citations Issued')
plt.xlabel('Day')
plt.savefig(home_directory + '/results/Citations_per_weekday.png')
plt.show()<jupyter_output><empty_output><jupyter_text>Creating an hourly time series for displaying common hours for traffic stops<jupyter_code>hourly_timeseries = df_seattle[df_seattle.citation_issued == True].groupby(pd.Grouper(key='time', freq='60min')).count()
hourly_timeseries.head()
ax = hourly_timeseries.citation_issued.plot(kind='bar', color='blue')
ticklabels = hourly_timeseries.index.strftime('%H:%Mh')
ax.xaxis.set_major_formatter(matplotlib.ticker.FixedFormatter(ticklabels))
plt.title('Citations per hour')
plt.ylabel('Citations')
plt.savefig(home_directory + '/results/Citations_per_hour.png')
plt.show()<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.1351634793337446, 'fraction_numerical': 0.017088217149907464, 'mean_word_length': 7.93168044077135, 'pattern_counts': {'":': 0, '<': 78, '<?xml version=': 0, '>': 88, 'https://': 1, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '20', 'dir': 'jupyter-structured-clean-dedup', 'id': '46666', 'license': 'MIT', 'n_tokens_mistral': 6030, 'n_tokens_neox': 5637, 'n_words': 932, 'path': 'src/final-project-notebook.ipynb', 'repo_name': 'ablew/data-512-final-project', 'size': '16210'}