data_cleanup.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# This notebook cleans and filters the data\n",
    "Parse/convert the times from human readable to machine readable (e.g. \"after school lets out\" --> \"3:00PM\")\n",
    "Parse the addresses to City, Zipcode (when available), State, Country\n",
    "Use external lookup and custom city/state corrections to determine County for each address"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Imports"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Import dependencies\n",
    "import pandas as pd\n",
    "import numpy as np\n",
    "import os,sys\n",
    "import json\n",
    "import re\n",
    "import matplotlib.pyplot as plt\n",
    "from datetime import date, datetime\n",
    "import dateparser\n",
    "from pprint import pprint\n",
    "\n",
    "DATA_DIR = \"data\""
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Lookup for City->County and County->CountyFIPS matchings that have mismatching city names"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "city_reps = {\"AR\": {\"Helena-West Helena\":\"West Helena\",\n",
    "                    \"DeWitt\":\"Dewitt\",\n",
    "                    \"Fargo\":\"Brinkley\"},\n",
    "            \"CA\": {\"South Whittier\":\"Whittier\",\n",
    "                   \"Florence-Graham\":\"Los Angeles\",\n",
    "                   \"San Antonio\":\"San Rafael\",\n",
    "                   \"Willowbrook\":\"Watts\",\n",
    "                   \"Florence-Graham\":\"Oakwood\"},\n",
    "            \"DE\": {\"Pennyhill\":\"Bear\"},\n",
    "            \"FL\": {\"Opa-locka\":\"Opa Locka\"},\n",
    "            \"GA\": {\"LaGrange\":\"Lagrange\"},\n",
    "            \"IA\": {\"DeWitt\":\"Dewitt\"},\n",
    "            \"IL\": {\"Glendale\":\"Cloverdale\",\n",
    "                   \"DeKalb\":\"Dekalb\"},\n",
    "            \"KS\": {\"Alta Mills\":\"Newton\"},\n",
    "            \"LA\": {\"Phoenix\":\"Buras\"},\n",
    "            \"MA\": {\"Manchester-by-the-Sea\":\"Manchester By The Sea\"},\n",
    "            \"MD\": {\"Greater Landover\":\"Landover\",\n",
    "                   \"Fair Hill\":\"Perryville\"},\n",
    "            \"MI\": {\"Mt Morris Township\":\"Mount Morris\",\n",
    "                   \"Calumet Twp\":\"Calumet\",\n",
    "                   \"Kalamazoo Twp\":\"Kalamazoo\",\n",
    "                   \"Charter Twp of Clinton\":\"Armada\",\n",
    "                   \"Bridgeport Charter Township\":\"Bridgeport\"},\n",
    "            \"MT\": {\"Bonner-West Riverside\":\"Bonner\"},\n",
    "            \"NC\": {\"McLeansville\":\"Mcleansville\",\n",
    "                   \"Feltonville\":\"Warrenton\"},\n",
    "            \"NJ\": {\"Lafayette Township\":\"Lafayette\",\n",
    "                   \"Bethlehem\":\"Califon\",\n",
    "                   \"Upper Deerfield Township\":\"Upper Deerfield Twp\",\n",
    "                   \"White\": \"Allamuchy\"},\n",
    "            \"NY\": {\"Hartland\":\"Burt\"},\n",
    "            \"OK\": {\"McAlester\":\"Mcalester\"},\n",
    "            \"PA\": {\"West Conshohocken\":\"Conshohocken\",\n",
    "                   \"McKeesport\":\"Mckeesport\",\n",
    "                   \"Bensalem Township\":\"Bensalem\",\n",
    "                   \"Wilkes-Barre\":\"Wilkes Barre\",\n",
    "                   \"Warwick Township\":\"Deodate\"},\n",
    "            \"SC\": {\"Dacusville\":\"Central\"},\n",
    "            \"TX\": {\"Falcon Mesa\":\"Falcon\"},\n",
    "            \"WI\": {\"DeForest\":\"De Forest\",\n",
    "                   \"Ashford\":\"Eden\",\n",
    "                   \"Mt Morris\":\"Coloma\"},\n",
    "            \"WV\": {\"Davis Creek\":\"Belle\"}\n",
    "            }\n",
    "\n",
    "county_reps = {\n",
    "            \"IL\": {\"DE WITT\": \"DEWITT\"}\n",
    "}"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Create a map from (County Name) to (State and County FIPS IDs) so we can uniquely and universally reference counties in outside datasets\n",
    "\n",
    "Source: https://www.census.gov/geo/reference/codes/cou.html"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "counties_fp_raw = pd.read_csv(os.path.join(DATA_DIR, \"national_county.txt\"), header=None, names=[\n",
    "    \"State\", \"StateFP\", \"CountyFP\", \"County\", \"ClassFP\"])\n",
    "counties_fp_map = {}\n",
    "for index, row in counties_fp_raw.iterrows():\n",
    "    state = row[\"State\"]\n",
    "    state_fp = str(row[\"StateFP\"])\n",
    "    county = row[\"County\"]\n",
    "    county_fp = str(row[\"CountyFP\"])\n",
    "    \n",
    "    state_fp = \"0\"*(2-len(state_fp)) + state_fp\n",
    "    county_fp = \"0\"*(3-len(county_fp)) + county_fp  \n",
    "    \n",
    "    county = county.upper() \\\n",
    "        .replace(\" COUNTY\", \"\") \\\n",
    "        .replace(\".\", \"\") \\\n",
    "        .replace(\"'\", \"\") \\\n",
    "        .replace(\"-\", \" \")\n",
    "    \n",
    "    if state in county_reps and county in county_reps[state]:\n",
    "        county = county_reps[state][county]\n",
    "    \n",
    "    if state not in counties_fp_map:\n",
    "        counties_fp_map[state] = {}\n",
    "        \n",
    "    counties_fp_map[state][county] = {\"state\":state_fp, \"county\":county_fp}"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Create map from State/City to County using external data\n",
    "Source: https://github.com/grammakov/USA-cities-and-states/blob/master/us_cities_states_counties.csv"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "counties_map_raw = pd.read_csv(os.path.join(DATA_DIR, \"us_cities_states_counties.csv\"), delimiter=\"|\")\n",
    "counties_map_raw = counties_map_raw.where((pd.notnull(counties_map_raw)), None)\n",
    "counties_fp_raw = counties_fp_raw.where((pd.notnull(counties_fp_raw)), None)\n",
    "\n",
    "counties_map = {}\n",
    "counties_map_alias = {}\n",
    "for index, row in counties_map_raw.iterrows():\n",
    "    state = row[\"State short\"]\n",
    "    city = row[\"City\"]\n",
    "    city_alias = row[\"City alias\"]\n",
    "    county_name = row[\"County\"]\n",
    "    \n",
    "    if not state or not county_name or state in [\"AE\", \"PW\", \"AS\", \"FM\", \"MP\", \"MH\"]:\n",
    "        continue\n",
    "    try:\n",
    "        q_county = [county_name, county_name.replace(\"-\", \" \"), \n",
    "                    county_name + \" CITY\", county_name.replace(\"SAINT \", \"ST \"),\n",
    "                    county_name.replace(\"SAINT \", \"ST \") + \" ISLAND\",\n",
    "                    county_name.replace(\"SAINT \", \"ST \") + \" PARISH\",\n",
    "                    county_name.replace(\"SAINTE \", \"STE \"),\n",
    "                    county_name + \" MUNICIPIO\", county_name + \" MUNICIPALITY\",\n",
    "                    county_name + \" CENSUS AREA\", county_name + \" BOROUGH\",\n",
    "                    county_name + \" CITY AND BOROUGH\",\n",
    "                    county_name + \" ISLAND\", county_name.replace(\" \", \"\")]\n",
    "    except:\n",
    "        print(county_name)\n",
    "    fp = None\n",
    "    for c in q_county:\n",
    "        if c in counties_fp_map[state]:\n",
    "            fp = counties_fp_map[state][c]\n",
    "            \n",
    "    if not fp:\n",
    "        print(state)\n",
    "        print(county_name)\n",
    "        continue\n",
    "        \n",
    "    state_fp = str(fp[\"state\"])\n",
    "    county_fp = str(fp[\"county\"])\n",
    "    state_fp = \"0\"*(2-len(state_fp)) + state_fp\n",
    "    county_fp = \"0\"*(3-len(county_fp)) + county_fp\n",
    "    \n",
    "    if state not in counties_map:\n",
    "        counties_map[state] = {}\n",
    "    if city not in counties_map[state]:\n",
    "        counties_map[state][city] = {\"county_name\":county_name, \"county_fp\":county_fp, \"state_fp\":state_fp}\n",
    "    \n",
    "    if state not in counties_map_alias:\n",
    "        counties_map_alias[state] = {}\n",
    "    if city_alias not in counties_map[state]:\n",
    "        counties_map_alias[state][city_alias] = {\"county_name\":county_name, \"county_fp\":county_fp, \"state_fp\":state_fp}\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Create Map for County Lookup"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "county_data_map = {}"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Population"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "pop_data_raw = pd.read_excel(os.path.join(DATA_DIR, \"Population_2010_2016.xls\"), skiprows=[0,1], header=0)\n",
    "for index, row in pop_data_raw.iterrows():\n",
    "    fp = str(row[\"FIPS\"])\n",
    "    fp = \"0\"*(5-len(fp)) + fp\n",
    "    state = row[\"State\"]\n",
    "    if state in [\"PR\"]:\n",
    "        continue\n",
    "    if county_fp not in county_data_map:\n",
    "        county_data_map[fp] = {}\n",
    "    if \"pop\" not in county_data_map[fp]:\n",
    "        county_data_map[fp][\"pop\"] = {}\n",
    "    pop_key = \"POP_ESTIMATE_\"\n",
    "    for y in range(2010, 2017):\n",
    "        pop = row[pop_key + str(y)]\n",
    "        county_data_map[fp][\"pop\"][y] = str(pop)\n",
    "    # Account for 2017 pop by adding half of the change from 2015-2016\n",
    "    pop_2016 = int(row[pop_key + \"2016\"])\n",
    "    est_chg_2017 = 0.5*int(row[\"N_POP_CHG_2016\"])\n",
    "    county_data_map[fp][\"pop\"][2017] = str(int(pop_2016 + est_chg_2017))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Land Size"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "land_data_raw = pd.read_csv(os.path.join(DATA_DIR, \"2010 Geo Data.csv\"), encoding=\"ISO-8859-1\", skiprows=[0], header=0)\n",
    "for index, row in land_data_raw.iterrows():\n",
    "    fp = row[\"Target Geo Id2\"]\n",
    "    if np.isnan(fp):\n",
    "        fp = \"00000\"\n",
    "    fp = str(int(fp))\n",
    "    if len(fp) < 3:\n",
    "        fp = fp + \"000\"\n",
    "    fp = \"0\"*(5-len(fp)) + fp\n",
    "    if fp not in county_data_map:\n",
    "        county_data_map[fp] = {}\n",
    "    land = row[\"Area in square miles - Land area\"]\n",
    "    county_data_map[fp][\"land\"] = land"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Population Density (Calculated)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Mean Density: 268\n",
      "# of Counties: 3192\n",
      "# of High Density Counties: 910\n"
     ]
    }
   ],
   "source": [
    "densities = []\n",
    "for fp, data in county_data_map.items():\n",
    "    if \"pop\" not in data or \"land\" not in data:\n",
    "        continue\n",
    "    pop_data = data[\"pop\"]\n",
    "    land = data[\"land\"]\n",
    "    pops = []\n",
    "    for y, pop in pop_data.items():\n",
    "        pops.append(int(pop))\n",
    "    mean_pop = np.mean(pops)\n",
    "    density = mean_pop / land\n",
    "    densities.append(density)\n",
    "    county_data_map[fp][\"pop_density\"] = density\n",
    "mean_density = np.mean(densities)\n",
    "high_density = [density for density in densities if density > 100]\n",
    "print(\"Mean Density: %d\" % int(mean_density))\n",
    "print(\"# of Counties: %d\" % len(densities))\n",
    "print(\"# of High Density Counties: %d\" % len(high_density))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Unemployment\n",
    "Represented as Unemployment Rate"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [],
   "source": [
    "unemp_data_raw = pd.read_csv(os.path.join(DATA_DIR, \"unemployment.csv\"), encoding=\"ISO-8859-1\", header=0)\n",
    "for index, row in unemp_data_raw.iterrows():\n",
    "    state_fp = str(int(row[\"State\"]))\n",
    "    state_fp = \"0\"*(2-len(state_fp)) + state_fp\n",
    "    county_fp = str(int(row[\"County\"]))\n",
    "    county_fp = \"0\"*(3-len(county_fp)) + county_fp\n",
    "    fp = state_fp + county_fp\n",
    "    if fp not in county_data_map:\n",
    "        county_data_map[fp] = {}\n",
    "    if \"unemp\" not in county_data_map[fp]:\n",
    "        county_data_map[fp][\"unemp\"] = {}\n",
    "    for y in range(2010, 2018):\n",
    "        y_i = str(y) + \" UR\"\n",
    "        unemp = row[y_i]\n",
    "        county_data_map[fp][\"unemp\"][y] = unemp"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Education\n",
    "Represented as (100 - Percent of adults with less than a high school diploma, 2012-2016)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [],
   "source": [
    "edu_data_raw = pd.read_csv(os.path.join(DATA_DIR, \"Education_12_16.csv\"), encoding=\"ISO-8859-1\", header=0)\n",
    "for index, row in edu_data_raw.iterrows():\n",
    "    fp = str(int(row[\"FIPS Code\"]))\n",
    "    fp = \"0\"*(5-len(fp)) + fp\n",
    "    edu = 100 - row[\"Percent of adults with less than a high school diploma, 2012-2016\"]\n",
    "    if fp not in county_data_map:\n",
    "        county_data_map[fp] = {}\n",
    "    county_data_map[fp][\"edu\"] = edu"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Poverty\n",
    "Assume 2016 poverty rate is representative of 2010-2017"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [],
   "source": [
    "edu_data_raw = pd.read_csv(os.path.join(DATA_DIR, \"poverty_2016.csv\"), encoding=\"ISO-8859-1\", header=0)\n",
    "for index, row in edu_data_raw.iterrows():\n",
    "    fp = str(int(row[\"FIPStxt\"]))\n",
    "    fp = \"0\"*(5-len(fp)) + fp\n",
    "    pov_str = row[\"POVALL_2016\"]\n",
    "    pov_str = str(pov_str)\n",
    "    if pov_str == 'nan':\n",
    "        continue\n",
    "    pov_str = pov_str.replace(\",\", \"\")\n",
    "    num_poverty = int(pov_str)\n",
    "    pop_2016 = county_data_map[fp][\"pop\"][2016]\n",
    "    poverty_rate = num_poverty/int(pop_2016)\n",
    "    county_data_map[fp][\"poverty\"] = poverty_rate"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Parse human-readable address and determine counties / FIPS IDs"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [],
   "source": [
    "def parse_address(address):\n",
    "    if not address:\n",
    "        return None\n",
    "    a = [_.strip() for _ in address.split(\",\")]\n",
    "    cur = -1\n",
    "    country = a[cur]\n",
    "    cur -= 1\n",
    "    zcode = None\n",
    "    state = None\n",
    "    city = None\n",
    "    # State and Zip have no comma\n",
    "    if -cur <= len(a):\n",
    "        if len(a[cur]) > 5:\n",
    "            state_zip = a[cur].split(\" \")\n",
    "            if len(state_zip) == 2:\n",
    "                zcode = state_zip[1]\n",
    "                state = state_zip[0]\n",
    "                cur -= 1\n",
    "        # Has zip\n",
    "        elif len(a[cur]) > 2: \n",
    "            if -cur <= len(a):\n",
    "                zcode = a[cur]\n",
    "                cur -= 1\n",
    "            if -cur <= len(a):\n",
    "                state = a[cur]\n",
    "                cur -= 1\n",
    "        # Doesn't have zip\n",
    "        else:\n",
    "            if -cur <= len(a):\n",
    "                state = a[cur]\n",
    "                cur -= 1\n",
    "        \n",
    "    if -cur <= len(a):\n",
    "        city = a[cur]\n",
    "    \n",
    "    #### Get County ####\n",
    "    county = None\n",
    "    city_orig = city\n",
    "    \n",
    "    if city:\n",
    "        # Cleanup city name\n",
    "        if state in city_reps and city in city_reps[state]:\n",
    "            city = city_reps[state][city]\n",
    "\n",
    "        city = city.replace(\"St.\", \"Saint\")\n",
    "        city = city.replace(\"Ste.\", \"Sainte\")\n",
    "        city = city.replace(\"ñ\", \"n\")\n",
    "        city = city.replace(\"-\", \" \")\n",
    "        if city[:3] == \"Mt \":\n",
    "            city = \"Mount \" + city[3:]\n",
    "        if city[:3] == \"St \":\n",
    "            city = \"Saint \" + city[3:]\n",
    "\n",
    "        if state in counties_map:\n",
    "            if city in counties_map[state]:\n",
    "                county = counties_map[state][city]\n",
    "            elif city in counties_map_alias[state]:\n",
    "                county = counties_map_alias[state][city]\n",
    "\n",
    "    return county, city_orig, state, zcode, country"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Parse time and convert to machine-readable format since the data is very unstructured and inconsistent"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [],
   "source": [
    "def parse_time(fts, cts, recurse=True):\n",
    "    fine_time = \"\"\n",
    "    coarse_time = \"\"\n",
    "    orig_fts = fts\n",
    "    orig_cts = cts\n",
    "    \n",
    "    stop_chars = ['.', ',']\n",
    "    days = [\"monday\", \"tuesday\", \"wednesday\", \n",
    "            \"thursday\", \"friday\", \"saturday\", \"sunday\"]\n",
    "    \n",
    "    numbers = [\"one\", \"two\", \"three\", \"four\", \"five\", \"six\", \n",
    "               \"seven\", \"eight\", \"nine\", \"ten\", \"eleven\", \"twelve\"]\n",
    "    \n",
    "    # Get fine time from clock-time\n",
    "    fts = fts.lower()\n",
    "    fts = fts.replace(\"p. m.\", \"p.m.\")\n",
    "    fts = fts.replace(\"a.m.\", \"a.m.\")\n",
    "    fts = fts.strip()\n",
    "    fts = re.sub(r\"([0-9]+)\\.([0-9]+)\", r\"\\1:\\2\", fts)\n",
    "    fts = re.sub(r\"([0-9]+)[ ]*\\:[ ]*([0-9]+)\", r\"\\1:\\2\", fts)\n",
    "    for i, number in enumerate(numbers):\n",
    "        fts = fts.replace(number, str(i+1))\n",
    "    \n",
    "    while fts and not fts[0].isdigit():\n",
    "        fts = fts[1:]\n",
    "    if fts:\n",
    "        try:\n",
    "            fine_time = dateparser.parse(fts).time()\n",
    "        except:\n",
    "            # Just keep going if you couldn't figure it out\n",
    "            _ = True\n",
    "            #print(\"Fine Time Parse Error 1: '\" + orig_fts + \"'  |  '\" + orig_cts + \"'\")\n",
    "        \n",
    "    # Get coarse time from time-day\n",
    "    cts = cts.lower()\n",
    "    coarse_time = datetime.now().replace(hour=0, minute=0, second=0, microsecond=0)\n",
    "    \n",
    "    for day in days:\n",
    "        cts = cts.replace(day, '')\n",
    "    for c in stop_chars:\n",
    "        cts = cts.replace(c, '')\n",
    "    cts = cts.replace('  ', ' ')\n",
    "    cts = cts.strip()\n",
    "    if cts == \"\":\n",
    "        coarse_time = \"\"\n",
    "    elif cts in [\"early morning\", \"dawn\"] or \\\n",
    "        (\"early\" in cts and \"morning\" in cts) or \"wake\" in cts:\n",
    "        coarse_time = coarse_time.replace(hour=7)\n",
    "    elif (\"early\" in cts and \"afternoon\" in \"cts\") or \"noon\" in cts or \"lunchtime\" in cts:\n",
    "        coarse_time = coarse_time.replace(hour=12)\n",
    "    elif cts in [\"morning\", \"am\", \"am\", \"early\"]:\n",
    "        coarse_time = coarse_time.replace(hour=10)\n",
    "    elif cts in [\"afternoon\", \"pm\", \"mid-day\", \"midday\"] or \"after school\" in cts:\n",
    "        coarse_time = coarse_time.replace(hour=15)\n",
    "    elif cts in [\"night\", \"evening\", \"late\"] or \"night\" in cts or \"nigh\" in cts \\\n",
    "        or \"after business\" in cts or \"dinnertime\" in cts:\n",
    "        coarse_time = coarse_time.replace(hour=21)\n",
    "    elif cts == \"midnight\" or \"midnight\" in cts:\n",
    "        coarse_time = coarse_time.replace(hour=0)\n",
    "    elif cts in [\"overnight\", \"before dawn\"]:\n",
    "        #TODO: Make this the next day (maybe that doesn't matter)\n",
    "        # .    Could be good enough to indicate \"3am\" or something\n",
    "        # .    since the day doesn't matter when examining time of day\n",
    "        # .    as an independent factor\n",
    "        coarse_time = coarse_time.replace(hour=3)\n",
    "    # Repeat this later in the priority list in case we hit \"early afternoon\" or something\n",
    "    elif \"early\" in cts:\n",
    "        coarse_time = coarse_time.replace(hour=10)\n",
    "    else:\n",
    "        coarse_time = \"\"\n",
    "        if not fine_time:\n",
    "            while cts and not cts[0].isdigit():\n",
    "                cts = cts[1:]\n",
    "            if cts:\n",
    "                try:\n",
    "                    fine_time = dateparser.parse(cts).time()\n",
    "                except:\n",
    "                    # Just keep going if you couldn't figure it out\n",
    "                    _ = True\n",
    "                    #print(\"Coarse Time Parse Error 1: '\" + orig_fts + \"'  |  '\" + orig_cts + \"'\")\n",
    "    # If no time was found, switch the two and try to parse\n",
    "    # This helps because sometimes the coarse time has something like \"2 a.m\" \n",
    "    # and the fine time has something like \"midnight\"\n",
    "    if coarse_time:\n",
    "        coarse_time = coarse_time.time()\n",
    "    if recurse and not fine_time and not coarse_time:\n",
    "        fine_time, coarse_time = parse_time(orig_cts, orig_fts, recurse=False)\n",
    "    return fine_time, coarse_time"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Read data and use parsers to clean it up"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [],
   "source": [
    "events_raw = pd.read_csv(os.path.join(DATA_DIR, \"Events.tsv\"), delimiter=\"\\t\")\n",
    "events_raw.columns = [\"Address\", \"Date\", \"Details\", \"Participants\"]\n",
    "events = pd.DataFrame(columns=[\"Address\", \"Details\", \"Country\", \"State\", \"StateFP\",\n",
    "                               \"County\",\"CountyFP\", \"City\", \"Zipcode\", \n",
    "                               #\"Population\", \"Land Area\", \"Population Density\",\n",
    "                               \"Date\", \"Fine Time\", \"Coarse Time\", \n",
    "                               \"Shots Fired\", \"Gun Type\", \"Participants\"])\n",
    "\n",
    "events[\"Date\"] = events.Date.astype(datetime)\n",
    "events[\"Zipcode\"] = events.Zipcode.astype(int)\n",
    "for index, event in events_raw.iterrows():\n",
    "    new_event = {}\n",
    "    new_event[\"Address\"] = event[\"Address\"]\n",
    "    county, city, state, zcode, country = parse_address(event[\"Address\"])\n",
    "    new_event[\"Country\"] = country\n",
    "    new_event[\"State\"] = state\n",
    "    new_event[\"City\"] = city\n",
    "    new_event[\"Zipcode\"] = zcode\n",
    "    \n",
    "    new_event[\"Date\"] = datetime.strptime(event[\"Date\"][:-6], \n",
    "                                          \"%Y-%m-%d %H:%M:%S\").date()\n",
    "    # Incorrect date entry\n",
    "    if new_event[\"Date\"].year == 2026:\n",
    "        new_event[\"Date\"] = new_event[\"Date\"].replace(year=2016)\n",
    "        \n",
    "    if county:\n",
    "        county_fp = county[\"county_fp\"]\n",
    "        state_fp = county[\"state_fp\"]\n",
    "        new_event[\"StateFP\"] = state_fp\n",
    "        new_event[\"County\"] = county[\"county_name\"]\n",
    "        new_event[\"CountyFP\"] = county_fp\n",
    "        year = new_event[\"Date\"].year\n",
    "    \n",
    "        '''\n",
    "        if year >= 2010 and year <= 2017:\n",
    "            fp = str(state_fp) + str(county_fp)\n",
    "            fp = \"0\"*(5-len(fp)) + fp\n",
    "            new_event[\"Population\"] = county_data_map[fp][\"pop\"][year]\n",
    "            new_event[\"Land Area\"] = county_data_map[fp][\"land\"]\n",
    "            new_event[\"Population Density\"] = county_data_map[fp][\"pop_density\"]            \n",
    "            new_event[\"Education\"]\n",
    "        '''\n",
    "        \n",
    "    event_details = json.loads(event[\"Details\"])\n",
    "    fine_time, coarse_time = parse_time(event_details[\"clock-time\"][\"value\"], \n",
    "                                    event_details[\"time-day\"][\"value\"])\n",
    "    new_event[\"Details\"] = event_details[\"details\"][\"value\"]\n",
    "    new_event[\"Fine Time\"] = fine_time\n",
    "    new_event[\"Coarse Time\"] = coarse_time\n",
    "    new_event[\"Shots Fired\"] = event_details[\"number-of-shots-fired\"][\"value\"]\n",
    "    new_event[\"Gun Type\"] = event_details[\"type-of-gun\"][\"value\"]\n",
    "    new_event[\"Participants\"] = event[\"Participants\"]\n",
    "    events = events.append(new_event, ignore_index=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Plot data to determine how to filter"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Collect data on dates and times of shootings\n",
    "# Prefer Fine Time to Coarse Time if the record has it\n",
    "event_dates = []\n",
    "event_times = []\n",
    "for index, event in events.iterrows():\n",
    "    event_dates.append(event[\"Date\"])\n",
    "    if event[\"Fine Time\"]:\n",
    "        event_times.append(event[\"Fine Time\"])\n",
    "    else:\n",
    "        event_times.append(event[\"Coarse Time\"])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "3524 total events\n",
      "942 events with no time\n"
     ]
    }
   ],
   "source": [
    "# Create bar charts of # events by year and # events by time of day\n",
    "num_events_by_year = []\n",
    "years = range(2000, 2018)\n",
    "for year in years:\n",
    "    start = date(year, 1, 1)\n",
    "    end = date(year, 12,31)\n",
    "    events_in_year = len([d for d in event_dates if d >= start and d <= end])\n",
    "    num_events_by_year.append(events_in_year)\n",
    "    \n",
    "num_events_by_tod = []\n",
    "time_ranges = [(0, 6), (6, 12), (12, 18), (18, 0)]\n",
    "for tr in time_ranges:\n",
    "    start = datetime.now().replace(hour=tr[0], minute=0, second=0, microsecond=0).time()\n",
    "    end = datetime.now().replace(hour=(tr[1]-1) % 24, minute=59, second=59, microsecond=0).time()\n",
    "    events_in_time_range = len([t for t in event_times if t and t >= start and t <= end])\n",
    "    num_events_by_tod.append(events_in_time_range)\n",
    "    \n",
    "num_empty_times = len([t for t in event_times if not t])\n",
    "\n",
    "print(\"%d total events\" % len(event_times))\n",
    "print(\"%d events with no time\" % num_empty_times)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Text(0,0.5,'Number of Gun Violence Events')"
      ]
     },
     "execution_count": 17,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAxUAAAFNCAYAAAB2X0K5AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvhp/UCwAAIABJREFUeJzs3XmcZFV9///XW0ZQZBmWCUG2AUUN7jhR3I1EBSRiXFGDoBiiUWNcohhNMC4RY1x/GggRFAyCSIwSxSgBkfhVkEVAQJEJojMjCMgiLoDo5/fHPS1F091T01XV3UW/no9HPerec7dzb926nzp1zzk3VYUkSZIkzdbd5jsDkiRJksabhQpJkiRJA7FQIUmSJGkgFiokSZIkDcRChSRJkqSBWKiQJEmSNBALFYtYkk8keec8bTtJPp7k+iTfmo88rKskX0qyfx/zLU9SSZbMRb7uipIckOTr850PaVx4PV83Xs/nzlxdz5O8KMlXRr2dPvKxVZIzktyU5H3znZ+5ZKFiAUlyRZKrk9yrJ+1lSU6fx2yNyuOApwDbVtUjp5ohydZJ/i3Jj5P8PMnlLXA+YNiZSfK9JC+dIv01Sc4BqKo9q+roYW97LrSg+It2HCdebxzxNue9YJDknUlOnZR2vyQ/S/Lg+cqX7vq8nt+R1/PhWcTX84t79vc3SW7uGf/bqjq2qp46n3lsDgKuBTapqtdPntjO+1tboeOmJBcleXeSTec+q8NloWLhWQ94zXxnYl0lWW8dF9kBuKKqfjHN+rYAvgFsCDwe2BjYFfgaXfAatqOBF0+Rvl+bdlfw0KraqOf1T/OdoTnwDuD3k/w5dP+oAv8GvL+qvjPMDflPpqbg9Ryv5yOy6K7nVfXAif0F/hd4Vc/+/+N856/HDsAlNfPTpf+pqjYGlgEvAXYD/l/vnxBjqap8LZAXcAVwMHAdsLSlvQw4vQ0vBwpY0rPM6cDL2vABwP8DPgDcAFwOPKalrwKuBvbvWfYTwOHAKcBNdBf4HXqmP6BNuw64FHjepGUPA04GfgH88RT7c2/gpLb8SuDPW/qBwM3Ab4CfA/8wxbLvBC4A7jbD8XoSsHqKY/jHbfhtwAnAMW3/LgZWTLOubYHbJu3/LsCtwJZTHOu7AW8FftiO6zHAplN9TsCmwJHAlcCatm/r9XxmXwf+Gbge+AGwZ08eNgc+Dvy4Tf9cz7S9gfPbZ/0N4CEzHKsC7jvNZ/QrYPOetIfT/cty9zb+UuC7bftfnnSMCng5cFnLx0eBAH8w6TO+oc2/F3BJ+zzWAG+YJr8H0J3LHwFuBL4H7N6mPRc4d9L8rwM+P826HgX8tO3rXwAXTuxbz3fse23/vgRs1zPtI8Bq4GfA2cBjJp2jnwaOa/tzwHxfQ3wtnBdez3uX9XpeXs8ZwvV8qu/KpO18fdL+/GXbn5vo/mS6Tzu+P2vn0/qz/AweQxcTbqQnNtB9l35Nd679nKm/S58A3jkpbeN2Tr2qjd8HOI0udl0LHMvt15G/Af5j0vIfBj407OvYur7m/cLrq+fDaBdQ4LMTJxzrHoRuoyv1rkd3sftRuzBsADy1fbE2avN/oo0/oU3/0MQXErgXXeB6CbCE2y9Mu/QseyPwWLoL8j2m2J8zgH8B7gE8DLgGeHJPXr8+w7E4E3jbWo7Xk1h7ELqZ7sK3HvBu4MwZ1ncK8Nae8Xdzx4t+77F+KV1g3QnYqH1mn5zqcwL+E/jXdkx/D/gW8Bc9x+HXwJ+3PL6CLuCkTf8i3Q/XzYC7A09s6Q+nC36Pasvt3/Z9g2n2bcog1KadRvuB0MbfCxzehvdp+/kH7Tx4K/CNSev9ArAU2L59xntM9xnTXTQf34Y3A3adJk8H0J3Lr237/fx2vm1Od65eB/xBz/zfBp49w2f7PuBUunN4RU/6s+l+YN2/7d/bgP/tmb5f2+YS4E10gXODNu2ddIHjT+i+A/ec72uIr4Xzwut577Jez8vrOUO6nk/+/CZtZ3Kh4vPAJsADgVvo4sBOdIXDS2gF83X5DFq+r6eLD0uAF7TxLXq+T++cIe9TTqcrzH66Dd+X7i7eBnR3M84APtimbU1X+J8oZCxpeX/EqK5n/b7m/cLrq+fDuD0IPah94Zax7kHosp5pD27zb9WT9lPgYW34E8DxPdM2ovsnYrv2pf/fSfn7V+CQnmWPmWFftmvr2rgn7d3AJ3ryOlMQWgm8vGf8GXT/HtwEfKWlPYm1B6H/6Zm2C/CrGbb5Z8ClbfhudAH8T6c51qcCf9kz7f50wWRJ7+cEbEV3Ibtnz7wvAL7acxxW9kzbsC37++3C8Vtgsynyehjwjklpl9KC1BTzF90/Mzf0vJ7Wpr0MOK0Nh+7HxxPa+JeAA3vWczfgl7R/t9p6H9cz/QTg4Ok+43ZM/4KurulM34UD6AnGLe1bwH49+/+uNvxAugv6lAG4zXNPun8hPzAp/RTu+G/vkvZ5bTPFOtLOvwe28XdOHDdfvia/8Hreu7zXc6/nQ7ueT/78Jm1ncqHisT3j5wJv6hl/H7f/UO/7M6ArTHxrUto3aXermX2h4lDglGmWeSbw7Z7xL3H73cK96apbzft1zzYVC1BVXUT3b8HBs1j8Jz3Dv2rrm5y2Uc/4qp7t/pzuH4N709UJfFSSGyZewIvoLo53WnYK9wauq6qbetJ+CGzT5378lO4iPJG3k6pqKd0/Hev3uQ6Aq3qGfwncY4a6758Ftk6yG12A25Dun6Wp3Jtufyb8kNuDTq8d6P6ZubLnOP4r3T9cd8pjVf2yDW5EF8ivq6rrp9j+DsDrJ30+27V8TWfXqlra8/pyS/8P4NFJtqb7l/O3dPVVJ7bzoZ5tXEcXqHo/x8nHuPf8muzZdP80/jDJ15I8eoZ511S7YjY/7Nm/o4EXtjYS+wEnVNUt062oqn5FVxXh4kmTdgA+2rN/19Lt/7YASd7YGn3eSBfo7gVs2bP8TN8Byet5x+u51/OhXc/X0eTvy3Tfn3X5DCafL7Bu34fpbEP3mUz0IHV8kjVJfgb8O3eMPUfTFZxp758ccNtDYaFi4TqE7hZq70k60Qhuw5603qAwG9tNDCTZiO623o/pAszXJl20NqqqV/QsW0zvx8DmSTbuSduervpIP04FnplkpnP0F/Qci9a4cFmf67+TFgBOpGvgtx/dv363TjP7j+kuQhO2p7u9+5NJ862i+2dry57juElVPbCPLK2iO4ZLp5n2rkmfz4ZVdVwf672DFuS+Qvdv5gvp9nvis11Fd2u/dzv3rKpv9LPqKbZ1dlXtQxeEP0f3T9h0tmlBZsL2dMedqjqTrurR41ueZ3tBXUX3z93k/TsryR/R1e19Nl11gM3o6sj25mmm74A0weu513Ov56O/ng9iXT6DyecLrNv34U7a9/WPub0A+I90x/zBVbUJXcGh9/h9DnhIkgfR3ak4drbbHiYLFQtUVa2kq3v5Vz1p19CdtH+WZL3WZd59BtzUXkkel2R9ukZMZ1bVKrp/1u6XZL8kd2+vP0zyB33mfxVdQ6d3J7lHkofQNej79z7z9X66H3GfTHKfdDamq8s74ft0/1Q9Pcnd6eqHbtDn+qdzNN3F+NnM3EvIccBrk+zYLgb/SFcX8rbemarqSroL/PuSbJLkbm1/nri2jLRlvwT8S5LN2mfwhDb534CXJ3lUOzb3asdh4+nXOKNP0QXf57ThCYcDb07yQIAkmyZ5bp/r/AmwbTu3SLJ+un7EN62qX9Pdvv/tDMv/HvBXbb+fS1cP+OSe6cfQNfz7dVXNtqvDw4G3TJzXSZYmeU6btjHdD4tr6f6dfBvdnQppnXg993ru9XxOrueDWJfP4GS679MLkyxJ8ny66nhfWNeNJtkgySPoCgnX0zXkhy7+/By4Mck2dI2zf6eqbqYrNH+KrirWj9Z126NgoWJhezt3/hHz53Qn10/p6h728w/DTD5F9y/adcAjaLfT2m3upwL70pXKrwLew7pd5F9AVx/1x3SN2w6pqv/pZ8Gqupaui7Wb6XrTuImuV4aN6Rq/UVU30vXs8DG64PwLup56BnEGXf3n1VV19gzzHUX3b8oZdNVqbgZePc28L6a7xX8J3UXjRHqqAqzFfnR1e79H1xDrrwGq6hy6c+EjbZ0r6eqTzuSC3LFf8w/2TDsJ2Bm4qqoumEisqv+k+9yPT3cL9iJgzz7zfhpddaOrklzbsz9XtHW9nK4KxnTOanm6FngX8Jyq+mnP9E/S1Vfv94fNnVTVZ+h+8Hym5elC4Glt8snA/9D1HHIFXdC8crbb0qLn9dzrudfzEV7PB7Eun0HL997A6+m+u28E9m7neb/emOSmtvwxdO09HlO3d8v8D3TdLt9IV23vs1Os42i6tlYLouoT3N4jgSSNlST3pAvMu1bVZfOdH0nS7Hg9X3dJtqcroP5+Vf1svvMD3qmQNL5eAZxtAJKksef1fB2ka5/0Oro2MwuiQAFd7waSNFaSXEHXaO2Z85wVSdIAvJ6vm3RP3f4JXY9Te8xzdu7A6k+SJEmSBmL1J0mSJEkDsVAhSZIkaSB3yTYVW265ZS1fvny+syFJC9a55557bVXN+uFidxXGC0maWb/x4i5ZqFi+fDnnnHPOfGdDkhasJD+c7zwsBMYLSZpZv/HC6k+SJEmSBmKhQpIkSdJALFRIkiRJGoiFCkmSJEkDsVAhSZIkaSAWKiRJkiQNxEKFJEmSpIFYqJAkSZI0EAsVkiRJkgZioUKSJEnSQCxUSJIkSRrIkvnOgCSNu+UHf3Gdl7ni0KePICeSpNmazbV8XMxFzPFOhSRJkqSBWKiQJEmSNBALFZIkSZIGYqFCkjSvkhyV5OokF/WkbZ7klCSXtffNWnqSfDjJyiQXJtm1Z5n92/yXJdl/PvZFkhYrCxWSpPn2CWCPSWkHA6dW1c7AqW0cYE9g5/Y6CDgMukIIcAjwKOCRwCETBRFJ0uhZqJAkzauqOgO4blLyPsDRbfho4Jk96cdU50xgaZKtgacBp1TVdVV1PXAKdy6oSJJGxEKFJGkh2qqqrmzDVwFbteFtgFU9861uadOlS5LmgIUKSdKCVlUF1LDWl+SgJOckOeeaa64Z1molaVEbWaFiqoZ3Lf3VSb6X5OIk/9ST/ubW8O7SJE/rSd+jpa1McjCSpMXgJ61aE+396pa+BtiuZ75tW9p06XdSVUdU1YqqWrFs2bKhZ1ySFqNR3qn4BJPqsyb5I7r6sA+tqgcC/9zSdwH2BR7YlvmXJOslWQ/4KF3DvF2AF7R5JUl3bScBEz047Q98vif9xa0XqN2AG1s1qS8DT02yWWug/dSWJkmaA0tGteKqOiPJ8knJrwAOrapb2jwT/zztAxzf0n+QZCVd7x0AK6vqcoAkx7d5LxlVviVJcyvJccCTgC2TrKbrxelQ4IQkBwI/BJ7XZj8Z2AtYCfwSeAlAVV2X5B3A2W2+t1fV5MbfkqQRGVmhYhr3Ax6f5F3AzcAbqupsusZ0Z/bM19vAbnLDu0fNRUYlSXOjql4wzaTdp5i3gFdOs56jgKOGmDVJUp/mulCxBNgc2A34Q7p/oXYaxoqTHETXZznbb7/9MFYpSZIkqQ9z3fvTauCzrX/xbwG/BbbEhneSJEnS2JrrQsXngD8CSHI/YH3gWrqGd/sm2SDJjnRPSv0WXd3YnZPsmGR9usbcJ81xniVJkiTNYGTVn6ZpeHcUcFTrZvZWYP9WP/biJCfQNcC+DXhlVf2mredVdD14rAccVVUXjyrPkiRJktbdKHt/mq7h3Z9NM/+7gHdNkX4yXW8fkiRJkhYgn6gtSZIkaSAWKiRJkiQNxEKFJEmSpIFYqJAkSZI0EAsVkiRJkgZioUKSJEnSQCxUSJIkSRqIhQpJkiRJA7FQIUmSJGkgFiokSZIkDcRChSRJkqSBWKiQJEmSNBALFZIkSZIGYqFCkiRJ0kAsVEiSJEkaiIUKSZIkSQOxUCFJkiRpIBYqJEmSJA3EQoUkSZKkgViokCRJkjQQCxWSJEmSBmKhQpIkSdJARlaoSHJUkquTXDTFtNcnqSRbtvEk+XCSlUkuTLJrz7z7J7msvfYfVX4lSZIkzc4o71R8AthjcmKS7YCnAj/qSd4T2Lm9DgIOa/NuDhwCPAp4JHBIks1GmGdJkiRJ62hkhYqqOgO4bopJHwDeCFRP2j7AMdU5E1iaZGvgacApVXVdVV0PnMIUBRVJkiRJ82dO21Qk2QdYU1UXTJq0DbCqZ3x1S5suXZIkSdICsWSuNpRkQ+Bv6ao+jWL9B9FVnWL77bcfxSYkSZIkTWGtdyqSPDbJvdrwnyV5f5IdZrGt+wA7AhckuQLYFjgvye8Da4DteubdtqVNl34nVXVEVa2oqhXLli2bRfYkSYMYYryQJI2Zfqo/HQb8MslDgdcD/wccs64bqqrvVNXvVdXyqlpOV5Vp16q6CjgJeHHrBWo34MaquhL4MvDUJJu1BtpPbWmSpIVnKPFCkjR++ilU3FZVRdeY+iNV9VFg47UtlOQ44JvA/ZOsTnLgDLOfDFwOrAT+DfhLgKq6DngHcHZ7vb2lSZIWnlnFC0nS+OunTcVNSd4M/BnwhCR3A+6+toWq6gVrmb68Z7iAV04z31HAUX3kU5I0v2YVLyRJ46+fOxXPB24BDmxVlbYF3jvSXEmSxpHxQpIWqX7uVLy2qt40MVJVP0rywBHmSZI0nowXkrRI9XOn4ilTpO057IxIksae8UKSFqlp71QkeQVdg+mdklzYM2lj4BujzpgkaTwYLyRJM1V/+hTwJeDdwME96TfZA5MkqYfxQpIWuWmrP1XVjVV1RevFaTXwa6CAjZL4yGpJEjDaeJHktUkuTnJRkuOS3CPJjknOSrIyyaeTrN/m3aCNr2zTlw+6b5Kk/vTzRO1XAT8BTgG+2F5fGHG+JEljZtjxIsk2wF8BK6rqQcB6wL7Ae4APVNV9geuBiecgHQhc39I/0OaTJM2Bfnp/+mvg/lX101FnRpI01kYRL5YA90zya2BD4ErgycAL2/SjgbfRPc17nzYMcCLwkSRpz0KSJI1QP70/rQJuHHVGJEljb6jxoqrWAP8M/IiuMHEjcC5wQ1Xd1mZbDWzThrdpeaBNvxHYYlj5kSRNr587FZcDpyf5It1DjQCoqvePLFeSpHE01HiRZDO6uw87AjcAnwH2GDSTSQ4CDgLYfnubCErSMPRTqPhRe63fXpIkTWXY8eKPgR9U1TUAST4LPBZYmmRJuxuxLbCmzb8G2A5YnWQJsClwp6pYVXUEcATAihUrrBolSUOw1kJFVf0DQJINq+qXo8+SJGkcjSBe/AjYLcmGwK+A3YFzgK8CzwGOB/YHPt/mP6mNf7NNP832FJI0N/rp/enRSS4BvtfGH5rkX0aeM0nSWBl2vKiqs+gaXJ8HfIcuZh0BvAl4XZKVdG0mjmyLHAls0dJfxx2fmSFJGqF+qj99EHga3T9AVNUFSZ4w0lxJksbR0ONFVR0CHDIp+XLgkVPMezPw3EG2J0manX56f6KqVk1K+s0I8iJJGnPGC0lanPq5U7EqyWOASnJ34DXAd0ebLUnSGDJeSNIi1c+dipcDr6Tr/3sN8LA2LklSL+OFJC1S/dypSFW9aOQ5kSSNO+OFJC1S/dyp+H9JvpLkwCRLR54jSdK4Ml5I0iK11kJFVd0PeCvwQOC8JF9I8mcjz5kkaawYLyRp8eq396dvVdXr6Lrwuw44eqS5kiSNJeOFJC1O/Tz8bpMk+yf5EvAN4Eqm6B9ckrS4GS8kafHq507FBXQ9eLy9qu5XVW+qqnPXtlCSo5JcneSinrT3JvlekguT/Gdvndskb06yMsmlSZ7Wk75HS1uZxKejStLCNat4IUkaf/0UKnaqqtdW1TcnEpL002vUJ4A9JqWdAjyoqh4CfB94c1vfLsC+dPVw9wD+Jcl6SdYDPgrsCewCvKDNK0laeGYbLyRJY27aQkWSrwNUVSX55KTJ31rbiqvqDLr6tL1pX6mq29romcC2bXgf4PiquqWqfgCspLtl/khgZVVdXlW3Ase3eSVJC8Sg8UKSNP5mulNxr57hB02aliFs+6XAl9rwNsCqnmmrW9p06XeS5KAk5yQ555prrhlC9iRJfRp1vJAkLXAzFSpqmuGpxtdJkrcAtwHHDrKeXlV1RFWtqKoVy5YtG9ZqJUlrN7J4IUkaDzPVdV2a5E/pCh5LkzyrpQfYdLYbTHIAsDewe1VNBJs1wHY9s23b0pghXZK0MIwkXkiSxsdMhYqvAc/oGf6TnmlnzGZjSfYA3gg8sap+2TPpJOBTSd4P3BvYma4eboCdk+xIV5jYF3jhbLYtSRqZoccLSdJ4mbZQUVUvGWTFSY4DngRsmWQ1cAhdb08bAKckATizql5eVRcnOQG4hK5a1Cur6jdtPa8CvgysBxxVVRcPki9J0nANGi8kSeNvZF39VdULpkg+cob53wW8a4r0k4GTh5g1SZIkSUPUz3MqJEmSJGlaFiokSZIkDWSthYokGyb5uyT/1sZ3TrL36LMmSRonxgtJWrz6uVPxceAW4NFtfA3wzpHlSJI0rowXkrRI9VOouE9V/RPwa4DWFaxPSJUkTWa8kKRFqp9Cxa1J7kl7KmqS+9D9EyVJUi/jhSQtUv10KXsI8N/AdkmOBR4LHDDKTEmSxpLxQpIWqbUWKqrqlCTnAbvR3cZ+TVVdO/KcSZLGivFCkhavfnp/+lPgtqr6YlV9AbgtyTNHnzVJ0jgxXkjS4tVPm4pDqurGiZGquoHuFrckSb2MF5K0SPVTqJhqnn7aYkiSFhfjhSQtUv0UKs5J8v4k92mv9wPnjjpjkqSxY7yQpEWqn0LFq4FbgU+31y3AK0eZKUnSWDJeSNIi1U/vT78ADp6DvEiSxpjxQpIWr7UWKpLcD3gDsLx3/qp68uiyJUkaN8YLSVq8+mlA9xngcOBjwG9Gmx1J0hgzXkjSItVPoeK2qjps5DmRJI0744UkLVL9NNT+ryR/mWTrJJtPvEaeM0nSuDFeSNIi1c+div3b+9/0pBWw0/CzI0kaY8YLSVqk+un9ace5yIgkabwZLyRp8Vpr9ackGyZ5a5Ij2vjOSfYefdYkSePEeCFJi1c/bSo+Tvcwo8e08TXAO9e2UJKjklyd5KKetM2TnJLksva+WUtPkg8nWZnkwiS79iyzf5v/siT7T7UtSdKCMKt4IUkaf/20qbhPVT0/yQsAquqXSdLHcp8APgIc05N2MHBqVR2a5OA2/iZgT2Dn9noUcBjwqNbA7xBgBV293HOTnFRV1/e1d5K0FssP/uKslrvi0KcPOSd3CbONF5KkMdfPnYpbk9yT7kc9Se4D3LK2harqDOC6Scn7AEe34aOBZ/akH1OdM4GlSbYGngacUlXXtYLEKcAefeRZkjT3ZhUvJEnjr587FYcA/w1sl+RY4LHAAbPc3lZVdWUbvgrYqg1vA6zqmW91S5suXZK08AwzXkiSxkg/vT+dkuQ8YDcgwGuq6tpBN1xVlaQGXc+EJAcBBwFsv/32w1qtJKlPo4gXSZbSPaH7QXR3QF4KXAp8GlgOXAE8r6qub1WtPgTsBfwSOKCqzhtk+5Kk/kxb/SnJrhMvYAfgSuDHwPa9DanX0U9atSba+9UtfQ2wXc9827a06dLvpKqOqKoVVbVi2bJls8yeJGldjSheTPgQ8N9V9QDgocB3ub193s7AqW0c7tg+7yC69nmSpDkw052K980wrYAnz2J7J9E9HOnQ9v75nvRXJTmerqH2jVV1ZZIvA/840UsU8FTgzbPYriRpdEYRL0iyKfAEWhWqqrqVrt3GPsCT2mxHA6fTdfrxu/Z5wJlJlibZuqfarSRpRKYtVFTVHw2y4iTH0V30t0yymq6u7aHACUkOBH4IPK/NfjLd7eqVdLesX9LycF2SdwBnt/neXlWTG39LkubRoPFiBjsC1wAfT/JQ4FzgNax7+7w7FCqsLitJw7fWNhVJ7g68gu7fIuj+EfrXqvr1TMtV1QummbT7FPMW8Mpp1nMUcNTa8ilJml+zjRczWALsCry6qs5K8iFur+oEzK59XlUdARwBsGLFiqG17ZOkxayfLmUPAx4B/Et7PQLrqUqS7mzY8WI1sLqqzmrjJ9IVMta1fZ4kacT66VL2D6vqoT3jpyW5YFQZkiSNraHGi6q6KsmqJPevqkvp7nRf0l59t8+b7fYlSf3rp1DxmyT3qar/A0iyE/Cb0WZLkjSGRhEvXg0cm2R94HK6Nnd3Yx3a50mSRq+fQsXfAF9Ncjldv+M74IVaknRnQ48XVXU+sGKKSevUPk+6q1t+8BfnOwsjccWhT5/vLKhP/Tz87tQkOwP3b0mXVtUto82WJGncGC8kafGatlCR5MlVdVqSZ02adN8kVNVnR5w3SdIYMF5Ikma6U/FE4DTgT6aYVoBBQpIExgtJWvRmKlQcBlBVtp+QJM1kUcaLu2oddrAeu6R1N9NzKs5P8j9JDkyydM5yJEkaN8YLSVrkZipUbAO8F3gccGmSzyfZN8k95yZrkqQxYbyQpEVu2kJFVf2mqr7cbmdvBxwF7AP8IMmxc5VBSdLCZryQJM10p+J3qupWuieYfhf4GfAHo8yUJGk8GS8kaXGasVCRZLskf5PkPOALbf5nVNWuc5I7SdJYMF5I0uI203MqvkFXT/YE4M+r6tw5y5UkaWwYLyRJM3UpezDwv1VVc5UZSdJYMl5I0iI3baGiqs6Yy4xIksaT8UKS1FdDbUmSJEmajoUKSZIkSQOZqU0FAEk2AJ4NLO+dv6rePrpsSZLGjfFCkhavtRYqgM8DNwLnAreMNjuSpDFmvJCkRaqfQsW2VbXHyHMiSRp3xgtJWqT6aVPxjSQPHnlOJEnjznghSYtUP4WKxwHnJrk0yYVJvpPkwkE2muS1SS5OclGS45LcI8mOSc5KsjLJp5Os3+bdoI2vbNOXD7JtSdLIDD1eSJLGQz/Vn/Yc5gaTbAP8FbBLVf0qyQnAvsBewAeq6vgkhwMHAoe19+ur6r5J9gXeAzx/mHmSJA3FUOOFJGl89HOnoqZ5DWIJcM8kS4ANgSuBJwMntulHA89sw/u0cdr03ZNkwO1LkoZvFPFCkjQG+rlT8UW6oBDgHsCOwKUcRGF7AAAgAElEQVTAA2ezwapak+SfgR8BvwK+QtdTyA1VdVubbTWwTRveBljVlr0tyY3AFsC1s9m+JGlkhhovJEnjY62Fiqq6Q6O7JLsCfznbDSbZjO7uw47ADcBngIF7C0lyEHAQwPbbbz/o6iRJ62jY8UKSND7W+YnaVXUe8KgBtvnHwA+q6pqq+jXwWeCxwNJWHQpgW2BNG14DbAfQpm8K/HSKfB1RVSuqasWyZcsGyJ4kaRiGEC8kSWOinydqv65n9G7ArsCPB9jmj4DdkmxIV/1pd+Ac4KvAc4Djgf3pHqIEcFIb/2abflpVWUdXkhaYEcQLSdKY6KdNxcY9w7fR1Zn9j9lusKrOSnIicF5b37eBI9p6j0/yzpZ2ZFvkSOCTSVYC19H1FCVJWniGGi8kSeOjnzYV/zDsjVbVIcAhk5IvBx45xbw3A88ddh4kScM1inghSRoP07apSPK4JC/uGT8xyWnt9eS5yZ4kaaEzXkiSZrpT8Q/Aq3vG7w8cANwL+FvgtNFlS5I0RowXkrTIzdT70yZVdUnP+GVVdW5VncEd681KkhY344UkLXIzFSqW9o5U1bN6RrcaTXYkSWPIeCFJi9xMhYrvJXn65MQke9M9IVWSJDBeSNKiN1ObitcCX0zyHLruXwEeATwG2HvUGZMkjQ3jhSQtctPeqaiqlcBDgP8FlrfXGcBDqur7c5E5SdLCZ7yQJM34nIqqugU4ao7yIkkaU8YLSVrcZmpTIUmSJElrZaFCkiRJ0kBmeqL2qe39PXOXHUnSuDFeSJJmalOxdZLHAM9IcjyQ3olVdd7Ui0mSFhnjhSQtcjMVKv4e+DtgW+D9k6YV8ORRZUqSNFaMF5K0yE1bqKiqE4ETk/xdVb1jDvMkSRojxgtJ0oxdygJU1TuSPAN4Qks6vaq+MNpsSZLGzajiRZL1gHOANVW1d5IdgeOBLYBzgf2q6tYkGwDH0D1476fA86vqikG3L0lau7X2/pTk3cBrgEva6zVJ/nHUGZMkjZcRxovXAN/tGX8P8IGqui9wPXBgSz8QuL6lf6DNJ0maA/10Kft04ClVdVRVHQXsAew92mxJksbQ0ONFkm3bej/WxkPXRuPENsvRwDPb8D5tnDZ99za/JGnE1lr9qVkKXNeGNx1RXiRJ42/Y8eKDwBuBjdv4FsANVXVbG18NbNOGtwFWAVTVbUlubPNf27vCJAcBBwFsv/32Q8iiBrH84C/OdxZG4opDnz7fWZDmVD+FincD307yVbpuAp8AHDzSXEmSxtFQ40WSvYGrq+rcJE8aThahqo4AjgBYsWJFDWu9krSY9dNQ+7gkpwN/2JLeVFVXjTRXkqSxM4J48Vi6Z1/sBdwD2AT4ELA0yZJ2t2JbYE2bfw2wHbA6yRK6OyU/HWD7kqQ+9dOmgqq6sqpOai8LFJKkKQ0zXlTVm6tq26paDuwLnFZVLwK+CjynzbY/8Pk2fFIbp00/raq8EyFJc6CvQoUkSQvIm4DXJVlJ12biyJZ+JLBFS38dVtWVpDnTb0PtoUqylK4njwfRPW31pcClwKeB5cAVwPOq6vrWc8eHgL2AXwIHVNV585BtSdI8qarTgdPb8OXAI6eY52bguXOaMUkSsJY7FUnWS/K9EWz3Q8B/V9UDgIfS9T9+MHBqVe0MnMrt/zDtCezcXgcBh40gP5KkAYwwXkiSxsCMhYqq+g1waZKh9bmXZFO6HkGObNu4tapu4I79i0/ud/yY6pxJ10Bv62HlR5I0uFHEC0nS+Oin+tNmwMVJvgX8YiKxqp4xy23uCFwDfDzJQ4Fz6Z6WulVVXdnmuQrYqg3/rt/xZqJP8iuRJC0kw44XkqQx0U+h4u9GsM1dgVdX1VlJPsSkxnRVVUnWqccOH2YkSfNu2PFCkjQm1tr7U1V9ja7h9N3b8NnAIA2lVwOrq+qsNn4iXSHjJxPVmtr71W36RL/jE3r7JO/N5xFVtaKqVixbtmyA7EmSZmME8UKSNCbWWqhI8ud0P/z/tSVtA3xuthts/ZavSnL/lrQ7cAl37F98cr/jL05nN+DGnmpSkqQFYtjxQpI0Pvqp/vRKuq77zgKoqsuS/N6A2301cGyS9YHLgZfQFXBOSHIg8EPgeW3ek+m6k11J16XsSwbctiRpNEYRLyRJY6CfQsUtVXVr97gISLKE7tkSs1ZV5wMrppi0+xTzFl2gkiQtbEOPF5Kk8dDPE7W/luRvgXsmeQrwGeC/RpstSdIYMl5I0iLVT6HiYLouYL8D/AVddaS3jjJTkqSxZLyQpEVqrdWfquq3SY6mqyNbwKWtSpIkSb9jvJCkxWuthYokTwcOB/4PCLBjkr+oqi+NOnOSpPFhvJCkxaufhtrvA/6oqlYCJLkP8EXAICFJ6mW8kKRFqp82FTdNBIjmcuCmEeVHkjS+jBeStEhNe6ciybPa4DlJTgZOoKsj+1y6p6RKkmS8kCTNWP3pT3qGfwI8sQ1fA9xzZDmSJI0b44UkLXLTFiqqyidXS5LWynghSeqn96cdgVcDy3vnr6pnjC5bkqRxY7yQpMWrn96fPgccSfdU1N+ONjuSpDFmvJCkRaqfQsXNVfXhkedEkjTujBeStEj1U6j4UJJDgK8At0wkVtV5I8uVJGkcGS8kaZHqp1DxYGA/4Mncfju72rgkSROMF5K0SPVTqHgusFNV3TrqzEiSxprxQpIWqX6eqH0RsHTUGZEkjT3jhSQtUv3cqVgKfC/J2dyxjqxdBEqSehkvJGmR6qdQccjIcyFJuiswXkjSIrXWQkVVfW0uMiJJGm/GC0lavPp5ovZNdL13AKwP3B34RVVtMsqMSZLGi/FCkhavfu5UbDwxnCTAPsBuo8yUJGn8GC8kafHqp/en36nO54CnjSg/kqS7AOOFJC0u/VR/elbP6N2AFcDNg244yXrAOcCaqto7yY7A8cAWwLnAflV1a5INgGOARwA/BZ5fVVcMun1J0nCNKl5Ikha+fnp/+pOe4duAK+huaQ/qNcB3gYm6tu8BPlBVxyc5HDgQOKy9X19V902yb5vv+UPYviRpuEYVLyRJC1w/bSpeMuyNJtkWeDrwLuB1re7tk4EXtlmOBt5GV6jYpw0DnAh8JEmqqpAkLRijiBeSpPEwbaEiyd/PsFxV1TsG2O4HgTcCE436tgBuqKrb2vhqYJs2vA2wqm30tiQ3tvmvnZTfg4CDALbffvsBsiZJWhcjjheSpDEwU0PtX0zxgq460ptmu8EkewNXV9W5s13HVKrqiKpaUVUrli1bNsxVS5JmNpJ4IUkaH9Peqaiq900MJ9mYrg3ES+gaU79vuuX68FjgGUn2Au5B16biQ8DSJEva3YptgTVt/jXAdsDqJEuATekabEuSFoARxgtJ0piYsUvZJJsneSdwIV0BZNeqelNVXT3bDVbVm6tq26paDuwLnFZVLwK+CjynzbY/8Pk2fFIbp00/zfYUkrSwjCJeSJLGx7SFiiTvBc4GbgIeXFVvq6rrR5iXN9E12l5J12biyJZ+JLBFS38dcPAI8yBJWkfzEC8kSQvMTL0/vR64BXgr8JaugyYAQtfwbpPpFuxXVZ0OnN6GLwceOcU8NwPPHXRbkqSRGXm8kCQtbDO1qVinp21LkhYn44UkyUAgSZIkaSAWKiRJC1KS7ZJ8NcklSS5O8pqWvnmSU5Jc1t43a+lJ8uEkK5NcmGTX+d0DSVo8LFRIkhaq24DXV9UuwG7AK5PsQtdhx6lVtTNwKrd34LEnsHN7HQQcNvdZlqTFyUKFJGlBqqorq+q8NnwT8F1gG2Af4Og229HAM9vwPsAx1TmT7vlHW89xtiVpUbJQIUla8JIsBx4OnAVsVVVXtklXAVu14W2AVT2LrW5pkqQRs1AhSVrQkmwE/Afw11X1s95p7WGo6/RA1CQHJTknyTnXXHPNEHMqSYuXhQpJ0oKV5O50BYpjq+qzLfknE9Wa2vvEU7vXANv1LL5tS7uDqjqiqlZU1Yply5aNLvOStIhYqJAkLUjpnqJ3JPDdqnp/z6STgP3b8P7A53vSX9x6gdoNuLGnmpQkaYRmeqK2JEnz6bHAfsB3kpzf0v4WOBQ4IcmBwA+B57VpJwN7ASuBXwIvmdvsStLiZaFCkrQgVdXXgUwzefcp5i/glSPNlCRpSlZ/kiRJkjQQCxWSJEmSBmKhQpIkSdJALFRIkiRJGoiFCkmSJEkDsVAhSZIkaSAWKiRJkiQNxEKFJEmSpIFYqJAkSZI0EAsVkiRJkgYy54WKJNsl+WqSS5JcnOQ1LX3zJKckuay9b9bSk+TDSVYmuTDJrnOdZ0mSJEnTm487FbcBr6+qXYDdgFcm2QU4GDi1qnYGTm3jAHsCO7fXQcBhc59lSZIkSdOZ80JFVV1ZVee14ZuA7wLbAPsAR7fZjgae2Yb3AY6pzpnA0iRbz3G2JUmSJE1jXttUJFkOPBw4C9iqqq5sk64CtmrD2wCrehZb3dIkSZIkLQDzVqhIshHwH8BfV9XPeqdVVQG1jus7KMk5Sc655pprhphTSZIkSTOZl0JFkrvTFSiOrarPtuSfTFRrau9Xt/Q1wHY9i2/b0u6gqo6oqhVVtWLZsmWjy7wkSZKkO5iP3p8CHAl8t6re3zPpJGD/Nrw/8Pme9Be3XqB2A27sqSYlSZIkaZ4tmYdtPhbYD/hOkvNb2t8ChwInJDkQ+CHwvDbtZGAvYCXwS+Alc5tdSQvV8oO/OKvlrjj06UPOiSRJi9ucFyqq6utAppm8+xTzF/DKkWZKkiRJ0qz5RG1JkiRJA7FQIUmSJGkg89GmQpJsDyFJ0l2IdyokSZIkDcRChSRJkqSBWKiQJEmSNBDbVEjSAjCbNia2L5EkLRTeqZAkSZI0EAsVkiRJkgZioUKSJEnSQCxUSJIkSRqIhQpJkiRJA7FQIUmSJGkgFiokSZIkDcRChSRJkqSBWKiQJEmSNBALFZIkSZIGYqFCkiRJ0kAsVEiSJEkayJL5zoBGZ/nBX1znZa449OkjyInuijy/JEnSBO9USJIkSRqIdyommc2/r+A/sBov3mWQJEnDNDaFiiR7AB8C1gM+VlWHznOWpjWMgslCKdwshB+f43wsRpEPSdMbp1ghSXclY1GoSLIe8FHgKcBq4OwkJ1XVJfObM/XDgslwLYTjKS1ExgpJmj9jUagAHgmsrKrLAZIcD+wDGCg0ViwQSCNlrJCkeTIuDbW3AVb1jK9uaZIkTTBWSNI8SVXNdx7WKslzgD2q6mVtfD/gUVX1qp55DgIOaqP3By4dQVa2BK51HUNbx0LIg+tYeHlwHcPPw1R2qKplI1jvvOknVrT0uYgXwzaq82Cx8TgOh8dxOMblOPYVL8al+tMaYLue8W1b2u9U1RHAEaPMRJJzqmqF6xjOOhZCHlzHwsuD6xh+HhaRtcYKmJt4MWyeB8PhcRwOj+Nw3NWO47hUfzob2DnJjknWB/YFTprnPEmSFhZjhSTNk7G4U1FVtyV5FfBlum4Cj6qqi+c5W5KkBcRYIUnzZywKFQBVdTJw8jxnYxi3y13HwsqD61h4eXAdw8/DorFAYsUoeB4Mh8dxODyOw3GXOo5j0VBbkiRJ0sI1Lm0qJEmSJC1Qi7pQkWS7JF9NckmSi5O8pqVvnuSUJJe1981aepJ8OMnKJBcm2bVnXfu3+S9Lsv+47kuSHZKcl+T8tp6Xj+N+tGnbJ/lKku+29S0f4315T5KL2uv5C3w/HpDkm0luSfKGta1nHPelTbsiyXfad+Wcud4XLUxJ9khyafseHzzF9A2SfLpNP6v3upTkzS390iRPm2b9T0pyYzvvzk/y99PM94kkP+iZ72EtfbMk/9muMd9K8qDh7PngkixNcmKS77Xr9qMnTZ/pGrnWGJxki/b9/3mSj0ya9q4kq5L8fIb8LU/yq55jenjPtOe3PF2c5D2zPwrrLslRSa5OclFP2nvbcbywfd5Lp1l2yvM1XWcHZ7X0T6fr+GCq5ac8bkk+0HOcvp/khmmWP71tf2Le32vpOyQ5teX/9CTbzubY9GuaY/iwJGdOXOOTPHKaZac895I8osWIle28zRTLThsXk7yj7f/56X7L3Hua7f+m5/id1JP+5HS/5y5KcnSS0Td5qKpF+wK2BnZtwxsD3wd2Af4JOLilHwy8pw3vBXwJCLAbcFZL3xy4vL1v1oY3G9N9WR/YoA1vBFwB3Hvc9qNNOx14Ss++bDimn8nTgVPo2kDdi66Hm00W8H78HvCHwLuAN6xtPQv8M5lyX9q0K4At5zL/vhb2i65x+P8BO7Vr6QWTz3HgL4HD2/C+wKfb8C5t/g2AHdt61ptiG08CvtBHXj4BPGeK9PcCh7ThBwCnzvdx68nb0cDL2vD6wNJJ0weKwe36+Tjg5cBHJk3brV0ffj5D/pYDF02RvgXwI2BZz37sPofH7QnArr15A54KLGnD75m4pvV7vgInAPu24cOBV0yz7X6O26vpOk2YatrpwIop0j8D7N+Gnwx8ch6O4VeAPXvOvdOnWG7acw/4Vjs+aeftnlMsP21cpCfOA39Fu25MsY47HXu6mwargPu18bcDB476XFzUdyqq6sqqOq8N3wR8l+7pq/vQXRRo789sw/sAx1TnTGBpkq2BpwGnVNV1VXU93Q/APeZwV4a2L1V1a1Xd0ubZgDm+mzWs/UiyC90F9ZS2rp9X1S/HcV/ofmycUVW3VdUvgAuZw/NrXfejqq6uqrOBX/e5njkzrH2RpvFIYGVVXV5VtwLH051bvXrPtROB3ds/mPsAx1fVLVX1A2BlW9+w7QKcBlBV3wOWJ9lqBNtZJ0k2pfthdyRAi0WT/90eKAZX1S+q6uvAzVNMO7Oqrpxl9ncCLquqa9r4/wDPnuW61llVnQFcNyntK1V1Wxs9k+6ZLZNNeb628/HJdOcn3DFOTd52P8ftBcBxfe3M7X53ngJf5c7fo6Ga6hgCBWzShjcFfjzFolOee+283KQdnwKOYYpjOFNcrKqf9cx6r5affm0B3FpV32/jpzAH5+SiLlT0SncL+uHAWcBWPV+Sq4CJC+42dCW/Catb2nTp82LAfZm4HXdhm/6eqprqizRyA+7H/YAbknw2ybfbreD15iTjUxhwXy6gu0htmGRL4I+44wO+5kyf+7Gu65kXQ9iXAr6S5Nx0T2iWpvwOJ3l7kmdMnqf96LuR7gfAusSRRye5IMmXkjxwIjHJyZOqSLyrVZ/4QJINWtoFwLPa/I8EdmDqH5xzbUfgGuDj7Zr9sST3SvLy3F4Nd85jcJJnJHl7bz5b/r6W5PEtbSVw/3TVo5bQ/Xicl2v0NF5K9085Se6dZKJ3tOmO2xbADT2FklkfzyQ70H22p/WknT9pto+3qjt/11NF6HfnKfCnwMZJtphNHgbw18B7k6wC/hl4M0CSFUk+1uaZ6ZxcPUX6tKaKixPVy4AXAX8/xfYB7tGqZ52ZZKLgci2wJMnEg/WewxyckxYqgCQbAf8B/PWkkiGthDk2XWQNY1+qalVVPQS4L7D/fPyLNYT9WAI8HngDXfWVnYADhp/TtRt0X6rqK3RdZH6D7t+ebwK/GU1upzes78lM65krQ9qXx1XVrsCewCuTPGH4OdVdQVX9fVUN6yF85wE7VNVDgf8P+FzPdvbq+RPozXTVm/6QrmrGm1r6oXT/8J9PVy3l28zD9WQKS+iqnxxWVQ8HfkFXJfHwqjp85kVHp6pOqqqJditXAtu3/L0O+FSSTdo/1K8APg38L13VyIVwTEnyFuA24FiAqvpxVe01h1nYFzixqn53PKrqYT3TX1RVD6aL148H9mvpbwCemOTbwBOBNcz9MX0F8Nqq2g54LbffRTunql42zA1NF5Oq6i1t+8cCr5pm+ztU91TuFwIfTHKfFsf2BT6Q5FvATczB8Vv0hYokd6f7II+tqs+25J+0W1e096tb+hruWNLbtqVNlz6nhrQvv9OC00V0X/Q5M6T9WA2c327r3kYXeHdljg3rM6mqd1XVw6rqKXT1M7/PHFrH/VjX9cypYe1LVU18NlcD/8loqqpovPQTC343T/tXe1Pgp9Mtm+SVPY0w711VP6uqn8Pvnslx93YH8w5atYqqrjrrx2nnZ1v+Je2H3YuBZXT1wOfbamB1VU38S3sid75mr1MMTvKnPcduBQNqVdN+2obPpWuPcL82/l9V9aiqejRwKXN8jZ5KkgOAvel+uE/1R8l0x/OndAXPJb3pSdbrOZ5vpz/7MkPVp57r6E3Ap7j9PP1xVT2rFeDe0tKmbOw9QvsDEzHiM0x9jZ/pnNx2cnqrCTJxDF8OfcfFY5mm+lLPMbycro3Kw9v4N6vq8VX1SOAM5uCcXNSFinab7Ujgu1X1/p5JJ9GdTLT3z/ekvzid3YAbW5WJLwNPTderxmZ0DaS+PCc70QxrX5Jsm+SebZ2b0TVqu3ROdoKhfiZn010Ul7X5ngxcMvId6DHEz2S9idu+SR4CPISuAdmcmMV+rOt65swQ9+VeSTaeGKb7zl800zJaFM4Gdk7Xc876dD+oJt+h6D3XngOc1n7wnQTsm653qB2BnYFvVdVH2x8KD6uqHyf5/YkqIumqL92N7kfgHfQUkkNXHeeiNr40t/fk8zK69lrzcsewV1VdBaxKcv+WtDt3vmavUwyuqv/sOXYD99CWZFlaNdokO9F9Rpe38Yleizaja4z/senWMxeS7AG8EXhGTd+ecMrztZ2PX6U7P6FdE6vqNz3Hc8pexybl4QF0jZe/Oc30JRMF4vbDem9uP0+3TDLxG/XNwFFr3+uh+zHdXRLofkNcNsU80517VwI/S7Jb+w6+mO4Yruo5hofPFBeT7Nwzug/wvckbb9vdoA1vCTyW9r3pOSc3oLtTOfo7fjXiluAL+UX3g7noGr6e31570dUnPJXuBPofYPM2f4CP0v078R16eiygq7O4sr1eMq77AjylreOC9n7QOO7HpH35Dl1PKOuP474A96C7SFxC1+DuYQt8P36f7l/HnwE3tOFNplvPmO7LTu07cgFwMfCWudwPXwv31c6n77fv8Vta2tvpftxNfJ8/QxcrvgXs1LPsW9pylzJFTzFtnle1c+6Cdj14TM+0k2m99dHVYf8O3Y+0fwc2aumPbvm7lO5f2DntqXAtx+5hwDnt+/k5uh+kLwde3qYPHIPpqiZdB/y8fZ8netr5pzb+2/b+tpb+DODtbfjZ7difT1cN7U961nsct1+n953j43YcXdWsX7e8H9iOw6qe69xEj2P3Bk6e6Xxt6Tu183NlO183mGbbUx63Nu1twKFTLHN+e78XcG77vC8GPkTr8YyuQHNZy9vHptv+iI/h41r+LqBr5/CINu8K4GNrO/fafBe1Y/sR6B44PWm708ZFursXF7Vp/wVsM3n7wGPad+GC9n5gz7rfS9fw+1K6alUjPxd9orYkSZKkgSzq6k+SJEmSBmehQpIkSdJALFRIkiRJGoiFCkmSJEkDsVAhSZIkaSAWKqRZan2lfz3Jnj1pz03y3/OZL0kSJNmi50FjVyVZ0zP+jTnOy3FJLkzy2knpb+vJ12VJPptklyFu94NJntCGj01yaZL/v727C7GqCsM4/n+siynEIJUyQ4zRoUQzY0ZCC1SsKLoRpVKTkoIiseiD9KbRbqKisigsSpQkqaGkuihQREiLIE2lGc0uJhCDQBPMlBqQebtY63i2x5POuMeRcZ4fHOZ8vGftdQ7MWfvda+39dkham2tDFGNbJJ2UNC8/HiupoyZmpaTn8/3XJc3qq77awOekwuw8Rboe8xPAm5IaJA0FXgaWlGlX1SqmZmZ2niLiSORCY6TCX6uiWnhsWn/1Q9K1QEtE3BwRq+qEVPo1HmgDthYKt5bZ7nDgtojYlp/aANwITAKuIBU/rMReBrxK7wqrvgMsL9tPu3Q4qTArISI6SEVplgGtwPqI6JT0sKQf89Gn1ZXKoJI+kLRT0l5JpyqSSvpd0iuSdgNzLsqHMTMbJCQdz39nSPpW0leSfsu/wwvz73e7pMYcN1LSRkk78m16nTYbJK3L79staWZ+aTMwOo8Hd5ytXxHRluMX5DZb8/Y68vghSY2SdhW2O774uGAucGrmPCK+iYxU2O76QuxSUrG1Q+f88qrtHQCG56TJzEmFWR94iTQA3AO8JmkiKTGYlo+QXQ48mGOXR0QzMBm4s2aa+1BETImIz/qx72Zmg91k0qzzTcAioCkippIqOS/NMW+TZhRaSDvra+q0s4Q0iT0JmA98JKmBVJW7M89GbO9Bf3aRZhQA3o2IloiYSJpduC8iOoG/JN2SYxYD6+q0M51UEfo0ednTInLCIWk0acx6r04bjYUlY3tI31NtX89IsGxw8jILs5Ii4oSkNuB4RHRJmg20ADslQRoIDubw+ZIeJf3vXQdMAPbl19r6t+dmZgbsiIg/ACR1Ul0C1A5UZhtmAxPybzrAMElDI+J4oZ3bSUuCiIj9kg4ATcCxXvZHhfszJb0AXAlcDewlzY6vARZLehZ4AJhap51RwOE6z68GthUSnLeAZRHRXfh8FZ354FjqmLSy5vVDpLHMzEmFWR/pzjdIA8LaiHixGCBpPPA0MDUijkr6GGgohJzol56amVlRV+F+d+FxN9X9pCGk8xP+7Yf+TCEdlGogJQDNEXEw79BXxoyNwApgK/BTRByp084/nD7GIGkFMBJ4vPB0M/BpTihGAPdKOgns6UFfG/J2zLz8yewC2ALcL2kEnLoCyRhgGPA3cEzSKODui9hHMzPruc1Ul0JRWHpUtB1YmF9vAsYAv/ZmI5LmAncBn1BNCP7MFwKZV4nLyc0m0pKlekufAH4BxhXafow07syPiMpBMCLihogYGxFjgc+BJyPiyx52uQnoOGeUDQpOKsz6WES0k86z2CLpZ9JgdA1p7ek+YD+wHvj+onXSzMx64ymgOV8Wdh9nnlsAaVZhiKR20nLWRyKiq05crWcql5QFHjaxQH4AAAC/SURBVAJmRcThiDgKfEjaad8E7Kh53wbSbMr/XbHpa2BG4fH7pLHoh7y91rrv6qF8bsY4YGeZduzSoXQRADMzMzMbKHK9iKtql9rWxHxHOrn76AXY/hzg1rNt3wYXn1NhZmZmNoBI+gJoBM5VfO450jKsPk8qSPuQb1yAdm2A8kyFmZmZmZmV4nMqzMzMzMysFCcVZmZmZmZWipMKMzMzMzMrxUmFmZmZmZmV4qTCzMzMzMxKcVJhZmZmZmal/AdgP6Qf7dYLNAAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 936x360 with 2 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# Show the graphs\n",
    "fig, ax = plt.subplots(1, 2, figsize=(13, 5))\n",
    "ax[0].bar(years, num_events_by_year, tick_label=[str(year) if i % 3 == 0 else \"\" for i, year in enumerate(years)])\n",
    "ax[0].set_title(\"Number of Gun Violence Events by Year\")\n",
    "ax[0].set_xlabel(\"Year\")\n",
    "ax[0].set_ylabel(\"Number of Gun Violence Events\")\n",
    "\n",
    "ax[1].bar([tr[0] for tr in time_ranges], num_events_by_tod, width=4,\n",
    "        tick_label=[\"0:00-5:59\", \"6:00-11:59\", \"12:00-17:59\", \"18:00-23:59\"])\n",
    "ax[1].set_title(\"Number of Gun Violence Events by Time of Day\")\n",
    "ax[1].set_xlabel(\"Time of Day (24H)\")\n",
    "ax[1].set_ylabel(\"Number of Gun Violence Events\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The first graph indicates that our data is very heavily biased towards more recent events, and specifically 2016. No other data source cites a climb in gun violence of 800% from 2014 to 2016, so this bias is likely from more recent reports being available, or the data collection process.\n",
    "\n",
    "The second graph indicates that between 6PM and midnight is the most likely time for gun violence to occur. It is roughly twice as likely in that time range than in any other."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Filter such that we take only data from US after 2010\n",
    "The data outside the US is extremely sparse and would not contribute to our work\n",
    "\n",
    "The most recent census is from 2010 and removing data from before 2010 does not significantly impact our sample size. In addition, data from too long ago introduces a variable that we don't have enough data to account for (e.g. changes across decades)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {},
   "outputs": [],
   "source": [
    "usa_indices = []\n",
    "date_2010 = date(2010, 1, 1)\n",
    "for index, event in events.iterrows():\n",
    "    address = event[\"Address\"]\n",
    "    if address.endswith(\"USA\") and event[\"Date\"] > date_2010:\n",
    "        usa_indices.append(index)\n",
    "usa_events = events.iloc[usa_indices]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Write data to files"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Write clean data to CSV\n",
    "usa_events.to_csv(os.path.join(DATA_DIR, \"Events_USA_2010.tsv\"), sep=\"\\t\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Write county data to JSON\n",
    "with open(os.path.join(DATA_DIR, \"county_data.json\"), 'w') as _:\n",
    "    json.dump(county_data_map, _)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "2966"
      ]
     },
     "execution_count": 21,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "len(usa_events)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.6.4"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}