4. Get city loans year by year.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 31,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Importing packages...\n",
      "Now in directory: /home/idies/workspace/21cc/raddick/community_reinvestment_act\n",
      "ok\n"
     ]
    }
   ],
   "source": [
    "# BALTIMORE 288500\n",
    "debug = 1\n",
    "city = \"Pittsburgh\" # Possible values: Baltimore, Washington DC, Detroit, Newark, St. Louis, Richmond, San Francisco, Cleveland, Philadelphia, Pittsburgh\n",
    "g = 0\n",
    "\n",
    "latest_year = 2022\n",
    "earliest_year = 2012\n",
    "\n",
    "required_overlap = 0.5\n",
    "\n",
    "check_tract_and_city_boundaries = False\n",
    "check_tract_consistency = False\n",
    "\n",
    "cross_check_acs5_year_columns = False\n",
    "\n",
    "show_water = False\n",
    "show_roads = False\n",
    "\n",
    "print('Importing packages...')\n",
    "import os\n",
    "import pandas\n",
    "import numpy as np\n",
    "import time\n",
    "import geopandas\n",
    "from pprint import pprint\n",
    "from matplotlib import pyplot as plt\n",
    "from shapely.ops import unary_union\n",
    "from shapely.geometry import Point, box\n",
    "import matplotlib.patches as mpatches\n",
    "import requests\n",
    "import json\n",
    "from pprint import pprint\n",
    "import re\n",
    "import io\n",
    "\n",
    "pandas.set_option('display.max_colwidth', None)\n",
    "\n",
    "\n",
    "# Directories to look in\n",
    "thisdir = os.getcwd() + '/' #'/home/idies/workspace/Storage/raddick/Baltimore/community_reinvestment_act/'\n",
    "data_dir = '/home/idies/workspace/21cc/Data/'\n",
    "loans_dir = data_dir + 'Community Reinvestment Act/'\n",
    "jobs_dir = data_dir + 'LODES8/'\n",
    "\n",
    "census_dir = data_dir + 'Census/ACS5/'\n",
    "\n",
    "output_data_dir = thisdir + 'city_data/'\n",
    "if (not(os.path.exists(output_data_dir))):\n",
    "    os.makedirs(output_data_dir)\n",
    "\n",
    "#baltimore_dir = thisdir + 'baltimore/'\n",
    "\n",
    "census_shapefile_tiger_basedir = data_dir + 'Census/Shapefiles/TIGER/'\n",
    "\n",
    "#shapefile_dir = census_shapefile_tiger_basedir + '{0:.0f}/TRACT/'.format(thisyear)\n",
    "\n",
    "#acs5_basedir = '/home/idies/workspace/Temporary/raddick/census_scratch/acs5/'\n",
    "\n",
    "code_lookup_dir = thisdir + 'code_guide_lookups/'\n",
    "#inflation_dir = '/home/idies/workspace/Storage/raddick/Baltimore/community_reinvestment_act/datasets/inflation/'\n",
    "#extrasdir = '/home/idies/workspace/Storage/raddick/census/extras/'\n",
    "\n",
    "\n",
    "scale = 1\n",
    "\n",
    "equal_area_epsg = 5070\n",
    "\n",
    "os.chdir(thisdir)\n",
    "g = 0 # global time\n",
    "\n",
    "print('Now in directory: {0:}'.format(os.getcwd()))\n",
    "print('ok')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Selected Pittsburgh!\n"
     ]
    }
   ],
   "source": [
    "s = time.time()\n",
    "years = np.arange(latest_year, earliest_year-1, -1)\n",
    "\n",
    "if (city == 'Baltimore'):\n",
    "    thestate = 24\n",
    "    state_abbrev = 'md'\n",
    "    citycode = 4000\n",
    "    baltimore_shapefile_dir = '/home/idies/workspace/Storage/raddick/Baltimore/shapefiles/'\n",
    "    plotlimits = {'N': 39.38, 'S': 39.195, 'W': -76.745, 'E': -76.505}\n",
    "\n",
    "elif (city == 'Washington DC'):\n",
    "    thestate = 11  \n",
    "    state_abbrev = 'dc'\n",
    "    citycode = 50000\n",
    "    plotlimits = {'N': 39.01, 'S': 38.79, 'E': -76.9, 'W': -77.13}\n",
    "\n",
    "elif (city == 'Detroit'):\n",
    "    thestate = 26   # state_codes_df[state_codes_df['STATE_NAME'] == 'Missouri']\n",
    "    state_abbrev = 'mi'\n",
    "    citycode = 22000\n",
    "    plotlimits = {'N': 42.45, 'S': 42.25, 'E': -82.9, 'W': -83.3}\n",
    "    \n",
    "elif (city == 'Newark'):\n",
    "    thestate = 34\n",
    "    state_abbrev = 'nj'\n",
    "    citycode = 51000\n",
    "    plotlimits = {'N': 40.8, 'S': 40.6, 'E': -74.1, 'W': -74.3}\n",
    "\n",
    "elif (city == 'St. Louis'):\n",
    "    thestate = 29   # state_codes_df[state_codes_df['STATE_NAME'] == 'Missouri']\n",
    "    state_abbrev = 'mo'\n",
    "    citycode = 65000\n",
    "    plotlimits = {'N': 38.8, 'S': 38.5, 'E': -90.1, 'W': -90.4}\n",
    "    \n",
    "elif (city == 'Richmond'):\n",
    "    thestate = 51   # state_codes_df[state_codes_df['STATE_NAME'] == 'Missouri']\n",
    "    state_abbrev = 'va'\n",
    "    citycode = 67000\n",
    "    plotlimits = {'N': 37.61, 'S': 37.43, 'E': -77.37, 'W': -77.61}\n",
    "elif (city == 'San Francisco'):\n",
    "    thestate = 6   # state_codes_df[state_codes_df['STATE_NAME'] == 'California']\n",
    "    state_abbrev = 'ca'\n",
    "    plotlimits = {'N': 37.9, 'S': 37.65, 'E': -122.2, 'W': -123.2}  # with Farallon Islands\n",
    "    #plotlimits = {'N': 37.82, 'S': 37.7, 'E': -122.35, 'W': -122.53}  # without Farallon Islands\n",
    "    citycode = 67000    \n",
    "elif (city == 'Pittsburgh'):\n",
    "    thestate = 42\n",
    "    state_abbrev = 'pa'\n",
    "    citycode = 61000\n",
    "#    plotlimits = {'N': 40.55, 'S': 40.3, 'E': -79.85, 'W': -80.15}\n",
    "    plotlimits = {'N': 40.5, 'S': 40.35, 'E': -79.85, 'W': -80.1}\n",
    "\n",
    "elif (city == 'Cleveland'):\n",
    "    thestate = 39\n",
    "    state_abbrev = 'oh'\n",
    "    citycode = 16000\n",
    "    plotlimits = {'N': 41.65, 'S': 41.35, 'E': -81.5, 'W': -81.9}\n",
    "    \n",
    "elif (city == 'Philadelphia'):\n",
    "    thestate = 42\n",
    "    state_abbrev = 'pa'\n",
    "    #thecounty = 101\n",
    "    citycode = 60000\n",
    "    plotlimits = {'N': 40.15, 'S': 39.85, 'E': -74.93, 'W': -75.3}\n",
    "else:\n",
    "    print('ERROR: Select city from list!')\n",
    "\n",
    "# Initiailze water and roads whether we need them or not\n",
    "water_gdf = None\n",
    "roads_gdf = None\n",
    "    \n",
    "cityname_file = city.lower().replace(' ','_')\n",
    "e = time.time()\n",
    "g = g + (e-s)\n",
    "print('Selected {0:}!'.format(city))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Get water and roads if needed"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "# show_year = 2022\n",
    "# water_gdf = None\n",
    "# roads_gdf = None\n",
    "# show_water = True\n",
    "# show_roads = True\n",
    "\n",
    "# if ((show_water) & (water_gdf is not None)):\n",
    "#     s = time.time()\n",
    "#     print('reading water files...')\n",
    "#     water_dir = census_shapefile_tiger_basedir+\"{0:.0f}/AREAWATER/{1:}/\".format(show_year, state_abbrev)\n",
    "#     #print(water_dir)\n",
    "#     water_files = [water_dir+x for x in os.listdir(water_dir) if ('tl_{0:.0f}_{1:02d}'.format(show_year,thestate) in x) and (x[-3:] == 'shp')]\n",
    "#     water_gdf = geopandas.GeoDataFrame()\n",
    "#     for i in range(0, len(water_files)):\n",
    "#         if ((np.mod(i,10) == 0) | (i == len(water_files)-1)):\n",
    "#             print('\\treading water file {0:,.0f} of {1:,.0f}...'.format(i+1, len(water_files)))\n",
    "#         water_gdf = pandas.concat((water_gdf, geopandas.read_file(water_files[i])))\n",
    "#     e = time.time()\n",
    "#     g += e-s\n",
    "#     print('Got {0:,.0f} bodies of water in {1:,.1f} seconds!'.format(len(water_gdf), e-s))\n",
    "    \n",
    "\n",
    "# if ((show_roads) & (roads_gdf is not None)):\n",
    "#     s = time.time()\n",
    "#     print('\\n')\n",
    "#     print('reading road files...')\n",
    "#     road_dir = census_shapefile_tiger_basedir+\"{0:.0f}/ROADS/{1:}/\".format(show_year, state_abbrev)\n",
    "#     road_files = [road_dir+x for x in os.listdir(road_dir) if ('tl_{0:.0f}_{1:02d}'.format(show_year,thestate) in x) and (x[-3:] == 'shp')]\n",
    "#     roads_gdf = geopandas.GeoDataFrame()\n",
    "#     for i in range(0, len(road_files)):\n",
    "#         if ((np.mod(i,10) == 0) | (i == len(road_files)-1)):\n",
    "#             print('\\treading road file {0:,.0f} of {1:,.0f}...'.format(i+1, len(road_files)))\n",
    "#         roads_gdf = pandas.concat((roads_gdf, geopandas.read_file(road_files[i])))\n",
    "#     e = time.time()\n",
    "#     g += e-s\n",
    "#     print('Got {0:,.0f} roads in {1:,.1f} seconds!'.format(len(roads_gdf), e-s))\n",
    "# print('Done!')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Do city boundaries line up with census tract boundaries?\n",
    "\n",
    "Baltimore: Yes\n",
    "\n",
    "Washington, DC: Yes\n",
    "\n",
    "Detroit: No\n",
    "\n",
    "Cleveland: Almost\n",
    "\n",
    "Newark: Yes\n",
    "\n",
    "St. Louis: Yes\n",
    "\n",
    "Richmond: No (the place shape changed slightly 2016 -> 2017, but tracts unchanged)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "2022: found 3,446 total census tracts and 1,888 total places in pa...\n",
      "\tMatching Pittsburgh 2022 data (city area = 151.1 square km)...\n",
      "\n",
      "\n",
      "2021: found 3,446 total census tracts and 1,888 total places in pa...\n",
      "\tMatching Pittsburgh 2021 data (city area = 151.1 square km)...\n",
      "\n",
      "\n",
      "2020: found 3,446 total census tracts and 1,888 total places in pa...\n",
      "\tMatching Pittsburgh 2020 data (city area = 151.1 square km)...\n",
      "\n",
      "\n",
      "2019: found 3,218 total census tracts and 1,761 total places in pa...\n",
      "\tMatching Pittsburgh 2019 data (city area = 151.1 square km)...\n",
      "\n",
      "\n",
      "2018: found 3,218 total census tracts and 1,761 total places in pa...\n",
      "\tMatching Pittsburgh 2018 data (city area = 151.1 square km)...\n",
      "\n",
      "\n",
      "2017: found 3,218 total census tracts and 1,761 total places in pa...\n",
      "\tMatching Pittsburgh 2017 data (city area = 151.1 square km)...\n",
      "\n",
      "\n",
      "2016: found 3,218 total census tracts and 1,761 total places in pa...\n",
      "\tMatching Pittsburgh 2016 data (city area = 151.1 square km)...\n",
      "\n",
      "\n",
      "2015: found 3,218 total census tracts and 1,762 total places in pa...\n",
      "\tMatching Pittsburgh 2015 data (city area = 151.1 square km)...\n",
      "\n",
      "\n",
      "2014: found 3,218 total census tracts and 1,762 total places in pa...\n",
      "\tMatching Pittsburgh 2014 data (city area = 151.1 square km)...\n",
      "\n",
      "\n",
      "2013: found 3,218 total census tracts and 1,762 total places in pa...\n",
      "\tMatching Pittsburgh 2013 data (city area = 151.1 square km)...\n",
      "\n",
      "\n",
      "2012: found 3,218 total census tracts and 1,764 total places in pa...\n",
      "\tMatching Pittsburgh 2012 data (city area = 151.1 square km)...\n",
      "\n",
      "\n",
      "Pittsburgh: Found 1,480 tract-years in 61 seconds!\n",
      "Done!\n"
     ]
    }
   ],
   "source": [
    "s = time.time()\n",
    "#required_overlap = 0.5\n",
    "#check_tract_and_city_boundaries = True\n",
    "\n",
    "show_year = 2022\n",
    "\n",
    "show_water = False\n",
    "show_roads = False\n",
    "\n",
    "xlim = [plotlimits['W'], plotlimits['E']]\n",
    "ylim = [plotlimits['S'], plotlimits['N']]\n",
    "\n",
    "# xlim = [plotlimits['W'], -77.47]\n",
    "# ylim = [37.51, plotlimits['N']]\n",
    "\n",
    "city_tracts_years_gdf = geopandas.GeoDataFrame()\n",
    "partial_overlaps = {}\n",
    "\n",
    "plot_geo = box(xlim[0], ylim[0], xlim[1], ylim[1])\n",
    "if (show_water):\n",
    "    plot_gdf = geopandas.GeoDataFrame(data=[[plot_geo]], columns=['geometry'], crs=water_gdf.crs)\n",
    "else:\n",
    "    plot_gdf = geopandas.GeoDataFrame(data=[[plot_geo]], columns=['geometry'], crs=city_tracts_years_gdf.crs)\n",
    "    \n",
    "if (show_water):\n",
    "    water_here_gdf = geopandas.overlay(water_gdf, plot_gdf, how='intersection')\n",
    "if (show_roads):\n",
    "    roads_here_gdf = geopandas.overlay(roads_gdf, plot_gdf, how='intersection')\n",
    "\n",
    "\n",
    "for thisyear in years:\n",
    "    # Get tracts this year\n",
    "    if (thisyear == 2010):\n",
    "        state_tracts_this_year_gdf = geopandas.read_file(census_shapefile_tiger_basedir +'{0:.0f}/TRACT/tl_{0:.0f}_{1:02d}_tract10.shp'.format(thisyear, thestate))#.set_index('GEOID')\n",
    "        for x in state_tracts_this_year_gdf.columns[:-1]:\n",
    "            state_tracts_this_year_gdf = state_tracts_this_year_gdf.rename(columns = {x: x[:-2]})\n",
    "        state_tracts_this_year_gdf = state_tracts_this_year_gdf.set_index('GEOID')\n",
    "    else:\n",
    "        state_tracts_this_year_gdf = geopandas.read_file(census_shapefile_tiger_basedir +'{0:.0f}/TRACT/tl_{0:.0f}_{1:02d}_tract.shp'.format(thisyear, thestate)).set_index('GEOID')\n",
    "    state_tracts_this_year_gdf = state_tracts_this_year_gdf.sort_index()\n",
    "    \n",
    "    # Get places this year\n",
    "    if (thisyear == 2010):\n",
    "        state_places_this_year_gdf = geopandas.read_file(census_shapefile_tiger_basedir + '{0:.0f}/PLACE/tl_{0:.0f}_{1:02d}_place10.shp'.format(thisyear, thestate))#.set_index('GEOID')\n",
    "        for x in state_places_this_year_gdf.columns[:-1]:\n",
    "            state_places_this_year_gdf = state_places_this_year_gdf.rename(columns = {x: x[:-2]})\n",
    "        state_places_this_year_gdf = state_places_this_year_gdf.set_index('GEOID')\n",
    "    else:\n",
    "        state_places_this_year_gdf = geopandas.read_file(census_shapefile_tiger_basedir + '{0:.0f}/PLACE/tl_{0:.0f}_{1:02d}_place.shp'.format(thisyear, thestate)).set_index('GEOID')\n",
    "    state_places_this_year_gdf = state_places_this_year_gdf.sort_index()\n",
    "    \n",
    "    # Prepare to make plot later\n",
    "    if (check_tract_and_city_boundaries):\n",
    "        if (thisyear == show_year):\n",
    "            state_tracts_show_year_gdf = state_tracts_this_year_gdf\n",
    "            state_places_show_year_gdf = state_places_this_year_gdf\n",
    "        \n",
    "        \n",
    "    print('{0:.0f}: found {1:,.0f} total census tracts and {2:,.0f} total places in {3:}...'.format(thisyear, len(state_tracts_this_year_gdf), len(state_places_this_year_gdf), state_abbrev))\n",
    "\n",
    "    # Get this city only\n",
    "    this_city_this_year_gdf = state_places_this_year_gdf[state_places_this_year_gdf['PLACEFP'] == '{0:05d}'.format(citycode)]\n",
    "    \n",
    "    if (len(this_city_this_year_gdf.geometry.values) == 1):\n",
    "        city_year_geo = this_city_this_year_gdf.geometry.values[0]\n",
    "        city_year_area_sq_m = this_city_this_year_gdf.to_crs(epsg=equal_area_epsg).geometry.area.values[0]\n",
    "        print('\\tMatching {0:} {1:.0f} data (city area = {2:,.1f} square km)...'.format(city, thisyear, city_year_area_sq_m/(1000**2)))\n",
    "    else:\n",
    "        print('wtf, this city has more than one geometry included?')\n",
    "        break\n",
    "    \n",
    "    cnt = 0\n",
    "    this_year_partial_overlap_list = []\n",
    "    city_tracts_this_year_gdf = geopandas.GeoDataFrame(data=None, columns=state_tracts_this_year_gdf.columns, crs=state_tracts_this_year_gdf.crs, geometry='geometry')\n",
    "    \n",
    "    for ix, thisrow in state_tracts_this_year_gdf.iterrows():\n",
    "# #         if (np.mod(cnt,100) == 0):\n",
    "# #             print('\\t\\tMatching tract {0:,.0f} of {1:,.0f}...'.format(cnt, ntracts_this_state_year))\n",
    "        if (thisrow.geometry.intersects(city_year_geo)):        \n",
    "            intersection_geo = state_tracts_this_year_gdf.loc[ix].geometry.intersection(city_year_geo)\n",
    "            pct_overlap = intersection_geo.area / thisrow.geometry.area\n",
    "            if (pct_overlap >= 0.995):\n",
    "                city_tracts_this_year_gdf = pandas.concat((city_tracts_this_year_gdf, state_tracts_this_year_gdf[state_tracts_this_year_gdf.index == ix]), axis=0)\n",
    "            elif (pct_overlap >= 0.01):\n",
    "                intersection_gdf = geopandas.GeoDataFrame(data=[[intersection_geo]], columns=['geometry'], crs=state_tracts_this_year_gdf.crs, geometry='geometry')\n",
    "                tract_total_area_sq_m = state_tracts_this_year_gdf.to_crs(epsg=equal_area_epsg).loc[ix].geometry.area\n",
    "                intersection_area_sq_m = intersection_gdf.to_crs(epsg=equal_area_epsg).geometry.apply(lambda x: x.area).values[0]\n",
    "                pct_overlap_sq_m = intersection_area_sq_m / tract_total_area_sq_m        \n",
    "                print('\\n\\t\\tPARTIAL overlap found: County {0:} Tract {1:}\\n\\t\\t\\tTotal area {2:,.2f} sq km, overlap area {3:,.2f} sq km ({4:.1%})'.format(thisrow['COUNTYFP'], thisrow['TRACTCE'], tract_total_area_sq_m/1000000, intersection_area_sq_m/1000000, pct_overlap_sq_m))\n",
    "                this_year_partial_overlap_list.append(ix)\n",
    "                if (pct_overlap >= required_overlap):\n",
    "                    city_tracts_this_year_gdf = pandas.concat((city_tracts_this_year_gdf, state_tracts_this_year_gdf[state_tracts_this_year_gdf.index == ix]), axis=0)\n",
    "        cnt += 1\n",
    "        \n",
    "    city_tracts_this_year_gdf = city_tracts_this_year_gdf.assign(year = thisyear)\n",
    "    city_tracts_years_gdf = pandas.concat((city_tracts_years_gdf, city_tracts_this_year_gdf), axis=0)\n",
    "    \n",
    "    #if (len(this_year_partial_overlap_list) > 0):\n",
    "    partial_overlaps[thisyear] = this_year_partial_overlap_list\n",
    "    \n",
    "    print('\\n')\n",
    "    \n",
    "city_tracts_years_gdf.loc[:, 'STATEFP'] = pandas.to_numeric(city_tracts_years_gdf['STATEFP'], errors='coerce')\n",
    "city_tracts_years_gdf.loc[:, 'COUNTYFP'] = pandas.to_numeric(city_tracts_years_gdf['COUNTYFP'], errors='coerce')\n",
    "city_tracts_years_gdf.loc[:, 'TRACTCE'] = pandas.to_numeric(city_tracts_years_gdf['TRACTCE'].apply(lambda x: '{0:}.{1:}'.format(x[0:-2],x[-2:])), errors='coerce')\n",
    "for x in ['ALAND', 'AWATER', 'INTPTLAT','INTPTLON']:\n",
    "    city_tracts_years_gdf.loc[:, x] = pandas.to_numeric(city_tracts_years_gdf[x], errors='coerce')\n",
    "\n",
    "city_tracts_years_gdf = city_tracts_years_gdf.rename(columns={'COUNTYFP': 'county', 'TRACTCE': 'census_tract'})\n",
    "city_tracts_years_gdf.index.name = 'GEOID'\n",
    "city_tracts_years_gdf = city_tracts_years_gdf.reset_index()\n",
    "city_tracts_years_gdf.loc[:, 'GEOID'] = city_tracts_years_gdf['GEOID'].apply(lambda x: '14000US'+x)\n",
    "city_tracts_years_gdf = city_tracts_years_gdf.set_index(['GEOID','year'])\n",
    "\n",
    "e = time.time()\n",
    "g = g + (e-s)\n",
    "print('{0:}: Found {1:,.0f} tract-years in {2:,.0f} seconds!'.format(city, len(city_tracts_years_gdf), e-s))\n",
    "\n",
    "# PLOT!\n",
    "if (check_tract_and_city_boundaries):\n",
    "    s = time.time()\n",
    "    print('\\n')\n",
    "    print('plotting tract vs city for {0:.0f}...'.format(show_year))\n",
    "    \n",
    "    partial_overlap_tracts_show_year_gdf = state_tracts_show_year_gdf.loc[partial_overlaps[show_year]]\n",
    "    \n",
    "    fig, ax = plt.subplots(1,1, figsize=(12,12))\n",
    "    state_tracts_show_year_gdf.plot(color='none', edgecolor='black', ax=ax, linewidth=0.25)    \n",
    "    state_places_show_year_gdf[state_places_show_year_gdf['PLACEFP'] == '{0:05d}'.format(citycode)].plot(ax=ax, color='none', edgecolor='yellow', lw=5)\n",
    "\n",
    "    city_tracts_years_gdf.xs(show_year, level='year').plot(ax=ax, color='red', alpha=0.5, edgecolor='black', linewidth=0.5)\n",
    "    #city_tracts_years_gdf.xs(2016, level='year').plot(ax=ax, color='green', alpha=0.5, edgecolor='black', linewidth=0.5)\n",
    "        \n",
    "    if (show_water):\n",
    "        print('\\tplotting water...')\n",
    "        water_here_gdf.plot(ax=ax, color='blue')\n",
    "\n",
    "    if (show_roads):\n",
    "        print('\\tplotting major roads...')\n",
    "        roads_here_gdf[roads_here_gdf['RTTYP'] == 'I'].plot(ax=ax, color='purple', linewidth=5)\n",
    "        roads_here_gdf[roads_here_gdf['RTTYP'] == 'U'].plot(ax=ax, color='purple', linewidth=3)\n",
    "        roads_here_gdf[roads_here_gdf['RTTYP'] == 'S'].plot(ax=ax, color='purple', linewidth=2)\n",
    "        print('\\tplotting minor roads...')\n",
    "        roads_here_gdf[roads_here_gdf['RTTYP'] == 'M'].plot(ax=ax, color='purple', linewidth=1)\n",
    "        roads_here_gdf[roads_here_gdf['RTTYP'].isnull()].plot(ax=ax, color='purple', linewidth=0.5)\n",
    "        \n",
    "    if (len(partial_overlap_tracts_show_year_gdf) > 0):\n",
    "        partial_overlap_tracts_show_year_gdf.plot(ax=ax, color='none', edgecolor='black', linewidth=3, alpha=0.5)\n",
    "        for ix, thisrow in partial_overlap_tracts_show_year_gdf.iterrows():\n",
    "            annotator = '{0:}'.format(np.float(thisrow['NAME']))\n",
    "            ax.annotate(annotator, (thisrow.geometry.centroid.x, thisrow.geometry.centroid.y), (thisrow.geometry.centroid.x, thisrow.geometry.centroid.y), \n",
    "                        color='black', bbox=dict(fc='white', ec='black'))\n",
    "\n",
    "    ax.set_xlim(xlim)\n",
    "    ax.set_ylim(ylim)\n",
    "        \n",
    "    plt.title('{0:.0f}'.format(show_year), size=24)\n",
    "    \n",
    "    plt.show()\n",
    "    e = time.time()\n",
    "    g += e-s \n",
    "\n",
    "print('Done!')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Which census tracts changed, when?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "How many tracts per year?\n",
      "\t2012: found 137 tracts!\n",
      "\t2013: found 137 tracts!\n",
      "\t2014: found 137 tracts!\n",
      "\t2015: found 137 tracts!\n",
      "\t2016: found 137 tracts!\n",
      "\t2017: found 137 tracts!\n",
      "\t2018: found 137 tracts!\n",
      "\t2019: found 137 tracts!\n",
      "\t2020: found 128 tracts!\n",
      "\t2021: found 128 tracts!\n",
      "\t2022: found 128 tracts!\n",
      "Done!\n"
     ]
    }
   ],
   "source": [
    "show_water = False\n",
    "show_roads = False\n",
    "\n",
    "xlim = [plotlimits['W'], plotlimits['E']]\n",
    "ylim = [plotlimits['S'], plotlimits['N']]\n",
    "\n",
    "xlim = [-79.97, -79.9]\n",
    "ylim = [40.37, 40.43]\n",
    "\n",
    "print('How many tracts per year?')\n",
    "theyears = city_tracts_years_gdf.index.get_level_values(level='year').drop_duplicates().sort_values().tolist()\n",
    "for thisyear in theyears:\n",
    "    print('\\t{0:}: found {1:,.0f} tracts!'.format(thisyear, len(city_tracts_years_gdf.xs(thisyear, level='year'))))\n",
    "\n",
    "    \n",
    "if (check_tract_consistency):\n",
    "    plot_geo = box(xlim[0], ylim[0], xlim[1], ylim[1])\n",
    "    plot_gdf = geopandas.GeoDataFrame(data=[[plot_geo]], columns=['geometry'], crs=city_tracts_years_gdf.crs)\n",
    "    if (show_water):\n",
    "        print('Finding water in plot area...')\n",
    "        water_here_gdf = geopandas.overlay(water_gdf, plot_gdf, how='intersection')\n",
    "    if (show_roads):\n",
    "        print('Finding roads in plot area...')\n",
    "        roads_here_gdf = geopandas.overlay(roads_gdf, plot_gdf, how='intersection')\n",
    "    \n",
    "    print('\\n')\n",
    "    print('Finding tracts that changed between 2019 and 2020...')\n",
    "    tracts_2019_gdf = city_tracts_years_gdf.xs(2019, level='year')\n",
    "    tracts_2020_gdf = city_tracts_years_gdf.xs(2020, level='year')\n",
    "\n",
    "\n",
    "    tracts_2019_but_not_2020 = tracts_2019_gdf[~tracts_2019_gdf.index.isin(tracts_2020_gdf.index)]['NAME'].tolist()\n",
    "    print('tracts yeeted: {0:,.0f}'.format(len(tracts_2019_but_not_2020)))\n",
    "    tracts_2020_but_not_2019 = tracts_2020_gdf[~tracts_2020_gdf.index.isin(tracts_2019_gdf.index)]['NAME'].tolist()\n",
    "    print('tracts gained: {0:,.0f}'.format(len(tracts_2020_but_not_2019)))\n",
    "\n",
    "    print('\\n')\n",
    "    print('YEETED TRACTS:')\n",
    "    print(sorted([np.float(x) for x in tracts_2019_but_not_2020]))\n",
    "    print('----------')\n",
    "    print('ADDED TRACTS:')\n",
    "    print(sorted([np.float(x) for x in tracts_2020_but_not_2019]))\n",
    "    print('\\n')\n",
    "    \n",
    "    \n",
    "    print('Plotting tracts that changed between 2019 and 2020...')\n",
    "    fig, ax = plt.subplots(1,1, figsize=(18,18))\n",
    "\n",
    "#     tracts_2019_gdf.plot(ax=ax, color='none', edgecolor='black')\n",
    "#     tracts_2019_gdf[(tracts_2019_gdf['NAME'].isin(tracts_2019_but_not_2020)) \n",
    "# #                     & \n",
    "# #                     (\n",
    "# #                         (pandas.to_numeric(tracts_2019_gdf['NAME'], errors='coerce', downcast='float') >= 1200)\n",
    "# #                            & \n",
    "# #                         (pandas.to_numeric(tracts_2019_gdf['NAME'], errors='coerce', downcast='float') <= 1300)\n",
    "# #                      )\n",
    "#                 ].plot(ax=ax, color='none', edgecolor='black', linewidth=4)\n",
    "#     for ix, thisrow in tracts_2019_gdf[tracts_2019_gdf.geometry.apply(lambda x: x.intersects(plot_geo))].iterrows():\n",
    "#        if (thisrow['NAME'] in tracts_2019_but_not_2020):\n",
    "#             if ((np.float(thisrow['NAME']) == 3103) | (np.float(thisrow['NAME']) == 5653)):\n",
    "#                 annotator = '{0:}'.format(np.float(thisrow['NAME']))\n",
    "#                 ax.annotate(annotator, (thisrow.geometry.centroid.x, thisrow.geometry.centroid.y), (thisrow.geometry.centroid.x, thisrow.geometry.centroid.y), \n",
    "#                             color='black')#, bbox=dict(fc='white', ec='black'))\n",
    "\n",
    "    tracts_2020_gdf.plot(ax=ax, color='none', edgecolor='red')\n",
    "    tracts_2020_gdf[(tracts_2020_gdf['NAME'].isin(tracts_2020_but_not_2019))\n",
    "#                    & \n",
    "#                     (\n",
    "#                         (pandas.to_numeric(tracts_2020_gdf['NAME'], errors='coerce', downcast='float') >= 1200)\n",
    "#                            & \n",
    "#                         (pandas.to_numeric(tracts_2020_gdf['NAME'], errors='coerce', downcast='float') <= 1300)\n",
    "#                      )\n",
    "                  ].plot(ax=ax, color='none', edgecolor='red', linewidth=2, alpha=0.25, linestyle='dashed')#'black', alpha=0.25)            \n",
    "    for ix, thisrow in tracts_2020_gdf[tracts_2020_gdf.geometry.apply(lambda x: x.intersects(plot_geo))].iterrows():\n",
    "       if (thisrow['NAME'] in tracts_2020_but_not_2019):\n",
    "            if ((np.float(thisrow['NAME']) == 5629.01) | (np.float(thisrow['NAME']) == 5653)):\n",
    "                annotator = '{0:}'.format(np.float(thisrow['NAME']))\n",
    "                ax.annotate(annotator, (thisrow.geometry.centroid.x, thisrow.geometry.centroid.y), (thisrow.geometry.centroid.x, thisrow.geometry.centroid.y), \n",
    "                            color='red', bbox=dict(fc='white', ec='red'))\n",
    "\n",
    "    if (show_water):\n",
    "        print('\\tplotting water...')\n",
    "        water_here_gdf.plot(ax=ax, color='blue', alpha=0.5)\n",
    "\n",
    "    if (show_roads):\n",
    "        print('\\tplotting major roads...')\n",
    "        roads_here_gdf[roads_here_gdf['RTTYP'] == 'I'].plot(ax=ax, color='purple', linewidth=5)\n",
    "        roads_here_gdf[roads_here_gdf['RTTYP'] == 'U'].plot(ax=ax, color='purple', linewidth=3)\n",
    "        roads_here_gdf[roads_here_gdf['RTTYP'] == 'S'].plot(ax=ax, color='purple', linewidth=2)\n",
    "        print('\\tplotting minor roads...')\n",
    "        roads_here_gdf[roads_here_gdf['RTTYP'] == 'M'].plot(ax=ax, color='purple', linewidth=1)\n",
    "        roads_here_gdf[roads_here_gdf['RTTYP'].isnull()].plot(ax=ax, color='purple', linewidth=0.5)\n",
    "\n",
    "    ax.set_xlim(xlim)\n",
    "    ax.set_ylim(ylim)\n",
    "    plt.title('BLACK = OLD (2019); RED = NEW (2020)', fontsize=18)\n",
    "\n",
    "    plt.show()\n",
    "print('Done!')\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "\n",
    "# BALTIMORE 2019 -> 2020\n",
    "\n",
    "# MULTIPLE 2019 TRACTS COMBINED INTO SINGLE 2020 TRACT\n",
    "## Tracts 1801 and 1802 were combined into 2806\n",
    "\n",
    "# -----------------------------------------------------------------\n",
    "\n",
    "# WASHINGTON 2019 -> 2020\n",
    "\n",
    "# 2019 TRACTS SPLIT INTO MULTIPLE 2020 TRACTS\n",
    "### Tract 1.0 split into 1.01 and 1.02\n",
    "### Tract 7.01 split into 7.03 and 7.04\n",
    "### Tract 8.01 split into 8.03 and 8.04\n",
    "### Tract 9.01 split into 9.03 and 9.04\n",
    "### Tract 10.01 split into 10.03 and 10.04\n",
    "### Tract 13.02 split into 13.03 and 13.04\n",
    "### Tract 25.02 split into 25.03 and 25.04\n",
    "### Tract 27.01 split into 27.03 and 27.04\n",
    "### Tract 37.0 split into 37.01 and 37.02\n",
    "### Tract 38.0 split into 38.01 and 38.02\n",
    "### Tract 39.0 split into 39.01 and 39.02\n",
    "### Tract 44.0 split into 44.01 and 44.02\n",
    "### Tract 47.01 split into 47.03 and 47.04\n",
    "### Tract 50.02 split into 50.03 and 50.04\n",
    "### Tract 52.01 split into 52.02 and 52.03\n",
    "### Tract 53.01 split into 53.02 and 53.03\n",
    "### Tract 55.0 split into 55.01 and 55.02\n",
    "### Tract 56.0 split into 56.01 and 56.02\n",
    "### Tract 58.0 split into 58.01 and 58.02\n",
    "### Tract 72.0 split into 72.01 and 72.02 and 72.03\n",
    "### Tract 95.01 split into 95.10 and 95.11\n",
    "### Tract 102.0 split into 102.01 and 102.02\n",
    "### Tract 106.0 split into 106.01 and 106.02 and 106.03\n",
    "### Tract 110.0 split into 110.01 and 110.02\n",
    "\n",
    "# 2019 TRACTS RENAMED IN 2020\n",
    "### Tract 62.01 renamed to 9800\n",
    "\n",
    "# MINOR THINGS TO IGNORE\n",
    "### Tract 39.02 includes a small uninhabited part of 27.02\n",
    "\n",
    "# -----------------------------------------------------------------\n",
    "\n",
    "# DETROIT 2019 -> 2020\n",
    "\n",
    "# TRACTS NEWLY CREATED IN 2020\n",
    "### Tract 9841 was created for a zero-population area (docks/water)\n",
    "### Tract 9842 was created for a zero-population area (Detroit Iron & Metal)\n",
    "### Tract 9858 was created for a zero-population area (Rouge Park)\n",
    "\n",
    "# 2019 TRACTS RENAMED WITH CHANGED BORDERS FOR 2020\n",
    "### Tract 5240.0 became 5240.01 after zero-pop tract created (9842)\n",
    "### Tract 5241.0 became 5241.01 after zero-pop tract created (9841)\n",
    "### Tract 5462.0 became 5462.01 after Rouge Park became a zero-pop tract (9858)\n",
    "### Tract 5463.0 became 5463.01 after Rouge Park became a zero-pop tract (9858)\n",
    "\n",
    "# 2019 TRACTS SPLIT INTO MULTIPLE 2020 TRACTS\n",
    "### Tract 5315.0 split into 5315.01 and 5315.02\n",
    "\n",
    "# MULTIPLE 209 TRACTS COMBINED INTO SINGLE 2020 TRACT\n",
    "### Tracts 5047.0 and 5048.0 and 5055.0 combined into 5057.0\n",
    "### Tracts 5049.0 and 5050.0 combined into 5056.0 \n",
    "### Tracts 5078.0 and 5079.0 combined into 5090.0\n",
    "### Tracts 5104.0 and 5105.0 combined into 5091.0\n",
    "### Tracts 5107.0 and 5110.0 combined into 5058.0\n",
    "### Tracts 5121.0 and 5122.0 combined into 5026.0\n",
    "### Tracts 5123.0 and 5124.0 combined into 5027.0\n",
    "### Tracts 5126.0 and 5129.0 combined into 5128.0\n",
    "### Tracts 5136.0 and 5156.0 combined into 5138.0\n",
    "### Tracts 5159.0 and 5161.0 combined into 5192.0\n",
    "### Tracts 5162.0 and 5163.0 combined into 5193.0\n",
    "### Tracts 5184.0 and 5185.0 combined into 5191.0\n",
    "### Tracts 5186.0 and 5188.0 combined into 5190.0\n",
    "### Tracts 5213.0 and 5221.0 combined into 5228.0 \n",
    "### Tracts 5222.0 and 5272.0 and 5273.0 combined into 5279.0 \n",
    "### Tracts 5249.0 and 5250.0 combined into 5246.0 (but subratcting zero-pop 9841.0)\n",
    "### Tracts 5255.0 and 5256.0 combined into 5259.0 \n",
    "### Tracts 5322.0 and 5323.0 combined into 5321.0 \n",
    "### Tracts 5335.0 and 5337.0 combined into 5338.0 \n",
    "### Tracts 5345.0 and 5346.0 combined into 5348.0 \n",
    "### Tracts 5354.0 and 5355.0 combined into 5358.0 \n",
    "### Tracts 5427.0 and 5428.0 combined into 5433.0 \n",
    "### Tracts 5436.0 and 5438.0 combined into 5448.0 \n",
    "### Tracts 5453.0 and 5454.0 combined into 5472.0 \n",
    "### Tracts 5464.0 and 5465.0 combined into 5470.0 \n",
    "### Tracts 5468.0 and 5469.0 combined into 5471.0 \n",
    "\n",
    "# MINOR THINGS TO IGNORE\n",
    "### Tract 5279.0 no longer includes Ken-Mac Metals plant that was in 5273.0\n",
    "### Holy Cross Cemetery: North part of what was 5241.0 now in 5260.0\n",
    "\n",
    "# -----------------------------------------------------------------\n",
    "\n",
    "# CLEVELAND 2019 -> 2020\n",
    "\n",
    "# TRACTS NEWLY CREATED IN 2020\n",
    "### Tract 9802 was created for a zero-population area (Edgwater Park)\n",
    "### Tract 9809 was created for a zero-population area (CSX railyard)\n",
    "\n",
    "# TRACTS RENAMED FROM 2019 TO 2020\n",
    "### Tract 1012.0 renamed 1012.01 after carving out zero-pop 9802    \n",
    "    \n",
    "# 2019 TRACTS SPLIT INTO MULTIPLE 2020 TRACTS\n",
    "### Tract 1172.02 split into 1172.03 and 9809 (CSX railyard)\n",
    "\n",
    "# MULTIPLE 2019 TRACTS COMBINED INTO SINGLE 2020 TRACT\n",
    "### Tracts 1031.0 and 1034.0 combined into 1975.0\n",
    "### Tracts 1039.0 and 1041.0 combined into 1977.0\n",
    "### Tracts 1042.0 and 1043.0 combined into 1978.0\n",
    "### Tracts 1046.0 and 1049.0 combined into 1976.0\n",
    "### Tracts 1063.0 and 1064.0 combined into 1974.0\n",
    "### Tracts 1105.01 and 1108.01 combined into 1979.0\n",
    "### Tracts 1115.0 and 1116.0 combined into 1989.0\n",
    "### Tracts 1118.0 and 1119.02 combined into 1990.0\n",
    "### Tracts 1124.0 and 1128.0 combined into 1987.0\n",
    "### Tracts 1125.0 and 1126.0 combined into 1988.0\n",
    "### Tracts 1131.01 and 1965.0 combined into 1972.0\n",
    "### Tracts 1135.0 and 1136.0 combined into 1985.0\n",
    "### Tracts 1138.01 and 1141.0 combined into 1984.0\n",
    "### Tracts 1143.0 and 1147.0 combined into 1148.0\n",
    "### Tracts 1149.0 and 1153.0 combined into 1980.0\n",
    "### Tracts 1151.0 and 1152.0 combined into 1981.0\n",
    "### Tracts 1161.00 and 1162.00 were combined into 1992.00\n",
    "### Tracts 1184.0 and 1185.0 combined into 1991.0\n",
    "### Tracts 1187.0 and 1191.0 combined into 1968.0\n",
    "### Tracts 1192.02 and 1193.0 combined into 1986.0\n",
    "\n",
    "\n",
    "\n",
    "# -----------------------------------------------------------------\n",
    "\n",
    "# NEWARK 2019 -> 2020\n",
    "\n",
    "# 2019 TRACTS SPLIT INTO MULTIPLE 2020 TRACTS\n",
    "### Tract 22.01 split into 22.03 and 22.04\n",
    "\n",
    "# -----------------------------------------------------------------\n",
    "\n",
    "# ST. LOUIS 2019 -> 2020\n",
    "# MULTIPLE 2019 TRACTS COMBINED INTO SINGLE 2020 TRACT\n",
    "### Tracts 1114.0 and 1115.0 combined into 1277.0\n",
    "### Tracts 1184.0 and 1211.0 combined into 1278.0\n",
    "\n",
    "\n",
    "# -----------------------------------------------------------------\n",
    "\n",
    "# RICHMOND 2019 -> 2020\n",
    "\n",
    "# 2019 TRACTS SPLIT INTO MULTIPLE 2020 TRACTS\n",
    "### Tract 102 split into 102.01 and 102.02\n",
    "### Tract 205 split into 205.01 and 205.02\n",
    "### Tract 305 split into 305.01 and 305.02\n",
    "### Tract 402 split into 402.01 and 402.02\n",
    "### Tract 605 split into 605.01 and 605.02\n",
    "### Tract 610 split into 610.01 and 610.02\n",
    "### Tract 708.01 split into 708.03 and 708.04\n",
    "### Tract 709 split into 709.01 and 709.02\n",
    "### Tract 710.01 split into 710.03 and 710.04\n",
    "\n",
    "\n",
    "\n",
    "# SAN FRANCISCO 2019 -> 2020\n",
    "\n",
    "\n",
    "# NEW TRACTS\n",
    "### 179.03 created for Treasure Island and Yerba Buena Island\n",
    "\n",
    "# RENAMED TRACTS\n",
    "### 179.02 renamed 9902.0 (zero-pop SF Bay)\n",
    "\n",
    "# 2019 TRACTS SPLIT INTO MULTIPLE 2020 TRACTS\n",
    "### 101.0 split into 101.01 and 101.02\n",
    "### 102.0 split into 102.01 and 102.02\n",
    "### 104.0 split into 104.01 and 104.02\n",
    "### 107.0 split into 107.01 and 107.02\n",
    "### 109.0 split into 109.01 and 109.02\n",
    "### 110.0 split into 110.01 and 110.02\n",
    "### 111.0 split into 111.01 and 111.02\n",
    "### 120.0 split into 120.01 and 120.02\n",
    "### 122.01 split into 122.03 and 122.04\n",
    "### 124.01 split into 124.03 and 124.04\n",
    "### 124.02 split into 124.05 and 124.06\n",
    "### 125.01 split into 125.03 and 124.04\n",
    "### 128.0 split into 128.01 and 128.02\n",
    "### 130.0 split into 130.01 and 130.02\n",
    "### 134.0 split into 134.01 and 134.02\n",
    "### 152.0 split into 152.01 and 152.02\n",
    "### 154.0 split into 154.01 and 154.02\n",
    "### 157.0 split into 157.01 and 157.02\n",
    "### 161.0 split into 161.01 and 161.02\n",
    "### 166.0 split into 166.01 and 166.02\n",
    "### 176.01 split into 176.02 and 176.03 and 176.04\n",
    "### 178.02 split into 178.03 and 178.04\n",
    "### 201.0 split into 201.01 and 201.02\n",
    "### 202.0 split into 202.01 and 202.02\n",
    "### 206.0 split into 206.01 and 206.02\n",
    "### 207.0 split into 207.01 and 207.02\n",
    "### 208.0 split into 208.01 and 208.02\n",
    "### 255.0 split into 255.01 and 255.02\n",
    "### 262.0 split into 262.01 and 262.02\n",
    "### 314.0 split into 314.01 and 314.02\n",
    "### 330.0 split into 330.01 and 330.02\n",
    "### 351.0 split into 351.01 and 351.02\n",
    "### 452.0 split into 452.01 and 452.02\n",
    "### 479.01 split into 479.03 and 479.04\n",
    "### 607.0 split into 607.01 and 607.02 and 607.03\n",
    "### 611.0 split into 611.01 and 611.02\n",
    "### 614.0 split into 614.01 and 614.02\n",
    "### 615.0 split into 615.01 and 615.02 and 615.03 and 615.04 and 615.05 and 615.06 and 615.07 and 615.08\n",
    "\n",
    "# MINOR THINGS TO IGNORE\n",
    "# Yacht berths by Marina Airfield moved from 126.02 to 127\n",
    "\n",
    "\n",
    "# -----------------------------------------------------------------\n",
    "\n",
    "# PITTSBURGH 2019 -> 2020\n",
    "\n",
    "# SPLIT\n",
    "# 103.0 split into 103.01 and 103.02\n",
    "# 5629.0 split into 5629.01 + 5653\n",
    "# 5632.0 split into 5632.01 and 5632.02\n",
    "\n",
    "# COMBINED\n",
    "# 1016 and 1017 combined into 1019\n",
    "# 1204 and 1208 combined into 1209\n",
    "# 1207 and 1301 combined into 1307\n",
    "# 1303 and 1304 combined into 1308\n",
    "# 1406 and 1410 combined into 1412\n",
    "# 2107 and 2507 combined into 5652\n",
    "# 2206 and 2503 combined into 5651\n",
    "# 2406 and 2412 combined into 2413\n",
    "# 2609 and 2612 combined into 2613\n",
    "# 2704 and 2715 combined into 2716\n",
    "# 5616.0 and 5617.0 combined into 1610\n",
    "\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 37,
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "File written!\n"
     ]
    }
   ],
   "source": [
    "s = time.time() \n",
    "tract_changes = {}\n",
    "tract_changes['Baltimore'] = {}\n",
    "tract_changes['Baltimore']['combined'] = []\n",
    "tract_changes['Baltimore']['combined'].append({2019: ['180100', '180200'], 2020: ['280600']})\n",
    "\n",
    "# -----------------------------------------------------------------\n",
    "\n",
    "\n",
    "tract_changes['Washington DC'] = {}\n",
    "tract_changes['Washington DC']['renamed'] = []\n",
    "tract_changes['Washington DC']['renamed'].append({2019: ['006200'], 2020: ['980000']})\n",
    "\n",
    "tract_changes['Washington DC']['split'] = []\n",
    "tract_changes['Washington DC']['split'].append({2019: ['000100'], 2020: ['000101', '000102']})\n",
    "tract_changes['Washington DC']['split'].append({2019: ['000701'], 2020: ['000703', '000704']})\n",
    "tract_changes['Washington DC']['split'].append({2019: ['000801'], 2020: ['000803', '000804']})\n",
    "tract_changes['Washington DC']['split'].append({2019: ['001001'], 2020: ['001003', '001004']})\n",
    "tract_changes['Washington DC']['split'].append({2019: ['001302'], 2020: ['001303', '001304']})\n",
    "tract_changes['Washington DC']['split'].append({2019: ['002502'], 2020: ['002503', '002504']})\n",
    "tract_changes['Washington DC']['split'].append({2019: ['002701'], 2020: ['002703', '002704']})\n",
    "tract_changes['Washington DC']['split'].append({2019: ['003700'], 2020: ['003701', '003702']})\n",
    "tract_changes['Washington DC']['split'].append({2019: ['003800'], 2020: ['003801', '003802']})\n",
    "tract_changes['Washington DC']['split'].append({2019: ['003900'], 2020: ['003901', '003902']})\n",
    "tract_changes['Washington DC']['split'].append({2019: ['004400'], 2020: ['004401', '004402']})\n",
    "tract_changes['Washington DC']['split'].append({2019: ['004701'], 2020: ['004703', '004704']})\n",
    "tract_changes['Washington DC']['split'].append({2019: ['005002'], 2020: ['005003', '005004']})\n",
    "tract_changes['Washington DC']['split'].append({2019: ['005201'], 2020: ['005202', '005203']})\n",
    "tract_changes['Washington DC']['split'].append({2019: ['005301'], 2020: ['005302', '005303']})\n",
    "tract_changes['Washington DC']['split'].append({2019: ['005500'], 2020: ['005501', '005502']})\n",
    "tract_changes['Washington DC']['split'].append({2019: ['005600'], 2020: ['005601', '005602']})\n",
    "tract_changes['Washington DC']['split'].append({2019: ['005800'], 2020: ['005801', '005802']})\n",
    "tract_changes['Washington DC']['split'].append({2019: ['007200'], 2020: ['007201', '007202']})\n",
    "tract_changes['Washington DC']['split'].append({2019: ['009501'], 2020: ['009510', '009511']})\n",
    "tract_changes['Washington DC']['split'].append({2019: ['010200'], 2020: ['010201', '010202']})\n",
    "tract_changes['Washington DC']['split'].append({2019: ['010600'], 2020: ['010601', '010602']})\n",
    "tract_changes['Washington DC']['split'].append({2019: ['011000'], 2020: ['011001', '011002']})\n",
    "\n",
    "\n",
    "# -----------------------------------------------------------------\n",
    "\n",
    "tract_changes['Detroit'] = {}\n",
    "tract_changes['Detroit']['new'] = []\n",
    "tract_changes['Detroit']['new'].append({2019: [''], 2020: ['984100']})\n",
    "tract_changes['Detroit']['new'].append({2019: [''], 2020: ['984200']})\n",
    "tract_changes['Detroit']['new'].append({2019: [''], 2020: ['985800']})\n",
    "\n",
    "tract_changes['Detroit']['renamed'] = []\n",
    "tract_changes['Detroit']['renamed'].append({2019: ['524000'], 2020: ['524001']})\n",
    "tract_changes['Detroit']['renamed'].append({2019: ['524100'], 2020: ['524101']})\n",
    "tract_changes['Detroit']['renamed'].append({2019: ['546200'], 2020: ['546201']})\n",
    "tract_changes['Detroit']['renamed'].append({2019: ['546300'], 2020: ['546301']})\n",
    "\n",
    "tract_changes['Detroit']['split'] = []\n",
    "tract_changes['Detroit']['split'].append({2019: ['531500'], 2020: ['531501', '531502']})\n",
    "\n",
    "tract_changes['Detroit']['combined'] = []\n",
    "tract_changes['Detroit']['combined'].append({2019: ['504700', '5048700'], 2020: ['505700']})\n",
    "tract_changes['Detroit']['combined'].append({2019: ['504900', '505000'], 2020: ['505600']})\n",
    "tract_changes['Detroit']['combined'].append({2019: ['507800', '507900'], 2020: ['509000']})\n",
    "tract_changes['Detroit']['combined'].append({2019: ['510400', '510500'], 2020: ['509100']})\n",
    "tract_changes['Detroit']['combined'].append({2019: ['510700', '511100'], 2020: ['510800']})\n",
    "tract_changes['Detroit']['combined'].append({2019: ['511200', '512100'], 2020: ['502600']})\n",
    "tract_changes['Detroit']['combined'].append({2019: ['512300', '512400'], 2020: ['502700']})\n",
    "tract_changes['Detroit']['combined'].append({2019: ['512600', '512900'], 2020: ['512800']})\n",
    "tract_changes['Detroit']['combined'].append({2019: ['513600', '515600'], 2020: ['513800']})\n",
    "tract_changes['Detroit']['combined'].append({2019: ['515900', '516100'], 2020: ['519200']})\n",
    "tract_changes['Detroit']['combined'].append({2019: ['516200', '516300'], 2020: ['519300']})\n",
    "tract_changes['Detroit']['combined'].append({2019: ['518400', '518500'], 2020: ['519100']})\n",
    "tract_changes['Detroit']['combined'].append({2019: ['518600', '518800'], 2020: ['519000']})\n",
    "tract_changes['Detroit']['combined'].append({2019: ['521300', '522100'], 2020: ['522800']})\n",
    "tract_changes['Detroit']['combined'].append({2019: ['522200', '527200', '527300'], 2020: ['527900']})\n",
    "tract_changes['Detroit']['combined'].append({2019: ['524900', '525000'], 2020: ['524600']})\n",
    "tract_changes['Detroit']['combined'].append({2019: ['525500', '525600'], 2020: ['525900']})\n",
    "tract_changes['Detroit']['combined'].append({2019: ['532200', '532300'], 2020: ['532100']})\n",
    "tract_changes['Detroit']['combined'].append({2019: ['533500', '533700'], 2020: ['533800']})\n",
    "tract_changes['Detroit']['combined'].append({2019: ['534500', '534600'], 2020: ['534800']})\n",
    "tract_changes['Detroit']['combined'].append({2019: ['542700', '542800'], 2020: ['543300']})\n",
    "tract_changes['Detroit']['combined'].append({2019: ['543600', '543800'], 2020: ['544800']})\n",
    "tract_changes['Detroit']['combined'].append({2019: ['545300', '545400'], 2020: ['547200']})\n",
    "tract_changes['Detroit']['combined'].append({2019: ['546400', '546500'], 2020: ['547000']})\n",
    "tract_changes['Detroit']['combined'].append({2019: ['546800', '546900'], 2020: ['547100']})\n",
    "\n",
    "# -----------------------------------------------------------------\n",
    "\n",
    "tract_changes['Cleveland'] = {}\n",
    "tract_changes['Cleveland']['new'] = []\n",
    "tract_changes['Cleveland']['new'].append({2019: [''], 2020: ['980200']})\n",
    "tract_changes['Cleveland']['new'].append({2019: [''], 2020: ['980900']})\n",
    "\n",
    "tract_changes['Cleveland']['renamed'] = []\n",
    "tract_changes['Cleveland']['renamed'].append({2019: ['101200'], 2020: ['101201']})\n",
    "\n",
    "tract_changes['Cleveland']['split'] = []\n",
    "tract_changes['Cleveland']['split'].append({2019: ['117202'], 2020: ['117203', '980900']})\n",
    "\n",
    "tract_changes['Cleveland']['combined'] = []\n",
    "tract_changes['Cleveland']['combined'].append({2019: ['103100', '103400'], 2020: ['197500']})\n",
    "tract_changes['Cleveland']['combined'].append({2019: ['103900', '104100'], 2020: ['197700']})\n",
    "tract_changes['Cleveland']['combined'].append({2019: ['104200', '104300'], 2020: ['197800']})\n",
    "tract_changes['Cleveland']['combined'].append({2019: ['104600', '104900'], 2020: ['197600']})\n",
    "tract_changes['Cleveland']['combined'].append({2019: ['106300', '106400'], 2020: ['197400']})\n",
    "tract_changes['Cleveland']['combined'].append({2019: ['110501', '110801'], 2020: ['197900']})\n",
    "tract_changes['Cleveland']['combined'].append({2019: ['111500', '111600'], 2020: ['198900']})\n",
    "tract_changes['Cleveland']['combined'].append({2019: ['111800', '111902'], 2020: ['199000']})\n",
    "tract_changes['Cleveland']['combined'].append({2019: ['112400', '112800'], 2020: ['198700']})\n",
    "tract_changes['Cleveland']['combined'].append({2019: ['112500', '112600'], 2020: ['198800']})\n",
    "tract_changes['Cleveland']['combined'].append({2019: ['113101', '196500'], 2020: ['197200']})\n",
    "tract_changes['Cleveland']['combined'].append({2019: ['113500', '113600'], 2020: ['198500']})\n",
    "tract_changes['Cleveland']['combined'].append({2019: ['113801', '114100'], 2020: ['198400']})\n",
    "tract_changes['Cleveland']['combined'].append({2019: ['114300', '114700'], 2020: ['114800']})\n",
    "tract_changes['Cleveland']['combined'].append({2019: ['114900', '115300'], 2020: ['198000']})\n",
    "tract_changes['Cleveland']['combined'].append({2019: ['115100', '115200'], 2020: ['198100']})\n",
    "tract_changes['Cleveland']['combined'].append({2019: ['116100', '116200'], 2020: ['199200']})\n",
    "tract_changes['Cleveland']['combined'].append({2019: ['118400', '118500'], 2020: ['199100']})\n",
    "tract_changes['Cleveland']['combined'].append({2019: ['118700', '119100'], 2020: ['196800']})\n",
    "tract_changes['Cleveland']['combined'].append({2019: ['119202', '119300'], 2020: ['198600']})\n",
    "\n",
    "\n",
    "# -----------------------------------------------------------------\n",
    "\n",
    "tract_changes['Newark'] = {}\n",
    "tract_changes['Newark']['split'] = []\n",
    "tract_changes['Newark']['split'].append({2019: ['002201'], 2020: ['002203', '002204']})\n",
    "\n",
    "\n",
    "# -----------------------------------------------------------------\n",
    "\n",
    "tract_changes['St. Louis'] = {}\n",
    "tract_changes['St. Louis']['combined'] = []\n",
    "tract_changes['St. Louis']['combined'].append({2019: ['111400', '111500'], 2020: ['127700']})\n",
    "tract_changes['St. Louis']['combined'].append({2019: ['118400', '121100'], 2020: ['127800']})\n",
    "\n",
    "# -----------------------------------------------------------------\n",
    "\n",
    "tract_changes['Richmond'] = {}\n",
    "tract_changes['Richmond']['split'] = []\n",
    "tract_changes['Richmond']['split'].append({2019: ['010200'], 2020: ['010201', '010202']})\n",
    "tract_changes['Richmond']['split'].append({2019: ['020500'], 2020: ['020501', '020502']})\n",
    "tract_changes['Richmond']['split'].append({2019: ['030500'], 2020: ['030501', '030502']})\n",
    "tract_changes['Richmond']['split'].append({2019: ['040200'], 2020: ['040201', '040202']})\n",
    "tract_changes['Richmond']['split'].append({2019: ['060500'], 2020: ['060501', '060502']})\n",
    "tract_changes['Richmond']['split'].append({2019: ['061000'], 2020: ['061001', '061002']})\n",
    "tract_changes['Richmond']['split'].append({2019: ['070801'], 2020: ['070803', '070804']})\n",
    "tract_changes['Richmond']['split'].append({2019: ['070900'], 2020: ['070901', '070902']})\n",
    "tract_changes['Richmond']['split'].append({2019: ['071001'], 2020: ['071003', '071004']})\n",
    "\n",
    "# -----------------------------------------------------------------\n",
    "\n",
    "tract_changes['San Francisco'] = {}\n",
    "tract_changes['San Francisco']['new'] = []\n",
    "tract_changes['San Francisco']['new'].append({2019: [''], 2020: ['017903']})\n",
    "\n",
    "tract_changes['San Francisco']['renamed'] = []\n",
    "tract_changes['San Francisco']['renamed'].append({2019: ['017903'], 2020: ['990200']})\n",
    "\n",
    "tract_changes['San Francisco']['split'] = []\n",
    "tract_changes['San Francisco']['split'].append({2019: ['010100'], 2020: ['010101', '010102']})\n",
    "tract_changes['San Francisco']['split'].append({2019: ['010200'], 2020: ['010201', '010202']})\n",
    "tract_changes['San Francisco']['split'].append({2019: ['010400'], 2020: ['010401', '010402']})\n",
    "tract_changes['San Francisco']['split'].append({2019: ['010700'], 2020: ['010701', '010702']})\n",
    "tract_changes['San Francisco']['split'].append({2019: ['010900'], 2020: ['010901', '010902']})\n",
    "tract_changes['San Francisco']['split'].append({2019: ['011000'], 2020: ['011001', '011002']})\n",
    "tract_changes['San Francisco']['split'].append({2019: ['011100'], 2020: ['011101', '011102']})\n",
    "tract_changes['San Francisco']['split'].append({2019: ['012000'], 2020: ['012001', '012002']})\n",
    "tract_changes['San Francisco']['split'].append({2019: ['012201'], 2020: ['012203', '012204']})\n",
    "tract_changes['San Francisco']['split'].append({2019: ['012401'], 2020: ['012403', '012404']})\n",
    "tract_changes['San Francisco']['split'].append({2019: ['012402'], 2020: ['012405', '012406']})\n",
    "tract_changes['San Francisco']['split'].append({2019: ['012501'], 2020: ['012503', '012504']})\n",
    "tract_changes['San Francisco']['split'].append({2019: ['012800'], 2020: ['012801', '012802']})\n",
    "tract_changes['San Francisco']['split'].append({2019: ['013000'], 2020: ['013001', '013002']})\n",
    "tract_changes['San Francisco']['split'].append({2019: ['013000'], 2020: ['013401', '013402']})\n",
    "tract_changes['San Francisco']['split'].append({2019: ['015200'], 2020: ['015201', '015202']})\n",
    "tract_changes['San Francisco']['split'].append({2019: ['015700'], 2020: ['015701', '015702']})\n",
    "tract_changes['San Francisco']['split'].append({2019: ['016100'], 2020: ['016100', '016100']})\n",
    "tract_changes['San Francisco']['split'].append({2019: ['016600'], 2020: ['016600', '016600']})\n",
    "tract_changes['San Francisco']['split'].append({2019: ['017601'], 2020: ['017602', '017603', '017604']})\n",
    "tract_changes['San Francisco']['split'].append({2019: ['017801'], 2020: ['017802', '017803']})\n",
    "tract_changes['San Francisco']['split'].append({2019: ['020100'], 2020: ['020101', '020102']})\n",
    "tract_changes['San Francisco']['split'].append({2019: ['020200'], 2020: ['020201', '020202']})\n",
    "tract_changes['San Francisco']['split'].append({2019: ['020600'], 2020: ['020601', '020602']})\n",
    "tract_changes['San Francisco']['split'].append({2019: ['020700'], 2020: ['020701', '020702']})\n",
    "tract_changes['San Francisco']['split'].append({2019: ['020800'], 2020: ['020801', '020802']})\n",
    "tract_changes['San Francisco']['split'].append({2019: ['025500'], 2020: ['025501', '025502']})\n",
    "tract_changes['San Francisco']['split'].append({2019: ['026200'], 2020: ['026201', '026202']})\n",
    "tract_changes['San Francisco']['split'].append({2019: ['031400'], 2020: ['031401', '031402']})\n",
    "tract_changes['San Francisco']['split'].append({2019: ['033000'], 2020: ['033001', '033002']})\n",
    "tract_changes['San Francisco']['split'].append({2019: ['035100'], 2020: ['035101', '035102']})\n",
    "tract_changes['San Francisco']['split'].append({2019: ['045200'], 2020: ['045201', '045202']})\n",
    "tract_changes['San Francisco']['split'].append({2019: ['047901'], 2020: ['047903', '047904']})\n",
    "tract_changes['San Francisco']['split'].append({2019: ['060700'], 2020: ['060701', '060702', '060703']})\n",
    "tract_changes['San Francisco']['split'].append({2019: ['061100'], 2020: ['061101', '061102']})\n",
    "tract_changes['San Francisco']['split'].append({2019: ['061100'], 2020: ['061401', '061402']})\n",
    "tract_changes['San Francisco']['split'].append({2019: ['061500'], 2020: ['061501', '061502', '061503', '061504', '061505', '061506', '061507', '061508']})\n",
    "\n",
    "# -----------------------------------------------------------------\n",
    "\n",
    "tract_changes['Philadelphia'] = {}\n",
    "\n",
    "tract_changes['Philadelphia']['new'] = []\n",
    "tract_changes['Philadelphia']['new'].append({2019: '', 2020: ['989300']})\n",
    "\n",
    "tract_changes['Philadelphia']['renamed'] = []\n",
    "tract_changes['Philadelphia']['renamed'].append({2019: ['005000'], 2020: ['989200']})\n",
    "tract_changes['Philadelphia']['renamed'].append({2019: ['006900'], 2020: ['980900']})\n",
    "\n",
    "tract_changes['Philadelphia']['split'] = []\n",
    "tract_changes['Philadelphia']['split'].append({2019: ['000100'], 2020: ['000101', '000102']})\n",
    "tract_changes['Philadelphia']['split'].append({2019: ['000402'], 2020: ['000403', '000404']})\n",
    "tract_changes['Philadelphia']['split'].append({2019: ['000700'], 2020: ['000701', '000702']})\n",
    "tract_changes['Philadelphia']['split'].append({2019: ['000804'], 2020: ['000805', '000806']})\n",
    "tract_changes['Philadelphia']['split'].append({2019: ['001202'], 2020: ['001203', '001204']})\n",
    "tract_changes['Philadelphia']['split'].append({2019: ['001300'], 2020: ['001301', '001302']})\n",
    "tract_changes['Philadelphia']['split'].append({2019: ['004102'], 2020: ['004103', '004104']})\n",
    "tract_changes['Philadelphia']['split'].append({2019: ['012500'], 2020: ['012501', '012502']})\n",
    "tract_changes['Philadelphia']['split'].append({2019: ['013700'], 2020: ['013701', '013702']})\n",
    "tract_changes['Philadelphia']['split'].append({2019: ['014200'], 2020: ['014201', '014202']})\n",
    "tract_changes['Philadelphia']['split'].append({2019: ['016000'], 2020: ['016001', '016002']})\n",
    "tract_changes['Philadelphia']['split'].append({2019: ['018800'], 2020: ['018801', '018802']})\n",
    "tract_changes['Philadelphia']['split'].append({2019: ['020700'], 2020: ['020701', '020702']})\n",
    "tract_changes['Philadelphia']['split'].append({2019: ['036900'], 2020: ['036901', '036902']})\n",
    "tract_changes['Philadelphia']['split'].append({2019: ['038300'], 2020: ['038301', '038302']})\n",
    "tract_changes['Philadelphia']['split'].append({2019: ['039000'], 2020: ['039001', '039002']})\n",
    "tract_changes['Philadelphia']['split'].append({2019: ['980000'], 2020: ['980001', '980002', '980003']})\n",
    "tract_changes['Philadelphia']['split'].append({2019: ['980700'], 2020: ['980701', '980702']})\n",
    "tract_changes['Philadelphia']['split'].append({2019: ['980900'], 2020: ['980901', '980902', '980903', '980904', '980905', '980906']})\n",
    "\n",
    "# -----------------------------------------------------------------\n",
    "\n",
    "tract_changes['Pittsburgh'] = {}\n",
    "tract_changes['Pittsburgh']['split'] = []\n",
    "tract_changes['Pittsburgh']['split'].append({2019: ['010300'], 2020: ['010301', '010302']})\n",
    "tract_changes['Pittsburgh']['split'].append({2019: ['562900'], 2020: ['562901', '565300']})\n",
    "tract_changes['Pittsburgh']['split'].append({2019: ['563200'], 2020: ['563201', '563202']})\n",
    "\n",
    "tract_changes['Pittsburgh']['combined'] = []\n",
    "tract_changes['Pittsburgh']['combined'].append({2019: ['101600', '101700'], 2020: ['101900']})\n",
    "tract_changes['Pittsburgh']['combined'].append({2019: ['120400', '120800'], 2020: ['120900']})\n",
    "tract_changes['Pittsburgh']['combined'].append({2019: ['120700', '130100'], 2020: ['130700']})\n",
    "tract_changes['Pittsburgh']['combined'].append({2019: ['130300', '130400'], 2020: ['130800']})\n",
    "tract_changes['Pittsburgh']['combined'].append({2019: ['140600', '141000'], 2020: ['141200']})\n",
    "tract_changes['Pittsburgh']['combined'].append({2019: ['210700', '250700'], 2020: ['565200']})\n",
    "tract_changes['Pittsburgh']['combined'].append({2019: ['220600', '250300'], 2020: ['565100']})\n",
    "tract_changes['Pittsburgh']['combined'].append({2019: ['260900', '261200'], 2020: ['261300']})\n",
    "tract_changes['Pittsburgh']['combined'].append({2019: ['270400', '271500'], 2020: ['271600']})\n",
    "tract_changes['Pittsburgh']['combined'].append({2019: ['561600', '561700'], 2020: ['161000']})\n",
    "\n",
    "#pprint(tract_changes)\n",
    "with open(output_data_dir+\"changes_to_tracts.json\", \"w\") as outfile: \n",
    "    json.dump(tract_changes, outfile)\n",
    "    \n",
    "with open(output_data_dir+\"changes_to_tracts.json\", \"r\") as readthefile: \n",
    "    y = json.load(readthefile)\n",
    "\n",
    "e = time.time()\n",
    "g += e-s\n",
    "print('File written!')\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Get geo-aggregated loans for this city\n",
    "\n",
    "and join the tract data onto the loans data\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "reading nationwide data...\n",
      "Read 3,686,478 nationwide tract-years in 24 seconds...\n",
      "converting columns to numeric...\n",
      "calculating total loans...\n",
      "\n",
      "\n",
      "Kept 1,480 tract-years in 35 seconds!\n"
     ]
    }
   ],
   "source": [
    "s = time.time()\n",
    "\n",
    "print('reading nationwide data...')\n",
    "agg_loans_all_df = pandas.read_csv(loans_dir+'agg_loans.csv', encoding='utf-8', low_memory=False, index_col='rownumber', keep_default_na=False)\n",
    "agg_loans_all_df = agg_loans_all_df.rename(columns={'activity_year': 'year'})\n",
    "e = time.time()\n",
    "print('Read {0:,.0f} nationwide tract-years in {1:,.0f} seconds...'.format(len(agg_loans_all_df), e-s))\n",
    "\n",
    "# Keep only business loans\n",
    "agg_loans_all_df = agg_loans_all_df[agg_loans_all_df['loan_type'] == 4]\n",
    "# Keep only loan originations\n",
    "agg_loans_all_df = agg_loans_all_df[agg_loans_all_df['action_taken_type'] == 1]\n",
    "# Keep only this state\n",
    "agg_loans_all_df = agg_loans_all_df[agg_loans_all_df['state'] == thestate]\n",
    "\n",
    "\n",
    "# Construct all GEOIDs\n",
    "agg_loans_all_df = agg_loans_all_df.assign(GEOID = np.nan)\n",
    "\n",
    "agg_loans_all_df.loc[agg_loans_all_df['census_tract'].apply(lambda x: np.floor(x) == x), 'GEOID'] = agg_loans_all_df[\n",
    "    agg_loans_all_df['census_tract'].apply(lambda x: np.floor(x) == x)\n",
    "].apply(lambda row: '14000US{0:02d}{1:03d}{2:04d}{3:02d}'.format(\n",
    "    row['state'], row['county'], \n",
    "    int(np.floor(row['census_tract'])),\n",
    "    int(np.round(row['census_tract'] - np.floor(row['census_tract']),0)),\n",
    "    row['census_tract']\n",
    "), axis=1)#[['state','county','census_tract']]#.drop_duplicates().sort_values(by='census_tract')\n",
    "\n",
    "agg_loans_all_df.loc[agg_loans_all_df['census_tract'].apply(lambda x: np.floor(10*x) != 10*x), 'GEOID'] = agg_loans_all_df[\n",
    "    agg_loans_all_df['census_tract'].apply(lambda x: np.floor(10*x) != 10*x)\n",
    "].apply(lambda row: '14000US{0:02d}{1:03d}{2:04d}{3:02d}'.format(\n",
    "    row['state'], row['county'], \n",
    "    int(np.floor(row['census_tract'])),\n",
    "    int(np.round(100*(row['census_tract'] - np.floor(row['census_tract'])),0))\n",
    "), axis=1)\n",
    "\n",
    "geoids_so_far_list = agg_loans_all_df['GEOID'].dropna().drop_duplicates().sort_values()\n",
    "\n",
    "agg_loans_all_df.loc[(agg_loans_all_df['GEOID'].isnull()), 'GEOID'] = agg_loans_all_df[\n",
    "    (agg_loans_all_df['GEOID'].isnull()) \n",
    "].apply(lambda row: \n",
    "        '14000US{0:02d}{1:03d}{2:04d}{3:02d}'.format(\n",
    "            row['state'], row['county'], \n",
    "            int(np.floor(row['census_tract'])),\n",
    "            int(np.round(100*(row['census_tract'] - np.floor(row['census_tract'])),0))\n",
    "            ), axis=1\n",
    "       )\n",
    "\n",
    "# Done building GEOIDs! Now let's join to city tracts we found earlier\n",
    "agg_loans_all_df = agg_loans_all_df.set_index(['GEOID','year'])\n",
    "\n",
    "data_gdf = city_tracts_years_gdf.join(agg_loans_all_df[[x for x in agg_loans_all_df.columns if x not in city_tracts_years_gdf.columns]], how='left')\n",
    "\n",
    "\n",
    "\n",
    "print('converting columns to numeric...')\n",
    "numeric_columns = []\n",
    "numeric_columns += ['nLoans1', 'amtLoans1', 'nLoans100k', 'amtLoans100k']\n",
    "numeric_columns += ['nLoans250k', 'amtLoans250k', 'nLoansToSmallest', 'amtLoansToSmallest']\n",
    "\n",
    "for x in numeric_columns:\n",
    "    data_gdf.loc[:, x] = pandas.to_numeric(data_gdf[x], errors='coerce')\n",
    "\n",
    "print('calculating total loans...')\n",
    "data_gdf = data_gdf.assign(nLoans = data_gdf['nLoans1'] + data_gdf['nLoans100k'] + data_gdf['nLoans250k'])\n",
    "data_gdf = data_gdf.assign(amtLoans = data_gdf['amtLoans1'] + data_gdf['amtLoans100k'] + data_gdf['amtLoans250k'])\n",
    "\n",
    "# print('backing up...')\n",
    "# data_gdf_bk = data_gdf\n",
    "\n",
    "e = time.time()\n",
    "g = g + (e-s)\n",
    "print('\\n')\n",
    "print('Kept {0:,.0f} tract-years in {1:,.0f} seconds!'.format(len(data_gdf), e-s))\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Add income groups, CRA levels, working loans for each tract-year"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "looking up income group names from income_group_total...\n",
      "Adding CRA income levels (low/moderate/middle/upper/unknown)...\n",
      "Getting CRA income levels for tracts where only CRA level was reported...\n",
      "calculating working loans...\n",
      "Kept 1,480 tract-years in Pittsburgh in 0.12 seconds!\n"
     ]
    }
   ],
   "source": [
    "s = time.time()\n",
    "\n",
    "# print('getting from backup...')\n",
    "# data_gdf = data_gdf_bk\n",
    "\n",
    "print('looking up income group names from income_group_total...')\n",
    "\n",
    "data_gdf = data_gdf.rename(columns = {'income_group_total': 'income_group_code'})\n",
    "data_gdf = data_gdf.assign(income_group = np.nan)\n",
    "\n",
    "data_gdf.loc[data_gdf['income_group_code'] == 1, 'income_group'] = '< 10% of Median Family Income (MFI)'\n",
    "data_gdf.loc[data_gdf['income_group_code'] == 2, 'income_group'] = '10% to 20% of MFI'\n",
    "data_gdf.loc[data_gdf['income_group_code'] == 3, 'income_group'] = '20% to 30% of MFI'\n",
    "data_gdf.loc[data_gdf['income_group_code'] == 4, 'income_group'] = '30% to 40% of MFI'\n",
    "data_gdf.loc[data_gdf['income_group_code'] == 5, 'income_group'] = '40% to 50% of MFI'\n",
    "data_gdf.loc[data_gdf['income_group_code'] == 6, 'income_group'] = '50% to 60% of MFI'\n",
    "data_gdf.loc[data_gdf['income_group_code'] == 7, 'income_group'] = '60% to 70% of MFI'\n",
    "data_gdf.loc[data_gdf['income_group_code'] == 8, 'income_group'] = '70% to 80% of MFI'\n",
    "data_gdf.loc[data_gdf['income_group_code'] == 9, 'income_group'] = '80% to 90% of MFI'\n",
    "data_gdf.loc[data_gdf['income_group_code'] == 10, 'income_group'] = '90% to 100% of MFI'\n",
    "data_gdf.loc[data_gdf['income_group_code'] == 11, 'income_group'] = '100% to 110% of MFI'\n",
    "data_gdf.loc[data_gdf['income_group_code'] == 12, 'income_group'] = '110% to 120% of MFI'\n",
    "data_gdf.loc[data_gdf['income_group_code'] == 13, 'income_group'] = '> 120% of MFI'\n",
    "data_gdf.loc[data_gdf['income_group_code'] == 14, 'income_group'] = 'unknown'\n",
    "\n",
    "print('Adding CRA income levels (low/moderate/middle/upper/unknown)...')\n",
    "# Get levels (low, moderate, middle, upper)\n",
    "data_gdf = data_gdf.assign(cra_level = np.nan)\n",
    "data_gdf.loc[(data_gdf['income_group_code'] >= 1) & (data_gdf['income_group_code'] <= 5), 'cra_level'] = 'low'\n",
    "data_gdf.loc[(data_gdf['income_group_code'] >= 6) & (data_gdf['income_group_code'] <= 8), 'cra_level'] = 'moderate'\n",
    "data_gdf.loc[(data_gdf['income_group_code'] >= 9) & (data_gdf['income_group_code'] <= 12), 'cra_level'] = 'middle'\n",
    "data_gdf.loc[(data_gdf['income_group_code'] == 13), 'cra_level'] = 'upper'\n",
    "data_gdf.loc[(data_gdf['income_group_code'] == 14), 'cra_level'] = 'unknown'\n",
    "\n",
    "print('Getting CRA income levels for tracts where only CRA level was reported...')\n",
    "data_gdf.loc[data_gdf['income_group_code'] == 101, 'cra_level'] = 'low'\n",
    "data_gdf.loc[data_gdf['income_group_code'] == 102, 'cra_level'] = 'moderate'\n",
    "data_gdf.loc[data_gdf['income_group_code'] == 103, 'cra_level'] = 'middle'\n",
    "data_gdf.loc[data_gdf['income_group_code'] == 104, 'cra_level'] = 'upper'\n",
    "data_gdf.loc[data_gdf['income_group_code'] == 105, 'cra_level'] = 'unknown'\n",
    "\n",
    "print('calculating working loans...')\n",
    "data_gdf = data_gdf.assign(avgSmallLoan = data_gdf['amtLoans1'] / data_gdf['nLoans1'])\n",
    "\n",
    "data_gdf = data_gdf.assign(nWorkingLoans = 0)\n",
    "data_gdf.loc[data_gdf['avgSmallLoan'] < 10000, \n",
    "                           'nWorkingLoans'] = data_gdf['nLoans'][data_gdf['avgSmallLoan'] < 10000] - data_gdf['nLoans1'][data_gdf['avgSmallLoan'] < 10000]\n",
    "data_gdf.loc[data_gdf['avgSmallLoan'] >= 10000, \n",
    "                           'nWorkingLoans'] = data_gdf['nLoans'][data_gdf['avgSmallLoan'] >= 10000]\n",
    "\n",
    "data_gdf = data_gdf.assign(amtWorkingLoans = 0)\n",
    "data_gdf.loc[data_gdf['avgSmallLoan'] < 10000, \n",
    "                           'amtWorkingLoans'] = data_gdf['amtLoans'][data_gdf['avgSmallLoan'] < 10000] - data_gdf['amtLoans1'][data_gdf['avgSmallLoan'] < 10000]\n",
    "data_gdf.loc[data_gdf['avgSmallLoan'] >= 10000, \n",
    "                           'amtWorkingLoans'] = data_gdf['amtLoans'][data_gdf['avgSmallLoan'] >= 10000]\n",
    "\n",
    "e = time.time()\n",
    "g = g + (e-s)\n",
    "\n",
    "# print('backing up...')\n",
    "# data_gdf_bk = data_gdf\n",
    "\n",
    "print('Kept {0:,.0f} tract-years in {1:} in {2:,.2f} seconds!'.format(len(data_gdf), city, e-s))\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Connect to jobs data"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Get raw jobs table and aggregate by year"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "getting raw jobs files...\n",
      "documenting job type and segment...\n",
      "getting geo crosswalk file...\n",
      "getting GEOIDs from crosswalk file...\n",
      "selecting jobs in city...\n",
      "finding columns to sum over...\n",
      "summing jobs values into tract-years...\n",
      "Extrapolating jobs for recent years needed...\n",
      "\t2022\n",
      "calculating total and small business jobs...\n",
      "joining to loan data...\n",
      "\n",
      "Added jobs data for 1,480 tract-years in 1 minutes 1 seconds!\n"
     ]
    }
   ],
   "source": [
    "s = time.time()\n",
    "print('getting raw jobs files...')\n",
    "raw_jobs_df = pandas.DataFrame()\n",
    "\n",
    "wac_file_list = sorted([\n",
    "    '{0:}{1:}/wac/'.format(jobs_dir,state_abbrev)+x for x \n",
    "    in os.listdir('{0:}{1:}/wac'.format(jobs_dir, state_abbrev)) if \n",
    "    (x[-4:] == '.csv') and (int(x[-8:-4]) in years) and ('_S000_JT02' in x)]\n",
    ")\n",
    "\n",
    "cnt = 0\n",
    "for this_file in wac_file_list:\n",
    "    #if (np.mod(cnt, 50) == 0):\n",
    "#    print('Reading file {0:,.0f} of {1:,.0f}...'.format(cnt, len(wac_file_list)))\n",
    "    raw_jobs_df_i = pandas.read_csv(this_file)\n",
    "    raw_jobs_df_i = raw_jobs_df_i.assign(job_segment = this_file[-13:-9])\n",
    "    raw_jobs_df_i = raw_jobs_df_i.assign(job_type = this_file[-18:-14])\n",
    "    raw_jobs_df_i = raw_jobs_df_i.assign(year = int(this_file[-8:-4]))\n",
    "    raw_jobs_df = pandas.concat((raw_jobs_df, raw_jobs_df_i), axis=0)\n",
    "    cnt += 1\n",
    "\n",
    "raw_jobs_df = raw_jobs_df.reset_index(drop=True)\n",
    "\n",
    "\n",
    "print('documenting job type and segment...')\n",
    "job_type_metadata_list = []\n",
    "job_type_metadata_list.append({'variable': 'w_geocode', 'description': 'Workplace Census Block Code'})\n",
    "job_type_metadata_list.append({'variable': 'h_geocode', 'description': 'Workplace Residence Block Code'})\n",
    "job_type_metadata_list.append({'variable': 'S000', 'description': 'Total number of jobs'})\n",
    "job_type_metadata_list.append({'variable': 'SA01', 'description': 'Number of jobs of workers age 29 or younger'})\n",
    "job_type_metadata_list.append({'variable': 'SA02', 'description': 'Number of jobs for workers age 30 to 54'})\n",
    "job_type_metadata_list.append({'variable': 'SA03', 'description': 'Number of jobs for workers age 55 or older'})\n",
    "job_type_metadata_list.append({'variable': 'SE01', 'description': 'Number of jobs with earnings $1250/month or less'})\n",
    "job_type_metadata_list.append({'variable': 'SE02', 'description': 'Number of jobs with earnings $1251/month to $3333/month'})\n",
    "job_type_metadata_list.append({'variable': 'SE03', 'description': 'Number of jobs with earnings greater than $3333/month'})\n",
    "job_type_metadata_list.append({'variable': 'SI01', 'description': 'Number of jobs in Goods Producing industry sectors'})\n",
    "job_type_metadata_list.append({'variable': 'SI02', 'description': \"Number of jobs in Trade, Transportation, and Utilities industry sectors\"})\n",
    "job_type_metadata_list.append({'variable': 'SI03', 'description': 'Number of jobs in All Other Services industry sectors'})\n",
    "job_type_metadata_list.append({'variable': 'createdate', 'description': \"Date on which data was created, formatted as YYYYMMDD\"})\n",
    "job_type_metadata_df = pandas.DataFrame.from_records(job_type_metadata_list).set_index('variable')\n",
    "\n",
    "job_segment_metadata_list = []\n",
    "job_segment_metadata_list.append({'variable': 'JT00', 'description': 'All Jobs'})\n",
    "job_segment_metadata_list.append({'variable': 'JT01', 'description': 'Primary Jobs'})\n",
    "job_segment_metadata_list.append({'variable': 'JT02', 'description': 'All Private Jobs'})\n",
    "job_segment_metadata_list.append({'variable': 'JT03', 'description': 'Private Primary Jobs'})\n",
    "job_segment_metadata_list.append({'variable': 'JT04', 'description': 'All Federal Jobs'})\n",
    "job_segment_metadata_list.append({'variable': 'JT05', 'description': 'Federal Primary Jobs'})\n",
    "job_segment_metadata_df = pandas.DataFrame.from_records(job_segment_metadata_list).set_index('variable')\n",
    "\n",
    "\n",
    "raw_jobs_df = raw_jobs_df.merge(job_type_metadata_df.reset_index(), left_on='job_type', right_on='variable').drop('variable', axis=1).rename(columns={'description': 'job_type_description'})\n",
    "raw_jobs_df = raw_jobs_df.merge(job_segment_metadata_df.reset_index(), left_on='job_segment', right_on='variable').drop('variable', axis=1).rename(columns={'description': 'job_segment_description'})\n",
    "\n",
    "print('getting geo crosswalk file...')\n",
    "geo_walk_df = pandas.read_csv('{0:}{1:}/{1:}_xwalk.csv'.format(jobs_dir,state_abbrev), low_memory=False)#, index_col='tabblk2010')\n",
    "\n",
    "print('getting GEOIDs from crosswalk file...')\n",
    "raw_jobs_df = raw_jobs_df.merge(geo_walk_df[['tabblk2020','cty','trct','ctyname','trctname']], how='left', left_on='w_geocode', right_on='tabblk2020')\n",
    "if (thestate <= 9):\n",
    "    raw_jobs_df = raw_jobs_df.assign(GEOID = raw_jobs_df['trct'].apply(lambda x: '14000US0'+str(x)))\n",
    "else:\n",
    "    raw_jobs_df = raw_jobs_df.assign(GEOID = raw_jobs_df['trct'].apply(lambda x: '14000US'+str(x)))\n",
    "\n",
    "print('selecting jobs in city...')\n",
    "raw_jobs_df = raw_jobs_df[raw_jobs_df['GEOID'].isin(data_gdf.index.get_level_values('GEOID').tolist())]\n",
    "\n",
    "print('finding columns to sum over...')\n",
    "jobs_metadata_df = pandas.read_csv(code_lookup_dir+'wac_jobs_metadata.csv', encoding='utf-8', index_col='varnum')\n",
    "\n",
    "jobs_total_columns = jobs_metadata_df[jobs_metadata_df['variable'] == 'C000']['variable'].tolist()\n",
    "jobs_age_columns = jobs_metadata_df[jobs_metadata_df['variable'].apply(lambda x: x[0:2] == 'CA')]['variable'].tolist()\n",
    "jobs_earnings_columns = jobs_metadata_df[jobs_metadata_df['variable'].apply(lambda x: x[0:2] == 'CE')]['variable'].tolist()\n",
    "jobs_sector_columns = jobs_metadata_df[jobs_metadata_df['variable'].apply(lambda x: x[0:3] == 'CNS')]['variable'].tolist()\n",
    "jobs_race_columns = jobs_metadata_df[jobs_metadata_df['variable'].apply(lambda x: x[0:2] == 'CR')]['variable'].tolist()\n",
    "jobs_ethnicity_columns = jobs_metadata_df[jobs_metadata_df['variable'].apply(lambda x: x[0:2] == 'CT')]['variable'].tolist()\n",
    "jobs_education_columns = jobs_metadata_df[jobs_metadata_df['variable'].apply(lambda x: x[0:2] == 'CD')]['variable'].tolist()\n",
    "jobs_sex_columns = jobs_metadata_df[jobs_metadata_df['variable'].apply(lambda x: x[0:2] == 'CS')]['variable'].tolist()\n",
    "jobs_firm_age_columns = jobs_metadata_df[jobs_metadata_df['variable'].apply(lambda x: x[0:3] == 'CFA')]['variable'].tolist()\n",
    "jobs_firm_size_columns = jobs_metadata_df[jobs_metadata_df['variable'].apply(lambda x: x[0:3] == 'CFS')]['variable'].tolist()\n",
    "\n",
    "jobs_columns = jobs_total_columns + jobs_age_columns + jobs_earnings_columns + jobs_sector_columns\n",
    "jobs_columns += jobs_race_columns + jobs_ethnicity_columns + jobs_education_columns + jobs_sex_columns\n",
    "jobs_columns += jobs_firm_age_columns + jobs_firm_size_columns\n",
    "\n",
    "jobs_metadata_df = jobs_metadata_df.set_index('variable')\n",
    "\n",
    "\n",
    "print('summing jobs values into tract-years...')\n",
    "jobs_df = pandas.DataFrame(data=None, columns=jobs_columns, index=pandas.MultiIndex.from_frame(raw_jobs_df[['GEOID','year']].drop_duplicates().sort_values(by=['GEOID','year'])))\n",
    "jobs_df.loc[:,:][jobs_columns] = raw_jobs_df.groupby(['GEOID','year'])[jobs_columns].sum()\n",
    "\n",
    "print('Extrapolating jobs for recent years needed...')\n",
    "latest_job_data_year = np.max(jobs_df.index.get_level_values('year'))\n",
    "for extrapolate_year in range(latest_job_data_year+1, latest_year+1):\n",
    "    print('\\t{0:.0f}'.format(extrapolate_year))\n",
    "    new_jobs_df = jobs_df.xs(latest_job_data_year, level='year').assign(year = extrapolate_year).reset_index().set_index(['GEOID','year'])\n",
    "    jobs_df = pandas.concat((jobs_df, new_jobs_df), axis=0)\n",
    "jobs_df = jobs_df.sort_index()\n",
    "\n",
    "print('calculating total and small business jobs...')\n",
    "jobs_df = jobs_df.assign(total_jobs = jobs_df['C000'])\n",
    "jobs_df = jobs_df.assign(sb_jobs = jobs_df['C000'] - jobs_df['CFS05'])\n",
    "\n",
    "print('joining to loan data...')\n",
    "data_gdf = data_gdf.sort_index().join(jobs_df)\n",
    "\n",
    "\n",
    "e = time.time()\n",
    "g = g + (e-s)\n",
    "print('\\nAdded jobs data for {0:,.0f} tract-years in {1:,.0f} minutes {2:,.0f} seconds!'.format(len(data_gdf), np.floor((e-s)/60), np.floor((e-s)%60)))\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Get loans per job"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Calulating loans per job (total and with firm size 0-19)...\n",
      "recoding zero/zero values to zero and infinite values to NaN...\n",
      "Done in 0.2 seconds!\n"
     ]
    }
   ],
   "source": [
    "s = time.time()\n",
    "\n",
    "print('Calulating loans per job (total and with firm size 0-19)...')\n",
    "\n",
    "data_gdf = data_gdf.assign(nLoans1_per_totaljob = data_gdf['nLoans1'] / data_gdf['total_jobs'])\n",
    "data_gdf = data_gdf.assign(amtLoans1_per_totaljob = data_gdf['amtLoans1'] / data_gdf['total_jobs'])\n",
    "data_gdf = data_gdf.assign(nLoans100k_per_totaljob = data_gdf['nLoans100k'] / data_gdf['total_jobs'])\n",
    "data_gdf = data_gdf.assign(amtLoans100k_per_totaljob = data_gdf['amtLoans100k'] / data_gdf['total_jobs'])\n",
    "data_gdf = data_gdf.assign(nLoans250k_per_totaljob = data_gdf['nLoans250k'] / data_gdf['total_jobs'])\n",
    "data_gdf = data_gdf.assign(amtLoans250k_per_totaljob = data_gdf['amtLoans250k'] / data_gdf['total_jobs'])\n",
    "data_gdf = data_gdf.assign(nLoansToSmallest_per_totaljob = data_gdf['nLoansToSmallest'] / data_gdf['total_jobs'])\n",
    "data_gdf = data_gdf.assign(amtLoansToSmallest_per_totaljob = data_gdf['amtLoansToSmallest'] / data_gdf['total_jobs'])\n",
    "data_gdf = data_gdf.assign(nLoans_per_totaljob = data_gdf['nLoans'] / data_gdf['total_jobs'])\n",
    "data_gdf = data_gdf.assign(amtLoans_per_totaljob = data_gdf['amtLoans'] / data_gdf['total_jobs'])\n",
    "data_gdf = data_gdf.assign(nWorkingLoans_per_totaljob = data_gdf['nWorkingLoans'] / data_gdf['total_jobs'])\n",
    "data_gdf = data_gdf.assign(amtWorkingLoans_per_totaljob = data_gdf['amtWorkingLoans'] / data_gdf['total_jobs'])\n",
    "\n",
    "data_gdf = data_gdf.assign(nLoans1_per_sbjob = data_gdf['nLoans1'] / data_gdf['sb_jobs'])\n",
    "data_gdf = data_gdf.assign(amtLoans1_per_sbjob = data_gdf['amtLoans1'] / data_gdf['sb_jobs'])\n",
    "data_gdf = data_gdf.assign(nLoans100k_per_sbjob = data_gdf['nLoans100k'] / data_gdf['sb_jobs'])\n",
    "data_gdf = data_gdf.assign(amtLoans100k_per_sbjob = data_gdf['amtLoans100k'] / data_gdf['sb_jobs'])\n",
    "data_gdf = data_gdf.assign(nLoans250k_per_sbjob = data_gdf['nLoans250k'] / data_gdf['sb_jobs'])\n",
    "data_gdf = data_gdf.assign(amtLoans250k_per_sbjob = data_gdf['amtLoans250k'] / data_gdf['sb_jobs'])\n",
    "data_gdf = data_gdf.assign(nLoansToSmallest_per_sbjob = data_gdf['nLoansToSmallest'] / data_gdf['sb_jobs'])\n",
    "data_gdf = data_gdf.assign(amtLoansToSmallest_per_sbjob = data_gdf['amtLoansToSmallest'] / data_gdf['sb_jobs'])\n",
    "data_gdf = data_gdf.assign(nLoans_per_sbjob = data_gdf['nLoans'] / data_gdf['sb_jobs'])\n",
    "data_gdf = data_gdf.assign(amtLoans_per_sbjob = data_gdf['amtLoans'] / data_gdf['sb_jobs'])\n",
    "data_gdf = data_gdf.assign(nWorkingLoans_per_sbjob = data_gdf['nWorkingLoans'] / data_gdf['sb_jobs'])\n",
    "data_gdf = data_gdf.assign(amtWorkingLoans_per_sbjob = data_gdf['amtWorkingLoans'] / data_gdf['sb_jobs'])\n",
    "\n",
    "\n",
    "per_job_columns = ['nLoans1_per_totaljob', 'amtLoans1_per_totaljob', 'nLoans100k_per_totaljob']\n",
    "per_job_columns += ['amtLoans100k_per_totaljob', 'nLoans250k_per_totaljob', 'amtLoans250k_per_totaljob']\n",
    "per_job_columns += ['nLoansToSmallest_per_totaljob', 'amtLoansToSmallest_per_totaljob']\n",
    "per_job_columns += ['nLoans_per_totaljob', 'amtLoans_per_totaljob', 'nWorkingLoans_per_totaljob']\n",
    "per_job_columns += ['amtWorkingLoans_per_totaljob', 'nLoans1_per_sbjob', 'amtLoans1_per_sbjob']\n",
    "per_job_columns += ['nLoans100k_per_sbjob', 'amtLoans100k_per_sbjob', 'nLoans250k_per_sbjob']\n",
    "per_job_columns += ['amtLoans250k_per_sbjob', 'nLoansToSmallest_per_sbjob', 'amtLoansToSmallest_per_sbjob']\n",
    "per_job_columns += ['nLoans_per_sbjob', 'amtLoans_per_sbjob', 'nWorkingLoans_per_sbjob']\n",
    "per_job_columns += ['amtWorkingLoans_per_sbjob']\n",
    "\n",
    "\n",
    "print('recoding zero/zero values to zero and infinite values to NaN...')\n",
    "for x in data_gdf[per_job_columns]:\n",
    "    data_gdf.loc[:, x] = data_gdf[x].fillna(0)\n",
    "    data_gdf.loc[data_gdf[x] == np.inf, x] = np.nan\n",
    "\n",
    "# print('backing up...')\n",
    "# data_gdf_bk = data_gdf\n",
    "\n",
    "e = time.time()\n",
    "g = g + (e-s)\n",
    "print('Done in {0:,.1f} seconds!'.format(e-s)) \n",
    "# data_gdf[(data_gdf['nWorkingLoans'] == 0) | (data_gdf['CFS01'] == 0)][\n",
    "#     ['nWorkingLoans', 'amtWorkingLoans', 'total_jobs', 'sb_jobs', 'nWorkingLoans_per_totaljob', 'amtWorkingLoans_per_totaljob', 'nWorkingLoans_per_sbjob', 'amtWorkingLoans_per_sbjob']\n",
    "# ].sample(4).T\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Get ACS data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Getting ACS 5-year census data...\n",
      "\t2022 estimates...\n",
      "\t2022 margins of error...\n",
      "\t2021 estimates...\n",
      "\t2021 margins of error...\n",
      "\t2020 estimates...\n",
      "\t2020 margins of error...\n",
      "\t2019 estimates...\n",
      "\t2019 margins of error...\n",
      "\t2018 estimates...\n",
      "\t2018 margins of error...\n",
      "\t2017 estimates...\n",
      "\t2017 margins of error...\n",
      "\t2016 estimates...\n",
      "\t2016 margins of error...\n",
      "\t2015 estimates...\n",
      "\t2015 margins of error...\n",
      "\t2014 estimates...\n",
      "\t2014 margins of error...\n",
      "\t2013 estimates...\n",
      "\t2013 margins of error...\n",
      "\t2012 estimates...\n",
      "\t2012 margins of error...\n",
      "setting index...\n",
      "renaming columns...\n",
      "nan-ifying missing value mask values in margins of error...\n",
      "\tFixing column 3 of 105: B01001_001_err...\n",
      "\tFixing column 4 of 105: B02001_001_err...\n",
      "\tFixing column 7 of 105: B08013_001_err...\n",
      "\tFixing column 14 of 105: B01001_011_err...\n",
      "\tFixing column 15 of 105: B01001_012_err...\n",
      "\tFixing column 18 of 105: B01001_015_err...\n",
      "\tFixing column 19 of 105: B01001_016_err...\n",
      "\tFixing column 29 of 105: B01001_026_err...\n",
      "\tFixing column 30 of 105: B01001_027_err...\n",
      "\tFixing column 33 of 105: B01001_030_err...\n",
      "\tFixing column 34 of 105: B01001_031_err...\n",
      "\tFixing column 38 of 105: B01001_035_err...\n",
      "\tFixing column 39 of 105: B01001_036_err...\n",
      "\tFixing column 42 of 105: B01001_039_err...\n",
      "\tFixing column 43 of 105: B01001_040_err...\n",
      "\tFixing column 53 of 105: B15003_001_err...\n",
      "\tFixing column 78 of 105: B17001_001_err...\n",
      "\tFixing column 79 of 105: B17001_002_err...\n",
      "\tFixing column 80 of 105: B17001A_001_err...\n",
      "\tFixing column 81 of 105: B17001A_002_err...\n",
      "\tFixing column 82 of 105: B17001B_001_err...\n",
      "\tFixing column 83 of 105: B17001B_002_err...\n",
      "\tFixing column 84 of 105: B19013_001_err...\n",
      "\tFixing column 85 of 105: B19013A_001_err...\n",
      "\tFixing column 86 of 105: B19013B_001_err...\n",
      "\tFixing column 87 of 105: B19113_001_err...\n",
      "\tFixing column 88 of 105: B23025_002_err...\n",
      "\tFixing column 89 of 105: B23025_005_err...\n",
      "\tFixing column 102 of 105: B25077_001_err...\n",
      "\tFixing column 104 of 105: B25035_001_err...\n",
      "\n",
      "\n",
      "\tmissing value mask values were [-666666666, -555555555, -333333333, -222222222]; now nan\n",
      "joining to rest of data...\n",
      "\tdeleting duplicate state and county columns from ACS5 data...\n",
      "backing up...\n",
      "Pittsburgh: Kept 1,480 tract-years in 6 minutes 4 seconds!\n"
     ]
    }
   ],
   "source": [
    "s = time.time()\n",
    "\n",
    "acs5_estimates_df = pandas.DataFrame()\n",
    "acs5_margins_of_error_df = pandas.DataFrame()\n",
    "print('Getting ACS 5-year census data...')\n",
    "\n",
    "\n",
    "for thisyear in years:\n",
    "    print('\\t{0:.0f} estimates...'.format(thisyear))\n",
    "    acs5_estimates_this_year_df = pandas.read_csv(census_dir+'{0:.0f}/data/estimates.csv'.format(thisyear), low_memory=False, encoding='utf-8')\n",
    "    acs5_estimates_this_year_df = acs5_estimates_this_year_df.drop([x for x in acs5_estimates_this_year_df.columns.tolist() if 'unnamed' in x.lower()], axis=1)\n",
    "\n",
    "    acs5_estimates_this_year_df = acs5_estimates_this_year_df.assign(year = thisyear)\n",
    "    acs5_estimates_df = pandas.concat((acs5_estimates_df, acs5_estimates_this_year_df), axis=0, sort=False)\n",
    "        \n",
    "    print('\\t{0:.0f} margins of error...'.format(thisyear))\n",
    "    acs5_margins_of_error_this_year_df = pandas.read_csv(census_dir+'{0:.0f}/data/margins_of_error.csv'.format(thisyear), low_memory=False, encoding='utf-8')\n",
    "    if (thisyear <= 2014):\n",
    "        acs5_margins_of_error_this_year_df = acs5_margins_of_error_this_year_df.drop([x for x in acs5_margins_of_error_this_year_df.columns.tolist() if 'unnamed' in x.lower()], axis=1)\n",
    "    acs5_margins_of_error_this_year_df = acs5_margins_of_error_this_year_df.assign(year = thisyear)\n",
    "    acs5_margins_of_error_df = pandas.concat((acs5_margins_of_error_df, acs5_margins_of_error_this_year_df), axis=0, sort=False)\n",
    "\n",
    "print('setting index...')\n",
    "acs5_estimates_df = acs5_estimates_df.set_index(['GEOID', 'year'])\n",
    "acs5_margins_of_error_df = acs5_margins_of_error_df.set_index(['GEOID', 'year'])\n",
    "\n",
    "print('renaming columns...')\n",
    "est_oldcols = acs5_estimates_df.columns[3:]\n",
    "est_newcols = [x[:-1] for x in est_oldcols]\n",
    "est_coldict = dict(zip(est_oldcols, est_newcols))\n",
    "acs5_estimates_df = acs5_estimates_df.rename(columns=est_coldict)\n",
    "\n",
    "merr_oldcols = acs5_margins_of_error_df.columns[3:]\n",
    "merr_newcols = [x[:-1]+'_err' for x in merr_oldcols]\n",
    "merr_coldict = dict(zip(merr_oldcols, merr_newcols))\n",
    "acs5_margins_of_error_df = acs5_margins_of_error_df.rename(columns=merr_coldict)\n",
    "\n",
    "\n",
    "print('nan-ifying missing value mask values in margins of error...')\n",
    "missing_val_flag_values = []\n",
    "for i in range(0, len(acs5_margins_of_error_df.columns)):\n",
    "    thiscol = acs5_margins_of_error_df.columns[i]\n",
    "    if (acs5_margins_of_error_df[thiscol].min() < 0):\n",
    "        print('\\tFixing column {0:.0f} of {1:.0f}: {2:}...'.format(i, len(acs5_margins_of_error_df.columns), thiscol))\n",
    "        missing_val_flag_values += [x for x in acs5_margins_of_error_df[acs5_margins_of_error_df[thiscol] < 0][thiscol].drop_duplicates().tolist() if x not in missing_val_flag_values]\n",
    "        acs5_margins_of_error_df.loc[acs5_margins_of_error_df[thiscol] < 0, thiscol] = np.nan\n",
    "missing_val_flag_values = sorted(missing_val_flag_values)\n",
    "print('\\n')\n",
    "print('\\tmissing value mask values were {0:}; now nan'.format(missing_val_flag_values))\n",
    "\n",
    "print('joining to rest of data...')\n",
    "data_gdf = data_gdf.join(acs5_estimates_df, how='left', lsuffix='_data', rsuffix='_est')\n",
    "data_gdf = data_gdf.rename(columns={'state_data': 'state', 'county_data': 'county'})\n",
    "data_gdf = data_gdf.join(acs5_margins_of_error_df, how='left', lsuffix='_data', rsuffix='_merr')\n",
    "data_gdf = data_gdf.rename(columns={'state_data': 'state', 'county_data': 'county'})\n",
    "print('\\tdeleting duplicate state and county columns from ACS5 data...')\n",
    "data_gdf = data_gdf.drop(['state_est', 'state_merr', 'county_est', 'county_merr'], axis=1)\n",
    "\n",
    "\n",
    "\n",
    "print('backing up...')\n",
    "data_gdf_bk = data_gdf\n",
    "\n",
    "e = time.time()\n",
    "g = g + (e-s)\n",
    "\n",
    "#print('{0:}: Kept {1:,.0f} tract-years in {2:,.0f} minutes {3:.0f} seconds!'.format(city, len(data_gdf), np.floor((e-s)/60), np.floor((e-s)%60)))\n",
    "print('{0:}: Kept {1:,.0f} tract-years in {2:,.0f} minutes {3:.0f} seconds!'.format(city, len(data_gdf), np.floor(g/60), np.floor(g%60)))\n",
    "\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Cross-check ACS columns\n",
    "\n",
    "What changed from year to year?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Not cross-checking all ACS5 columns\n"
     ]
    }
   ],
   "source": [
    "\n",
    "#cross_check_acs5_year_columns = False\n",
    "# earlier = 2020\n",
    "# later = 2021\n",
    "# print('Getting metadata for all ACS5 columns...')\n",
    "# acs5_columns = [x for x in data_gdf.columns.tolist() if (x[0] == 'B') and (x[-4:] != '_err')]\n",
    "\n",
    "# print('\\t{0:.0f}...'.format(earlier))\n",
    "# acs5_metadata_earlier_df = pandas.read_csv(data_dir+'acs5/{0:.0f}/variables.csv'.format(earlier), low_memory=False, encoding='utf-8')\n",
    "# # acs5_metadata_earlier_df.loc[:, 'Unique ID'] = acs5_metadata_earlier_df['Unique ID'].apply(lambda x: '{0:}_{1:}'.format(x[:-3], x[-3:]))\n",
    "\n",
    "# print('\\t{0:.0f}...'.format(later))\n",
    "# acs5_metadata_later_df = pandas.read_csv(data_dir+'acs5/{0:.0f}/variables.csv'.format(later), low_memory=False, encoding='utf-8')\n",
    "\n",
    "# print('Joining to compare...')\n",
    "# df = acs5_metadata_earlier_df.merge(acs5_metadata_later_df, how='outer', on='Unique ID', suffixes=(' {0:.0f}'.format(earlier), ' {0:.0f}'.format(later)))\n",
    "\n",
    "# print('\\n')\n",
    "# print('In both {0:.0f} and {1:.0f}: {2:,.0f}'.format(earlier, later, len(df[(df['Table ID {0:.0f}'.format(earlier)].notnull()) & (df['Table ID {0:.0f}'.format(later)].notnull())])))\n",
    "# print('Removed from {0:.0f} and not found in {1:.0f}: {2:,.0f}'.format(earlier, later, len(df[(df['Table ID {0:.0f}'.format(earlier)].notnull()) & (df['Table ID {0:.0f}'.format(later)].isnull())])))\n",
    "# print('Added in {1:.0f} was not found in {0:.0f}: {2:,.0f}'.format(earlier, later, len(df[(df['Table ID {0:.0f}'.format(earlier)].isnull()) & (df['Table ID {0:.0f}'.format(later)].notnull())])))\n",
    "\n",
    "\n",
    "# df[\n",
    "#     (\n",
    "#         (df['Table ID {0:.0f}'.format(earlier)].isnull()) \n",
    "#         & (df['Table ID {0:.0f}'.format(later)].notnull())\n",
    "#         #& (df['Unique ID'].apply(lambda x: '.5' not in x))\n",
    "#      )\n",
    "#    & \n",
    "#      (\n",
    "#          (df['Table ID {0:.0f}'.format(later)].apply(lambda x: 'B01002' not in str(x)))\n",
    "#          & (df['Table ID {0:.0f}'.format(later)].apply(lambda x: 'B17' not in str(x)))\n",
    "#      )\n",
    "# ].sort_values(by='Unique ID')\n",
    "\n",
    "# #df[df['Unique ID'].apply(lambda x: 'B23' in str(x))]\n",
    "\n",
    "if (not(cross_check_acs5_year_columns)):\n",
    "    print('Not cross-checking all ACS5 columns')\n",
    "else:\n",
    "    print('CHANGES from 2012 to 2013:')\n",
    "    print('\\tremoved:')\n",
    "    print('\\t\\tunweighted counts and response rates (tables starting with B98)...')   #df[(df['Label 2012'].isnull()) & (df['Label {0:.0f}'.format(later)].notnull())]['Table ID {0:.0f}'.format(later)].drop_duplicates()\n",
    "    print('\\tadded:')\n",
    "    print('\\t\\tDetailed race (B02003, B02014-19)')   #df[((df['Label 2012'].isnull()) & (df['Label {0:.0f}'.format(later)].notnull())) & (df['Table ID {0:.0f}'.format(later)].apply(lambda x: 'B020' in str(x)))]\n",
    "    print('\\t\\tPlace of birth excluding born at sea (B05006161)')   #df[df['Table ID {0:.0f}'.format(later)].apply(lambda x: ('B0500' in str(x)) & (re.search('\\d',str(x)[-1]) != None))][['Table ID {0:.0f}'.format(later), 'Title {0:.0f}'.format(later), 'Universe {0:.0f}'.format(later)]].drop_duplicates()\n",
    "    print('\\t\\tSex by educational attainment by race (B15002A-I)')   #df[(df['Label 2012'].isnull()) & (df['Label {0:.0f}'.format(later)].notnull()) & df['Table ID {0:.0f}'.format(later)].apply(lambda x: 'B15' in str(x))]\n",
    "    print('\\t\\tDetailed field of bachelors degree (B15010-12)')   #df[(df['Label 2012'].isnull()) & (df['Label {0:.0f}'.format(later)].notnull()) & df['Table ID {0:.0f}'.format(later)].apply(lambda x: 'B15' in str(x))]\n",
    "    print('\\t\\tRatio of income to poverty level (B17002) and poverty by age (B17020 main level)')   #df[(df['Label 2012'].isnull()) & (df['Label {0:.0f}'.format(later)].notnull()) & df['Table ID {0:.0f}'.format(later)].apply(lambda x: 'B17' in str(x))]\n",
    "    print('\\t\\tDisability by age, poverty, and employment (B181080-31)')   #df[(df['Label 2012'].isnull()) & (df['Label {0:.0f}'.format(later)].notnull()) & df['Table ID {0:.0f}'.format(later)].apply(lambda x: 'B18' in str(x))]\n",
    "    print('\\t\\tVeteran status by race+disability (B21001A-I) and poverty (B21007)') \n",
    "    print('\\t\\tSex by employment by race+disability (B23002A-I) and with hours (B23023)') \n",
    "    print('\\t\\tSex by occupation (B2401-07)')\n",
    "    print('\\t\\tHousehold owner costs as percentage of income (B25095)')\n",
    "    print('\\t\\tHealth insurance by race; health insurance source (B27001A-I)')\n",
    "    print('\\t\\tSource of health insurance (B27004-18)')\n",
    "    print('\\t\\tComputers (B28)')\n",
    "    print('\\t\\tImputed fields - compensate for nonresponse bias (B99)')\n",
    "    print('\\t\\tComparisons (C)')\n",
    "\n",
    "    print('-------')\n",
    "\n",
    "    print('CHANGES from 2013 to 2014:')\n",
    "    print('\\trenamed:')\n",
    "    print('\\t\\tmedian family income: was B19113 001, now B19113001...')\n",
    "    print('\\tremoved:')\n",
    "    print('\\t\\tfirst ancestry reported (B04)...')\n",
    "    print('\\t\\tplace of birth, foreign-born, Puerto Rico (B05006PR)...')\n",
    "    print('\\t\\twork status and disability (B23004-23)...')\n",
    "    print('\\t\\tsex by occupation (B24)...')\n",
    "    print('\\t\\thousehold rooms and rent (B25002-09)...')\n",
    "    print('\\t\\tcomparisions for all above, and also...')\n",
    "    print('\\t\\t\\tcomparisons of place of birth, geo mobility (C06, C07)...')\n",
    "    print('\\t\\t\\tworkers by sex and place (C08)...')\n",
    "    print('\\t\\t\\twomen who had a birth (C13)...')\n",
    "    print('\\tadded:')\n",
    "    print('\\t\\tsex by age for foreign-born population (B05013-14)...')\n",
    "    print('\\t\\thours worked per week, full-time status (B23026-27)...')\n",
    "    print('\\t\\thousehold income by household costs (B2509)...')\n",
    "    print('\\t\\thandheld computers (vars B28001005-009)')\n",
    "    print('\\t\\tunweighted counts and imputations (B98-99)...')\n",
    "\n",
    "    print('-------')\n",
    "\n",
    "    print('CHANGES from 2014 to 2015:')\n",
    "    #print('\\trenamed:')\n",
    "    #print('\\t\\tmedian family income: was B19113 001, now B19113001...')\n",
    "    print('\\tremoved:')\n",
    "    print('\\t\\tAmerican Indian, Asian, Pacific Islander by subgroup (B02005-07)...')\n",
    "    print('\\t\\tGrandparents living with grandchildren by age (B10060-61)...')\n",
    "    print('\\t\\tHousehold type by tenure (B11012)...')\n",
    "    print('\\t\\tMarital status by nativity (B12005)...')\n",
    "    print('\\t\\tMarital status by age for women (B13001)...')\n",
    "    print('\\t\\tYear moved into house (B25129)...')\n",
    "    print('\\tadded:')\n",
    "    print('\\t\\tMore options for place of birth by nativity (B05002_016 to 027)...')\n",
    "    print('\\t\\tPlace of birth: other northern America, Puerto Rico (B05006_162 and B05006PR)...')\n",
    "    print('\\t\\tMore options for place of birth by year of entry (B05007_095 to 106)...')\n",
    "    print('\\t\\tLiving arrangements of adults (B09021)...')\n",
    "    print('\\t\\tMore occupation categories for sex by occupation (B24010/20/60)...')\n",
    "    print('\\t\\tComparison industry by occupation (C24050/60)...')\n",
    "    print('\\t\\tRent (B25031/56/61/63/74)...')\n",
    "    print('\\t\\tYear structure built (B25034/36/111/112/126/127)...')\n",
    "    print('\\t\\tHome value (B25075/85)...')\n",
    "    print('\\t\\tHomeowner costs (B25087/94/95/96/104/107/108)...')\n",
    "\n",
    "    print('-------')\n",
    "\n",
    "    print('CHANGES from 2015 to 2016:')\n",
    "    print('\\tremoved:')\n",
    "    print('\\t\\tInternet subscriptions (B28002/B99283)...')\n",
    "    print('\\tadded:')\n",
    "    print('\\t\\tMore options for language spoken at home (B16001_120-128 and B16002_015-038)...')\n",
    "    print('\\t\\tComparisons on language spoken at home (C16001-02)...')\n",
    "    print('\\t\\tOptions for no computer in household (B28001_010-11)...')\n",
    "    print('\\t\\tTypes of computers in household (B28010)...')\n",
    "    print('\\t\\tInternet subscriptions in household (B28011)...')\n",
    "    print('\\t\\tHousing unit nonresponse rate with reasons (B98021_010)...')\n",
    "    print('\\t\\tImputed computer type (B99282_008-009)...')\n",
    "\n",
    "\n",
    "    print('-------')\n",
    "\n",
    "    print('CHANGES from 2016 to 2017:')\n",
    "    print('\\tremoved:')\n",
    "    print('\\t\\tsome options for health insurance by type (B27011_018-35 / B27012_047-98 / B27013_023-46 / B27014_023-46)...')\n",
    "    print('\\tremoved:')\n",
    "    print('\\t\\tearners in family (B19122)...')\n",
    "    print('\\t\\tgroup quarters type (B26101-220)...')\n",
    "    print('\\t\\tcomparisons for health insurance coverage (C27)...')\n",
    "    print('\\t\\tvoting age citizens (B29)...')\n",
    "\n",
    "\n",
    "    print('-------')\n",
    "\n",
    "    print('CHANGES from 2017 to 2018:')\n",
    "\n",
    "    print('\\tremoved:')\n",
    "    print('Option for naturalized before 1980 (B05011_011)...')\n",
    "    print('Options for group quarters by sex and age (B26101_215-232 / B26108_083-085 / B26201_254-298 / B26206_042-46 / B26208_104-109)...')\n",
    "\n",
    "    print('\\tadded:')\n",
    "    print('\\t\\tfamily size by cash assistance (B19123)...')\n",
    "    print('\\t\\tdetailed occupation (B24114-16)...')\n",
    "    print('\\t\\tmore choices for detailed occupation (B24121/22/23/24/25/26_527-566)...')\n",
    "    print('\\t\\tOptions for group quarters (B26114_023-35 / B26214_030-44)')\n",
    "\n",
    "\n",
    "    print('-------')\n",
    "\n",
    "    print('CHANGES from 2018 to 2019:')\n",
    "    print('\\tremoved:')\n",
    "    print('\\t\\tUnweighted sample count (B00001/02)...')\n",
    "    print('\\t\\tChildren in nonfamily households (B09005_006)...')\n",
    "    print('\\t\\tPresence of unmarried partner (B09008)...')\n",
    "    print('\\t\\tOptions for living situation (B09019_027-38)...')\n",
    "    print('\\t\\tPartner in unmarried household (B13004_010-11...')\n",
    "\n",
    "    print('\\tadded:')\n",
    "    print('\\t\\tOptions for foreign-born place of birth (B05006_163-68)...')\n",
    "    print('\\t\\tPlace of birth by year of entry (B05015)...')\n",
    "    print('\\t\\tCohabiting couples (B11008/09/12)...')\n",
    "    print('\\t\\tOptions for median duration of current marriage (B12504_010-13)...')\n",
    "    print('\\t\\tFinal interview count (B98003_001)...')\n",
    "\n",
    "    print('-------')\n",
    "\n",
    "    print('CHANGES from 2019 to 2020:')\n",
    "    print('\\tremoved: NONE')\n",
    "    print('\\tadded:')\n",
    "    print('\\tdetailed industry (B24134-6)...')\n",
    "    print('\\thousing units without a mortgage (B25081_009)...')\n",
    "    print('\\thousing unit vacancies (B25130-31)...')\n",
    "    print('\\t\\tage by computer ownership (B28012_021)...')\n",
    "    print('\\t\\tallocation of vacancy (B992524-5)...')\n",
    "\n",
    "\n",
    "    print('-------')\n",
    "\n",
    "    print('CHANGES from 2020 to 2021:')\n",
    "    print('removed:')\n",
    "    print('\\tDetailed race (B02003)...')\n",
    "    print('\\tEducational attainment by race (B15002A-I / B15010)...')\n",
    "    print('\\tDetailed language (B16002-3)...')\n",
    "    print('\\tIncome/poverty ratio (B17002)...')\n",
    "    print('\\tNumber of disabilies / employment by disability (B18108/12)...')\n",
    "    print('\\tSex by age by employment/occupation by race (B23-4A-I)...')\n",
    "    print('\\tReasons for vacanting houses (B25)...')\n",
    "    print('\\tHealth insurance by age by race (B27A-I)...')\n",
    "    print('\\tVacancy allocation (B9925I)...')\n",
    "    print('\\tComparison for sex by age (C01001)...')\n",
    "    print('\\tComparison for American Indian / Asian / Pacific islander subgroup (C02014-6)...')\n",
    "    print('\\tComparison for ancsetry/nativity/birth/mobility (C03-7)...')\n",
    "    print('\\tComparison for transportation (C08-11)...')\n",
    "    print('\\tComparison for education/language/poverty/disability (C14-18/21/30-1)...')\n",
    "    print('\\tComparison for employment and disability (C21)...')\n",
    "    print('\\tComparison for public assistance (C22)...')\n",
    "    print('\\tComparison for ex by age by employment/occupation by race (C23-24)...')\n",
    "    print('\\tComparison for reasons for vacanting houses (C25)...')\n",
    "    print('\\tComparison for health insurance by age by race (C27)...')\n",
    "    print('added:')\n",
    "    print('\\tPoverty status by housing unit (B17)...')\n",
    "    print('\\tElectricity costs (B25)...')\n",
    "    print('\\tMedian age by sex (B01002 plus A-I)...')\n",
    "    print('\\tPacific Islander by category (B02016/19)...')\n",
    "    print('\\tMedian age (B06002/07002)...')\n",
    "    print('\\tMoves within USA (B07101)...')\n",
    "    print('\\tMedian age by transportation (B08103/B08503)...')\n",
    "    print('\\tMFI for grandparent households (B10010)...')\n",
    "    print('\\tMedian age at first marriage (B12007)...')\n",
    "\n",
    "    print('\\tTotal fields of bachelors degree (B15012)...')\n",
    "    print('\\tQuintiles of household income (B19080-2)...')\n",
    "    print('\\tMedian age (B23/26)')\n",
    "    print('\\tCoverage rates (B98-9)')\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Count missing values in estimates"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Counting missing values in Pittsburgh...\n",
      "2012\n",
      "-------------------\n",
      "-------------------\n",
      "\n",
      "\n",
      "2013\n",
      "-------------------\n",
      "-------------------\n",
      "\n",
      "\n",
      "2014\n",
      "-------------------\n",
      "-------------------\n",
      "\n",
      "\n",
      "2015\n",
      "-------------------\n",
      "B08013_001: 11 missing values...\n",
      "B19013_001: 13 missing values...\n",
      "B19013A_001: 40 missing values...\n",
      "B19013B_001: 120 missing values...\n",
      "B19113_001: 18 missing values...\n",
      "B25035_001: 10 missing values...\n",
      "B25077_001: 25 missing values...\n",
      "-------------------\n",
      "\n",
      "\n",
      "2016\n",
      "-------------------\n",
      "B08013_001: 13 missing values...\n",
      "B19013_001: 11 missing values...\n",
      "B19013A_001: 21 missing values...\n",
      "B19013B_001: 52 missing values...\n",
      "B19113_001: 15 missing values...\n",
      "B25035_001: 10 missing values...\n",
      "B25077_001: 15 missing values...\n",
      "-------------------\n",
      "\n",
      "\n",
      "2017\n",
      "-------------------\n",
      "B08013_001: 17 missing values...\n",
      "B19013_001: 12 missing values...\n",
      "B19013A_001: 22 missing values...\n",
      "B19013B_001: 58 missing values...\n",
      "B19113_001: 16 missing values...\n",
      "B25035_001: 12 missing values...\n",
      "B25077_001: 17 missing values...\n",
      "-------------------\n",
      "\n",
      "\n",
      "2018\n",
      "-------------------\n",
      "B08013_001: 17 missing values...\n",
      "B19013_001: 13 missing values...\n",
      "B19013A_001: 26 missing values...\n",
      "B19013B_001: 54 missing values...\n",
      "B19113_001: 15 missing values...\n",
      "B25035_001: 13 missing values...\n",
      "B25077_001: 17 missing values...\n",
      "-------------------\n",
      "\n",
      "\n",
      "2019\n",
      "-------------------\n",
      "B08013_001: 22 missing values...\n",
      "B19013_001: 14 missing values...\n",
      "B19013A_001: 25 missing values...\n",
      "B19013B_001: 52 missing values...\n",
      "B19113_001: 19 missing values...\n",
      "B25035_001: 12 missing values...\n",
      "B25077_001: 17 missing values...\n",
      "-------------------\n",
      "\n",
      "\n",
      "2020\n",
      "-------------------\n",
      "B08013_001: 17 missing values...\n",
      "B19013_001: 14 missing values...\n",
      "B19013A_001: 22 missing values...\n",
      "B19013B_001: 53 missing values...\n",
      "B19113_001: 18 missing values...\n",
      "B25035_001: 13 missing values...\n",
      "B25077_001: 19 missing values...\n",
      "-------------------\n",
      "\n",
      "\n",
      "2021\n",
      "-------------------\n",
      "B08013_001: 30 missing values...\n",
      "B19013_001: 14 missing values...\n",
      "B19013A_001: 22 missing values...\n",
      "B19013B_001: 55 missing values...\n",
      "B19113_001: 20 missing values...\n",
      "B25035_001: 13 missing values...\n",
      "B25077_001: 20 missing values...\n",
      "-------------------\n",
      "\n",
      "\n",
      "2022\n",
      "-------------------\n",
      "B08013_001: 28 missing values...\n",
      "B19013_001: 13 missing values...\n",
      "B19013A_001: 20 missing values...\n",
      "B19013B_001: 56 missing values...\n",
      "B19113_001: 20 missing values...\n",
      "B25035_001: 12 missing values...\n",
      "B25077_001: 19 missing values...\n",
      "-------------------\n",
      "\n",
      "\n",
      "'k'\n"
     ]
    }
   ],
   "source": [
    "print('Counting missing values in {0:}...'.format(city))\n",
    "this_city_geoids = data_gdf.index.get_level_values(level='GEOID').drop_duplicates().tolist()\n",
    "\n",
    "for oneyear in years[-1::-1]:\n",
    "    #print(data_gdf.xs(oneyear, level='year').columns.tolist())\n",
    "    acs5_thisyear_metadata_df = pandas.read_csv('{0:}{1:.0f}/metadata/variables.csv'.format(census_dir, oneyear))\n",
    "    if (oneyear <= 2014):\n",
    "        acs5_thisyear_metadata_df.loc[:, 'Unique ID'] = acs5_thisyear_metadata_df['Unique ID'].apply(lambda x: '{0:}_{1:}'.format(x[:-3], x[-3:]))\n",
    "    #print(acs5_thisyear_metadata_df[acs5_thisyear_metadata_df['Unique ID'].apply(lambda x: '_' not in x)])\n",
    "    acs5_thisyear_metadata_df = acs5_thisyear_metadata_df.set_index('Unique ID')\n",
    "    print(oneyear)\n",
    "    print('-------------------')    \n",
    "    acs5_thisyear_columns = acs5_thisyear_metadata_df[acs5_thisyear_metadata_df.index.isin(data_gdf.xs(oneyear, level='year').columns.tolist())].index.tolist()\n",
    "    #print(acs5_thisyear_columns)\n",
    "    \n",
    "    \n",
    "    for thiscol in acs5_thisyear_columns:\n",
    "        \n",
    "    \n",
    "        if (len(data_gdf.xs(oneyear, level='year')[data_gdf.xs(oneyear, level='year')[thiscol].isnull()]) > 0):\n",
    "            print('{0:}: {1:,.0f} missing values...'.format(thiscol, len(\n",
    "                data_gdf.xs(oneyear, level='year')[(data_gdf.xs(oneyear, level='year')[thiscol].isnull()) & (data_gdf.xs(oneyear, level='year').index.isin(this_city_geoids))]\n",
    "            )))\n",
    "            \n",
    "            \n",
    "    print('-------------------')\n",
    "    print('\\n')\n",
    "#     for thiscol in data_gdf.xs(oneyear, level='year').columns.tolist():\n",
    "#         print(thiscol)\n",
    "    \n",
    "    \n",
    "#     print(oneyear)\n",
    "\n",
    "#     print(acs5_thisyear_metadata_df[acs5_thisyear_metadata_df.index.isin(data_gdf.xs(oneyear, level='year').columns.tolist())])\n",
    "#     print('-------------------')\n",
    "#     print('\\n')\n",
    "    \n",
    "    \n",
    "    \n",
    "#data_gdf.columns.tolist()\n",
    "pprint('k')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Calculate composite demographic columns"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "calculating and renaming estimates columns for IVs...\n",
      "population 25 and older with bachelors degree or higher (B15003_022 t.e.m. B15003_025...\n",
      "...population 25 years and older...\n",
      "...householder sex & race, unempoyment, poverty, household type, home value, home age, travel time...\n",
      "...race, owner-occupied units, mfi, vacants...\n",
      "....comparison variables: total population, total households, poverty status, rentage...\n",
      "MFI & median home value: substituting \".\" with np.nan, converting to numeric...\n",
      "median home value: substituting \".\" with np.nan, converting to numeric...\n",
      "Done in 0 seconds!\n"
     ]
    }
   ],
   "source": [
    "s = time.time()\n",
    "# print('getting from backup...')\n",
    "# data_gdf = data_gdf_bk\n",
    "\n",
    "print('\\ncalculating and renaming estimates columns for IVs...')\n",
    "\n",
    "if (debug >= 1):\n",
    "    print('population 25 and older with bachelors degree or higher (B15003_022 t.e.m. B15003_025...')\n",
    "h = data_gdf['B15003_022'] + data_gdf['B15003_023'] + data_gdf['B15003_024']  + data_gdf['B15003_025']\n",
    "data_gdf = data_gdf.assign(educated = pandas.to_numeric(h, errors='coerce'))\n",
    "\n",
    "if (debug >= 1):\n",
    "    print('...population 25 years and older...')\n",
    "data_gdf = data_gdf.rename(columns={'B15003_001': 'pop_25plus'})\n",
    "\n",
    "\n",
    "if (debug >= 1):\n",
    "    print('...householder sex & race, unempoyment, poverty, household type, home value, home age, travel time...')\n",
    "data_gdf = data_gdf.rename(columns = {  \n",
    "    'B11001_006': 'female_householders',\n",
    "    'B11001A_001': 'white_householders',\n",
    "    'B11001B_001': 'black_householders',\n",
    "    'B23025_005': 'unemployed_16plus',\n",
    "    'B17001_002': 'poverty_past_12_months',\n",
    "    'B25077_001': 'median_home_value',\n",
    "    'B25035_001': 'median_year_built',\n",
    "})\n",
    "\n",
    "if (debug >= 1):\n",
    "    print('...race, owner-occupied units, mfi, vacants...')\n",
    "data_gdf = data_gdf.rename(columns = {    \n",
    "    'B02001_002': 'pop_white',\n",
    "    'B02001_003': 'pop_black',\n",
    "    'B25003_002': 'owner_occ_housing_units',\n",
    "    'B19113_001': 'mfi',\n",
    "    'B25002_002': 'occupied_housing_units',\n",
    "    'B25002_003': 'vacant_housing_units',\n",
    "    'B11001_002': 'household_type_family',\n",
    "    'B11001_007': 'household_type_nonfamily',\n",
    "    'B25003_003': 'tenure_rent',\n",
    "    'B25034_001': 'year_built'\n",
    "})\n",
    "\n",
    "if (debug >= 1):\n",
    "    print('....comparison variables: total population, total households, poverty status, rentage...')\n",
    "data_gdf = data_gdf.rename(columns = {\n",
    "    'B01001_001': 'pop_total',\n",
    "    'B02001_001': 'pop_by_race_total',\n",
    "    'B11001_001': 'total_households',\n",
    "    'B25002_001': 'total_housing_units',\n",
    "    'B23025_002': 'labor_force_16plus',\n",
    "    'B17001_001': 'poverty_status_known',\n",
    "    'B25003_001': 'tenure_total'\n",
    "})\n",
    "\n",
    "if (debug >= 1):\n",
    "    print('MFI & median home value: substituting \".\" with np.nan, converting to numeric...')\n",
    "    print('median home value: substituting \".\" with np.nan, converting to numeric...')\n",
    "data_gdf.loc[data_gdf['mfi'] == '.', 'mfi'] = pandas.to_numeric(data_gdf['mfi'][data_gdf['mfi'] == '.'], errors='coerce')\n",
    "data_gdf.loc[data_gdf['median_home_value'] == '.', 'median_home_value'] = pandas.to_numeric(data_gdf['median_home_value'][data_gdf['median_home_value'] == '.'], errors='coerce')\n",
    "\n",
    "\n",
    "# print('backing up...')\n",
    "# data_gdf_bk = data_gdf\n",
    "\n",
    "e = time.time()\n",
    "g = g + (e-s)\n",
    "\n",
    "\n",
    "print('Done in {0:,.0f} seconds!'.format(e-s))\n",
    "\n",
    "\n",
    "# unrenamed_columns = []\n",
    "# for x in data_gdf.columns:\n",
    "#     if ((x[0] == 'B') and ('_err' not in x)):\n",
    "#         unrenamed_columns.append(x)\n",
    "# acs5_metadata_df[acs5_metadata_df['variable'].isin(unrenamed_columns)][-9:]\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Create error calculating functions"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Defined standard-error-calculating functions!\n",
      "ok\n"
     ]
    }
   ],
   "source": [
    "### Guide on how to calculate errors in percentages:\n",
    "# https://www.census.gov/content/dam/Census/library/publications/2018/acs/acs_general_handbook_2018_ch08.pdf\n",
    "    \n",
    "## Aggregating Data Across Population Subgroups: add error for each group in quadrature, divide by 1.645 for serr\n",
    "s = time.time()\n",
    "\n",
    "def find_serr_educated(row):\n",
    "\n",
    "    return pandas.to_numeric(np.sqrt(row['B15003_022_err']**2 + row['B15003_023_err']**2 + row['B15003_024_err']**2 + row['B15003_025_err']**2\n",
    "                                ) / 1.645, errors='coerce')\n",
    "\n",
    "# def find_serr_householders(row):\n",
    "#     return pandas.to_numeric(np.sqrt(row['B11001_002_err']**2 + row['B11001_007_err']**2 \n",
    "#                                 ) / 1.645, errors='coerce')\n",
    "\n",
    "e = time.time()\n",
    "g = g + (e-s)\n",
    "print('Defined standard-error-calculating functions!')\n",
    "print('ok')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Calculate errors"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "...standard errors for hs graduates 25 and older (using custom serr-finding function...\n",
      "...margins of error for householder sex & race, unempoyment, poverty, home value, home age...\n",
      "MFI & median home value: substituting \".\" with np.nan, converting to numeric...\n",
      "\n",
      "calculating and renaming margins of error for comparison variables...\n",
      "...race, owner-occupied units, mfi...\n",
      "Done in 0 seconds!\n"
     ]
    }
   ],
   "source": [
    "s = time.time()\n",
    "\n",
    "# print('getting from backup..')\n",
    "# data_gdf = data_gdf_bk\n",
    "\n",
    "if (debug >= 1):\n",
    "    print('...standard errors for hs graduates 25 and older (using custom serr-finding function...')\n",
    "data_gdf = data_gdf.assign(educated_serr = pandas.to_numeric(data_gdf.apply(lambda row: find_serr_educated(row), axis=1), errors='coerce'))\n",
    "\n",
    "if (debug >= 1):\n",
    "    print('...margins of error for householder sex & race, unempoyment, poverty, home value, home age...')\n",
    "data_gdf = data_gdf.rename(columns = {     \n",
    "    'B11001_001_err': 'total_households_err',\n",
    "    'B11001_006_err': 'female_householders_err',\n",
    "    'B11001A_001_err': 'black_householders_err',\n",
    "    'B11001B_001_err': 'white_householders_err',\n",
    "    'B23025_005_err': 'unemployed_16plus_err',\n",
    "    'B17001_002_err': 'poverty_past_12_months_err',\n",
    "    'B25077_001_err': 'median_home_value_err',\n",
    "    'B25035_001_err': 'median_year_built_err',\n",
    "    'B25003_002_err': 'owner_occ_housing_units_err',\n",
    "    'B19113_001_err': 'mfi_err',\n",
    "    'B25002_002_err': 'occupied_housing_units_err',\n",
    "    'B25002_003_err': 'vacant_housing_units_err',\n",
    "    'B11001_002_err': 'household_type_family_err',\n",
    "    'B11001_007_err': 'household_type_nonfamily_err',\n",
    "    'B25003_003_err': 'tenure_rent_err',\n",
    "    'B15003_001_err': 'pop_25plus_err'\n",
    "})\n",
    "\n",
    "data_gdf = data_gdf.assign(household_type_total = data_gdf['household_type_family'] + data_gdf['household_type_nonfamily'])\n",
    "data_gdf = data_gdf.assign(household_type_total_err = data_gdf['household_type_family_err'] + data_gdf['household_type_nonfamily_err'])\n",
    "\n",
    "if (debug >= 1):\n",
    "    print('MFI & median home value: substituting \".\" with np.nan, converting to numeric...')\n",
    "data_gdf.loc[data_gdf['mfi_err'] == '.', 'mfi_err'] = pandas.to_numeric(data_gdf['mfi_err'][data_gdf['mfi_err'] == '.'], errors='coerce')\n",
    "data_gdf.loc[data_gdf['median_home_value_err'] == '.', 'median_home_value_err'] = pandas.to_numeric(data_gdf['median_home_value_err'][data_gdf['median_home_value_err'] == '.'], errors='coerce')\n",
    "\n",
    "\n",
    "print('\\ncalculating and renaming margins of error for comparison variables...')\n",
    "if (debug >= 1):\n",
    "    print('...race, owner-occupied units, mfi...')\n",
    "data_gdf = data_gdf.rename(columns = {\n",
    "    'B01001_001_err': 'pop_total_err',\n",
    "    'B02001_002_err': 'pop_white_err',\n",
    "    'B02001_003_err': 'pop_black_err',\n",
    "    'B02001_001_err': 'pop_by_race_total',\n",
    "    'B17001_001_err': 'poverty_status_known_err',\n",
    "    'B02001_001_err': 'pop_by_race_total_err',\n",
    "    'B25002_001_err': 'total_housing_units_err',\n",
    "    'B23025_002_err': 'labor_force_16plus_err',\n",
    "    'B25003_001_err': 'tenure_total_err',\n",
    "    \n",
    "    \n",
    "})\n",
    "\n",
    "# if (debug >= 1):\n",
    "#     print('...total householders...')\n",
    "# data_gdf = data_gdf.assign(total_householders_serr = pandas.to_numeric(data_gdf.apply(lambda row: find_serr_householders(row), axis=1), errors='coerce'))\n",
    "\n",
    "\n",
    "#print('dropping columns we do not care about...')\n",
    "# columns_do_not_care = []\n",
    "# columns_do_not_care += ['B11001_002', 'B11001_007']\n",
    "# columns_do_not_care += ['B11001_002_err', 'B11001_007_err']\n",
    "# columns_do_not_care += ['B01001_048_err','B01001_049_err']#,'STATE']\n",
    "#data_gdf = data_gdf.drop(columns_do_not_care, axis=1)\n",
    "\n",
    "\n",
    "# print('Calculated errors for all columns!')\n",
    "# print('backing up...')\n",
    "# data_gdf_bk = data_gdf\n",
    "\n",
    "e = time.time()\n",
    "g = g + (e-s)\n",
    "\n",
    "unrenamed_columns_with_errs = []\n",
    "\n",
    "for x in data_gdf.columns:\n",
    "#    print(x)\n",
    "    if ((x[0] == 'B') and ('_err' in x)):\n",
    "        unrenamed_columns_with_errs.append(x)\n",
    "\n",
    "unrenamed_columns = []        \n",
    "for y in unrenamed_columns_with_errs:\n",
    "    unrenamed_columns.append(y[:-4])\n",
    "#acs5_metadata_df[acs5_metadata_df['variable'].isin(unrenamed_columns)][80:]\n",
    "#print(unrenamed_columns)\n",
    "\n",
    "e = time.time()\n",
    "print('Done in {0:,.0f} seconds!'.format(e-s))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Prepare to percentify"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "ok\n"
     ]
    }
   ],
   "source": [
    "s = time.time()\n",
    "vars_for_percentification = ['pop_white', 'pop_black', 'black_householders', 'white_householders']\n",
    "vars_for_percentification += ['owner_occ_housing_units', 'educated', 'female_householders']\n",
    "vars_for_percentification += ['unemployed_16plus', 'poverty_past_12_months']\n",
    "vars_for_percentification += ['tenure_rent', 'household_type_family']\n",
    "\n",
    "vars_for_percentification += ['pop_white_err', 'pop_black_err', 'black_householders_err', 'white_householders_err']\n",
    "vars_for_percentification += ['owner_occ_housing_units_err', 'educated_serr', 'female_householders_err']\n",
    "vars_for_percentification += ['unemployed_16plus_err', 'poverty_past_12_months_err']\n",
    "vars_for_percentification += ['tenure_rent_err', 'household_type_family_err']\n",
    "\n",
    "vars_for_percentification += ['pop_total', 'total_householders', 'pop_by_race_total', 'pop_25plus', 'labor_force_16plus']\n",
    "vars_for_percentification += ['poverty_status_known', 'vacant_housing_units', 'total_housing_units']\n",
    "vars_for_percentification += ['tenure_total', 'household_type_total']\n",
    "\n",
    "vars_for_percentification += ['pop_totatl_err', 'total_householders_serr', 'pop_by_race_total_err', 'pop_25plus_err', 'labor_force_16plus_err']\n",
    "vars_for_percentification += ['poverty_status_known_err', 'vacant_housing_units_err', 'total_housing_units_err']\n",
    "vars_for_percentification += ['tenure_total_err', 'household_type_total_err']\n",
    "#vars_for_percentification\n",
    "#city_tracts_years_df[vars_for_percentification].columns.tolist()\n",
    "e = time.time()\n",
    "g = g + (e-s)\n",
    "\n",
    "print('ok')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Calculate percentages for needed demographic variables"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "ok\n"
     ]
    },
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead tr th {\n",
       "        text-align: left;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr>\n",
       "      <th>GEOID</th>\n",
       "      <th>14000US42003262000</th>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>year</th>\n",
       "      <th>2022</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>STATEFP</th>\n",
       "      <td>42</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>county</th>\n",
       "      <td>3</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>census_tract</th>\n",
       "      <td>2620</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>NAME</th>\n",
       "      <td>2620</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>NAMELSAD</th>\n",
       "      <td>Census Tract 2620</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>...</th>\n",
       "      <td>...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>pct_unemployed</th>\n",
       "      <td>0.0773639</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>pct_poverty</th>\n",
       "      <td>0.296055</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>pct_vacant</th>\n",
       "      <td>0.126631</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>pct_rent</th>\n",
       "      <td>0.434974</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>pct_household_family</th>\n",
       "      <td>0.487698</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "<p>332 rows × 1 columns</p>\n",
       "</div>"
      ],
      "text/plain": [
       "GEOID                14000US42003262000\n",
       "year                               2022\n",
       "STATEFP                              42\n",
       "county                                3\n",
       "census_tract                       2620\n",
       "NAME                               2620\n",
       "NAMELSAD              Census Tract 2620\n",
       "...                                 ...\n",
       "pct_unemployed                0.0773639\n",
       "pct_poverty                    0.296055\n",
       "pct_vacant                     0.126631\n",
       "pct_rent                       0.434974\n",
       "pct_household_family           0.487698\n",
       "\n",
       "[332 rows x 1 columns]"
      ]
     },
     "execution_count": 19,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "s = time.time()\n",
    "\n",
    "data_gdf = data_gdf.assign(pct_white = pandas.to_numeric((data_gdf['pop_white'] / data_gdf['pop_by_race_total']), errors='coerce'))\n",
    "data_gdf = data_gdf.assign(pct_black = pandas.to_numeric((data_gdf['pop_black'] / data_gdf['pop_by_race_total']), errors='coerce'))\n",
    "\n",
    "data_gdf = data_gdf.assign(pct_white_householders = pandas.to_numeric((data_gdf['white_householders'] / data_gdf['total_households']), errors='coerce'))\n",
    "data_gdf = data_gdf.assign(pct_black_householders = pandas.to_numeric((data_gdf['black_householders'] / data_gdf['total_households']), errors='coerce'))\n",
    "data_gdf = data_gdf.assign(pct_female_householders = pandas.to_numeric((data_gdf['female_householders'] / data_gdf['total_households']), errors='coerce'))\n",
    "\n",
    "data_gdf = data_gdf.assign(peducated = pandas.to_numeric(data_gdf['educated'], errors='coerce') / pandas.to_numeric(data_gdf['pop_25plus'], errors='coerce'))\n",
    "data_gdf = data_gdf.assign(pct_unemployed = pandas.to_numeric(data_gdf['unemployed_16plus'], errors='coerce') / pandas.to_numeric(data_gdf['labor_force_16plus'], errors='coerce'))\n",
    "data_gdf = data_gdf.assign(pct_poverty = pandas.to_numeric(data_gdf['poverty_past_12_months'], errors='coerce') / pandas.to_numeric(data_gdf['poverty_status_known'], errors='coerce'))\n",
    "data_gdf = data_gdf.assign(pct_vacant = pandas.to_numeric(data_gdf['vacant_housing_units'], errors='coerce') / pandas.to_numeric(data_gdf['total_housing_units'], errors='coerce'))\n",
    "\n",
    "data_gdf = data_gdf.assign(pct_rent = pandas.to_numeric(data_gdf['tenure_rent'], errors='coerce') / pandas.to_numeric(data_gdf['tenure_total'], errors='coerce'))\n",
    "data_gdf = data_gdf.assign(pct_household_family = pandas.to_numeric(data_gdf['household_type_family'], errors='coerce') / pandas.to_numeric(data_gdf['household_type_total'], errors='coerce'))\n",
    "\n",
    "e = time.time()\n",
    "g = g + (e-s)\n",
    "print('ok')\n",
    "#data_gdf.columns.tolist()\n",
    "\n",
    "data_gdf.sample(1).T\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Functions to calculate errors in percentages"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Defined functions to calculate standard errors in percentages!\n",
      "ok\n"
     ]
    }
   ],
   "source": [
    "#Guide on how to do this:\n",
    "#### https://www.census.gov/content/dam/Census/library/publications/2018/acs/acs_general_handbook_2018_ch08.pdf\n",
    "\n",
    "# X and Y are the measured values (not the errors) - X for the subsgroup and Y for the whole sample\n",
    "# Let P = X/Y  (the proportion we calculated in the last step)\n",
    "# dX and dY are the measured errors\n",
    "# dP = (1/Y) * np.sqrt(dX**2 - (P**2 * dY**2))\n",
    "# Standard error of P is dP/1.645\n",
    "#### this calculation is done verbosely in fnid_pop_white_serr, quickly in other functions\n",
    "\n",
    "s = time.time()\n",
    "def find_errors_in_pct(X, Y, dX, dY, verboselevel = 0):\n",
    "    try:\n",
    "        P = X / Y\n",
    "        oneoverY = 1 / Y\n",
    "        dXsq = dX**2\n",
    "        dYsq = dY**2\n",
    "        Psq = P**2\n",
    "        PsqdYsq = Psq * dYsq\n",
    "        if (PsqdYsq <= dXsq):\n",
    "            underroot = dXsq - PsqdYsq\n",
    "        else:\n",
    "            underroot = dXsq + PsqdYsq\n",
    "        rooty = np.sqrt(underroot)\n",
    "        dP = oneoverY * rooty\n",
    "        SE = dP / 1.645\n",
    "        if (verboselevel >= 2):\n",
    "#            print('X = pop_white, Y = pop_total')\n",
    "            print('X = {0:.0f}, dX = {1:.0f} ({2:.1%} error)'.format(X, dX, dX/X))\n",
    "            print('Y = {0:.0f}, dY = {1:.0f} ({2:.1%} error)'.format(Y, dY, dY/Y))\n",
    "        if (verboselevel >= 3):\n",
    "            print('P = {0:.3f}'.format(P))\n",
    "            print('dXsq = {0:.0f}, dYsq = {1:.0f}, Psq = {2:.3f}'.format(dXsq, dYsq, Psq))\n",
    "            print('PsqdYsq = {0:.0f}, underroot = {1:.0f}, rooty = {2:.3f}'.format(PsqdYsq, underroot, rooty))\n",
    "            print('dP = {0:.3f}'.format(dP))\n",
    "            print('SE = {0:.3f}'.format(SE))\n",
    "        if (verboselevel >= 2):\n",
    "            print('RESULT: {0:.2%} +/- {1:.2%}'.format(P, SE)) \n",
    "            print('\\n')\n",
    "        return pandas.to_numeric(SE, errors='coerce')\n",
    "    except ZeroDivisionError:\n",
    "        return np.nan\n",
    "    \n",
    "e = time.time()\n",
    "g = g + (e-s)    \n",
    "print('Defined functions to calculate standard errors in percentages!')\n",
    "print('ok')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Calculate errors in percntages"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Calculating errors in percentages...\n",
      "ok\n"
     ]
    }
   ],
   "source": [
    "verboselevel = 0\n",
    "s = time.time()\n",
    "\n",
    "print('Calculating errors in percentages...')\n",
    "data_gdf = data_gdf.assign(pct_white_serr = np.nan)\n",
    "data_gdf = data_gdf.assign(pct_black_serr = np.nan)\n",
    "data_gdf = data_gdf.assign(pct_white_householders_serr = np.nan)\n",
    "data_gdf = data_gdf.assign(pct_black_householders_serr = np.nan)\n",
    "data_gdf = data_gdf.assign(pct_female_householders_serr = np.nan)\n",
    "data_gdf = data_gdf.assign(peducated_serr = np.nan)\n",
    "data_gdf = data_gdf.assign(pct_unemployed_serr = np.nan)\n",
    "data_gdf = data_gdf.assign(pct_poverty_serr = np.nan)\n",
    "data_gdf = data_gdf.assign(pct_vacant_serr = np.nan)\n",
    "\n",
    "\n",
    "data_gdf = data_gdf.assign(pct_rent_serr = np.nan)\n",
    "data_gdf = data_gdf.assign(pct_household_family_serr = np.nan)\n",
    "\n",
    "\n",
    "#data_gdf.loc[:, \n",
    "#              'poverty_status_known_last12months_total_err'] = pandas.to_numeric(data_gdf['poverty_status_known_last12months_total_err'], errors='coerce')\n",
    "\n",
    "\n",
    "for ix, thisrow in data_gdf.iterrows():\n",
    "    if (verboselevel >= 2):\n",
    "        print('Census tract {0:}...'.format(ix))\n",
    "    #print('pct_white_serr...')\n",
    "    data_gdf.loc[ix, 'pct_white_serr'] = find_errors_in_pct(thisrow['pop_white'], thisrow['pop_by_race_total'], thisrow['pop_white_err'], thisrow['pop_by_race_total_err'], verboselevel)\n",
    "    #print('pct_black_serr...')\n",
    "    data_gdf.loc[ix, 'pct_black_serr'] = find_errors_in_pct(thisrow['pop_black'], thisrow['pop_by_race_total'], thisrow['pop_black_err'], thisrow['pop_by_race_total_err'], verboselevel)\n",
    "    #print('pct_white_householders_serr...')\n",
    "    data_gdf.loc[ix, 'pct_white_householders_serr'] = find_errors_in_pct(thisrow['white_householders'], thisrow['total_households'], thisrow['white_householders_err'], thisrow['total_households_err'], verboselevel)\n",
    "    data_gdf.loc[ix, 'pct_black_householders_serr'] = find_errors_in_pct(thisrow['white_householders'], thisrow['total_households'], thisrow['white_householders_err'], thisrow['total_households_err'], verboselevel)\n",
    "    data_gdf.loc[ix, 'pct_female_householders_serr'] = find_errors_in_pct(thisrow['white_householders'], thisrow['total_households'], thisrow['white_householders_err'], thisrow['total_households_err'], verboselevel)\n",
    "    \n",
    "    data_gdf.loc[ix, 'peducated_serr'] = find_errors_in_pct(thisrow['educated'], thisrow['pop_25plus'], thisrow['educated_serr'], thisrow['pop_25plus_err'], verboselevel)\n",
    "    data_gdf.loc[ix, 'pct_unemployed_serr'] = find_errors_in_pct(thisrow['unemployed_16plus'], thisrow['unemployed_16plus_err'], thisrow['labor_force_16plus_err'], thisrow['labor_force_16plus_err'], verboselevel)\n",
    "    data_gdf.loc[ix, 'pct_poverty_serr'] = find_errors_in_pct(thisrow['poverty_past_12_months'], thisrow['poverty_status_known'], thisrow['poverty_past_12_months_err'], thisrow['poverty_status_known_err'], verboselevel)\n",
    "    data_gdf.loc[ix, 'pct_vacant_serr'] = find_errors_in_pct(thisrow['vacant_housing_units'], thisrow['total_housing_units'], thisrow['vacant_housing_units_err'], thisrow['total_housing_units_err'], verboselevel)\n",
    "\n",
    "    data_gdf.loc[ix, 'pct_rent_serr'] = find_errors_in_pct(thisrow['tenure_rent'], thisrow['tenure_total'], thisrow['tenure_rent_err'], thisrow['tenure_total_err'], verboselevel)\n",
    "    data_gdf.loc[ix, 'pct_household_family_serr'] = find_errors_in_pct(thisrow['household_type_family'], thisrow['household_type_total'], thisrow['household_type_family_err'], thisrow['household_type_total_err'], verboselevel)\n",
    "\n",
    "if (verboselevel >= 1):\n",
    "    for ix, thisrow in data_gdf.iterrows():\n",
    "        print('Census tract {0:,.0f}'.format(ix))\n",
    "        print('{0:,.0f} +/- {1:,.0f} white'.format(\n",
    "            thisrow['pop_white'], thisrow['pop_white_err']\n",
    "        ))\n",
    "        print('{0:,.0f} +/- {1:,.0f} total'.format(\n",
    "            thisrow['pop_total'], thisrow['pop_total_err']\n",
    "        ))\n",
    "        print('{0:.1%} +/- {1:.1%}'.format(\n",
    "            thisrow['pct_white'], thisrow['pct_white_serr']\n",
    "        ))\n",
    "        print('\\n')\n",
    "\n",
    "# print('backing up...')\n",
    "# data_gdf_bk = data_gdf\n",
    "    \n",
    "e = time.time()\n",
    "g = g + (e-s)\n",
    "print('ok')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Correct for inflation"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Get inflation values from Board of Labor Statistics Consumer Price Index data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>themonth</th>\n",
       "      <th>thevalue</th>\n",
       "      <th>thefactor</th>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>theyear</th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>2022</th>\n",
       "      <td>1</td>\n",
       "      <td>282.39</td>\n",
       "      <td>1.05815</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2021</th>\n",
       "      <td>1</td>\n",
       "      <td>262.518</td>\n",
       "      <td>1.13825</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2020</th>\n",
       "      <td>1</td>\n",
       "      <td>258.906</td>\n",
       "      <td>1.15413</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2019</th>\n",
       "      <td>1</td>\n",
       "      <td>252.561</td>\n",
       "      <td>1.18313</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2018</th>\n",
       "      <td>1</td>\n",
       "      <td>248.859</td>\n",
       "      <td>1.20073</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2017</th>\n",
       "      <td>1</td>\n",
       "      <td>243.618</td>\n",
       "      <td>1.22656</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2016</th>\n",
       "      <td>1</td>\n",
       "      <td>237.652</td>\n",
       "      <td>1.25735</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2015</th>\n",
       "      <td>1</td>\n",
       "      <td>234.747</td>\n",
       "      <td>1.27291</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2014</th>\n",
       "      <td>1</td>\n",
       "      <td>235.288</td>\n",
       "      <td>1.26998</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2013</th>\n",
       "      <td>1</td>\n",
       "      <td>231.679</td>\n",
       "      <td>1.28977</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2012</th>\n",
       "      <td>1</td>\n",
       "      <td>227.842</td>\n",
       "      <td>1.31149</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "        themonth thevalue thefactor\n",
       "theyear                            \n",
       "2022           1   282.39   1.05815\n",
       "2021           1  262.518   1.13825\n",
       "2020           1  258.906   1.15413\n",
       "2019           1  252.561   1.18313\n",
       "2018           1  248.859   1.20073\n",
       "2017           1  243.618   1.22656\n",
       "2016           1  237.652   1.25735\n",
       "2015           1  234.747   1.27291\n",
       "2014           1  235.288   1.26998\n",
       "2013           1  231.679   1.28977\n",
       "2012           1  227.842   1.31149"
      ]
     },
     "execution_count": 22,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "s = time.time()\n",
    "baseline_year = latest_year\n",
    "\n",
    "inflator_df = pandas.DataFrame(data=None, columns=['theyear','themonth','thevalue'])\n",
    "\n",
    "headers = {'Content-type': 'application/json'}\n",
    "data = json.dumps({\"seriesid\": ['CUSR0000SA0'],\"startyear\":earliest_year, \"endyear\":latest_year})  # All items in U.S. city average, all urban consumers, seasonally adjusted: https://www.bls.gov/cpi/tables/supplemental-files/historical-cpi-u-202007.pdf\n",
    "p = requests.post('https://api.bls.gov/publicAPI/v2/timeseries/data/', data=data, headers=headers)\n",
    "json_data = json.loads(p.text)\n",
    "\n",
    "if (json_data['status'] == 'REQUEST_NOT_PROCESSED'):\n",
    "    print('Request not processed, setting manually...')\n",
    "    inflator_df.loc[0, ['theyear', 'themonth', 'thevalue']] = [2021, 1, 262.518]\n",
    "    inflator_df.loc[1, ['theyear', 'themonth', 'thevalue']] = [2020, 1, 258.906]\n",
    "    inflator_df.loc[2, ['theyear', 'themonth', 'thevalue']] = [2019, 1, 251.712]\n",
    "    inflator_df.loc[3, ['theyear', 'themonth', 'thevalue']] = [2018, 1, 247.867]\n",
    "    inflator_df.loc[4, ['theyear', 'themonth', 'thevalue']] = [2017, 1, 242.839]\n",
    "    inflator_df.loc[5, ['theyear', 'themonth', 'thevalue']] = [2016, 1, 237.652]\n",
    "    inflator_df.loc[6, ['theyear', 'themonth', 'thevalue']] = [2015, 1, 234.747]\n",
    "    inflator_df.loc[7, ['theyear', 'themonth', 'thevalue']] = [2014, 1, 235.288]\n",
    "    inflator_df.loc[8, ['theyear', 'themonth', 'thevalue']] = [2013, 1, 231.679]\n",
    "    inflator_df.loc[9, ['theyear', 'themonth', 'thevalue']] = [2012, 1, 227.842]\n",
    "#     inflator_df.loc[10, ['theyear', 'themonth', 'thevalue']] = [2011, 1, 221.187]\n",
    "#     inflator_df.loc[11, ['theyear', 'themonth', 'thevalue']] = [2010, 1, 217.488]\n",
    "else:\n",
    "    cnt = 0\n",
    "    for series in json_data['Results']['series']:\n",
    "        seriesId = series['seriesID']\n",
    "        for item in series['data']:\n",
    "            inflator_df.loc[cnt, ['theyear', 'themonth', 'thevalue']] = [int(item['year']), int(item['period'][1:]), float(item['value'])]\n",
    "            cnt += 1 \n",
    "\n",
    "## 2022 not included because reasons, add it manually here        \n",
    "if (json_data['status'] == 'REQUEST_NOT_PROCESSED'):\n",
    "    inflator_df.loc[inflator_df.index.max()+1, ['theyear', 'themonth', 'thevalue']] = [2022, 1, 282.390]\n",
    "else:\n",
    "    inflator_df.loc[inflator_df.index.max()+1, ['theyear', 'themonth', 'thevalue']] = [2022, 1, 282.390]\n",
    "    inflator_df.loc[inflator_df.index.max()+1, ['theyear', 'themonth', 'thevalue']] = [2022, 2, 284.535]\n",
    "    inflator_df.loc[inflator_df.index.max()+1, ['theyear', 'themonth', 'thevalue']] = [2022, 3, 287.553] \n",
    "    inflator_df.loc[inflator_df.index.max()+1, ['theyear', 'themonth', 'thevalue']] = [2022, 4, 288.764]\n",
    "    inflator_df.loc[inflator_df.index.max()+1, ['theyear', 'themonth', 'thevalue']] = [2022, 5, 291.359]\n",
    "    inflator_df.loc[inflator_df.index.max()+1, ['theyear', 'themonth', 'thevalue']] = [2022, 6, 294.996]\n",
    "    inflator_df.loc[inflator_df.index.max()+1, ['theyear', 'themonth', 'thevalue']] = [2022, 7, 294.977]\n",
    "    inflator_df.loc[inflator_df.index.max()+1, ['theyear', 'themonth', 'thevalue']] = [2022, 8, 295.209]\n",
    "    inflator_df.loc[inflator_df.index.max()+1, ['theyear', 'themonth', 'thevalue']] = [2022, 2, 284.535]\n",
    "    inflator_df.loc[inflator_df.index.max()+1, ['theyear', 'themonth', 'thevalue']] = [2022, 9, 296.341]\n",
    "    inflator_df.loc[inflator_df.index.max()+1, ['theyear', 'themonth', 'thevalue']] = [2022, 10, 297.863]\n",
    "    inflator_df.loc[inflator_df.index.max()+1, ['theyear', 'themonth', 'thevalue']] = [2022, 11, 298.648]\n",
    "    inflator_df.loc[inflator_df.index.max()+1, ['theyear', 'themonth', 'thevalue']] = [2022, 12, 298.812]\n",
    "    \n",
    "inflator_df = inflator_df.sort_values(by=['theyear', 'themonth'], ascending=False).reset_index(drop=True)\n",
    "baseline_value = inflator_df[inflator_df['theyear'] == baseline_year]['thevalue'].values[0]\n",
    "inflator_df = inflator_df.assign(thefactor = baseline_value / inflator_df['thevalue'])\n",
    "inflator_df = inflator_df.set_index('theyear')\n",
    "inflator_df = inflator_df[inflator_df['themonth'] == 1]   # inflate by year by choosing only January\n",
    "\n",
    "e = time.time()\n",
    "g = g + (e-s)\n",
    "# #inflator_df[inflator_df['themonth'] == 1]\n",
    "\n",
    "inflator_df"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Pittsburgh: Kept 1,480 tract-years in 0 minutes 0 seconds!\n"
     ]
    }
   ],
   "source": [
    "s = time.time()\n",
    "\n",
    "# print('getting from backup...')\n",
    "# data_gdf = data_gdf_bk\n",
    "\n",
    "money_columns = []\n",
    "money_columns += ['amtLoans1', 'amtLoans100k', 'amtLoans250k', 'amtLoansToSmallest']\n",
    "money_columns += ['avgSmallLoan', 'amtLoans', 'amtWorkingLoans']\n",
    "money_columns += ['amtLoans1_per_totaljob', 'amtLoans100k_per_totaljob', 'amtLoans250k_per_totaljob']\n",
    "money_columns += ['amtLoansToSmallest_per_totaljob', 'amtLoans_per_totaljob', 'amtWorkingLoans_per_totaljob']\n",
    "money_columns += ['amtLoans1_per_sbjob', 'amtLoans100k_per_sbjob', 'amtLoans250k_per_sbjob']\n",
    "money_columns += ['amtLoansToSmallest_per_sbjob', 'amtLoans_per_sbjob', 'amtWorkingLoans_per_sbjob']\n",
    "money_columns += ['mfi', 'median_home_value', 'mfi_err', 'median_home_value_err']\n",
    "\n",
    "\n",
    "\n",
    "adjusted_money_df = pandas.DataFrame(data=None, columns=['GEOID']+money_columns).set_index('GEOID')\n",
    "\n",
    "for thisyear in years:\n",
    "    thefactor = inflator_df.loc[thisyear]['thefactor']\n",
    "    adjusted_money_df_i = data_gdf.xs(thisyear, level='year')[money_columns].apply(lambda x: x*thefactor)\n",
    "    adjusted_money_df_i = adjusted_money_df_i.assign(year = thisyear).reset_index()\n",
    "    adjusted_money_df = pandas.concat((adjusted_money_df, adjusted_money_df_i), axis=0)\n",
    "\n",
    "for x in money_columns:\n",
    "    adjusted_money_df = adjusted_money_df.rename(columns={x: x+'_adj'})\n",
    "# adjusted_money_df = adjusted_money_df.reset.set_index(['GEOID', 'year'])\n",
    "    \n",
    "# # data_gdf = pandas.concat((data_gdf, adjusted_money_df), axis=1)\n",
    "#.set_index(['GEOID', 'year'])\n",
    "adjusted_money_df = adjusted_money_df.set_index(['GEOID', 'year'])\n",
    "data_gdf = pandas.concat((data_gdf,  adjusted_money_df), axis=1)\n",
    "\n",
    "e = time.time()\n",
    "g = g + (e-s)\n",
    "\n",
    "print('{0:}: Kept {1:,.0f} tract-years in {2:,.0f} minutes {3:.0f} seconds!'.format(city, len(data_gdf), np.floor((e-s)/60), np.floor((e-s)%60)))\n",
    "#data_gdf.reset_index().groupby(['GEOID', 'year']).size().sort_values()\n",
    "#adjusted_money_df.xs(2012, level='year')\n",
    "#adjusted_money_df.reset_index().groupby(['GEOID', 'year']).size().sort_values(ascending=False)\n",
    "\n",
    "#data_gdf.reset_index().groupby(['GEOID', 'year']).size().sort_values(ascending=False)\n",
    "\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Add community statistical areas (Baltimore only)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 24,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "backing up...\n",
      "Done!\n"
     ]
    }
   ],
   "source": [
    "s = time.time()\n",
    "if (city == 'Baltimore'):\n",
    "    tract_to_csa_df = pandas.read_csv(code_lookup_dir+'census_tract_to_neighborhood.csv')\n",
    "    tract_to_csa_df.loc[:, 'GEOID10'] = tract_to_csa_df['GEOID10'].apply(lambda x: '14000US'+str(x))\n",
    "    data_gdf = data_gdf.reset_index().merge(tract_to_csa_df[['GEOID10','TRACTCE10','NAME10','CSA2010']], left_on='GEOID', right_on='GEOID10').set_index(['GEOID', 'year'])\n",
    "\n",
    "#data_gdf = data_gdf.drop(['GEOID10', 'TRACTCE10'], axis=1)\n",
    "\n",
    "print('backing up...')\n",
    "data_gdf_bk = data_gdf\n",
    "\n",
    "e = time.time()\n",
    "g = g + (e-s)\n",
    "print('Done!')\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Get only the columns we need, in the right order"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 25,
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "getting from backup...\n",
      "selecting column order...\n",
      "Got 365 columns for 1,480 tract-years in 0 minutes 0 seconds!\n"
     ]
    }
   ],
   "source": [
    "s = time.time()\n",
    "\n",
    "print('getting from backup...')\n",
    "data_gdf = geopandas.GeoDataFrame()\n",
    "data_gdf = data_gdf_bk\n",
    "\n",
    "data_gdf = data_gdf.assign(city_name = city)\n",
    "data_gdf = data_gdf.drop('state', axis=1)\n",
    "data_gdf = data_gdf.rename(columns={'STATEFP': 'state'})\n",
    "#data_gdf = data_gdf.drop(['tract_data', 'tract_merr', 'GEOID10', 'TRACTCE10'], axis=1)\n",
    "\n",
    "# print('renaming variables describing this tract\\'s state...')\n",
    "# data_gdf = data_gdf.rename(columns = {'STATE': 'state_acs5'})#, 'state': 'state_loans'})#, 'STATEFP': 'state'})\n",
    "# #data_gdf = data_gdf.rename(columns = {'state': 'state_loans', 'STATEFP': 'state'})\n",
    "\n",
    "# data_gdf = data_gdf.assign(city_name = city)\n",
    "\n",
    "print('selecting column order...')\n",
    "new_columns = []\n",
    "new_columns += ['state', 'county', 'census_tract']\n",
    "new_columns += ['city_name', 'msa'] # city name, metropolitan statistical area (MSA) number\n",
    "new_columns += ['loan_type', 'action_taken_type'] # loan type (always 4 for business) and action taken (always 1 for origination)\n",
    "new_columns += ['income_group_code'] # income group code (1 to 14 or 101 to 105)\n",
    "new_columns += ['income_group', 'cra_level'] # income group and CRA level (human-readable)\n",
    "new_columns += ['nLoans1','amtLoans1','nLoans100k','amtLoans100k','nLoans250k','amtLoans250k','nLoansToSmallest','amtLoansToSmallest'] # number and amount of loans (directly from CRA data)\n",
    "new_columns += ['nLoans', 'amtLoans', 'avgSmallLoan', 'nWorkingLoans', 'amtWorkingLoans'] # calculated total and working loans\n",
    "new_columns += ['total_jobs', 'sb_jobs'] # total and small business jobs\n",
    "\n",
    "# all jobs columns under original names\n",
    "new_columns += ['C000','CA01','CA02','CA03','CE01','CE02','CE03']\n",
    "new_columns += ['CNS01','CNS02','CNS03','CNS04','CNS05','CNS06','CNS07','CNS08','CNS09','CNS10','CNS11','CNS12','CNS13','CNS14','CNS15','CNS16','CNS17','CNS18','CNS19','CNS20']\n",
    "new_columns += ['CR01','CR02','CR03','CR04','CR05','CR07','CT01','CT02','CD01','CD02','CD03','CD04','CS01','CS02']\n",
    "new_columns += ['CFA01','CFA02','CFA03','CFA04','CFA05','CFS01','CFS02','CFS03','CFS04','CFS05']\n",
    "\n",
    "# loans per job\n",
    "new_columns += ['nLoans1_per_totaljob','amtLoans1_per_totaljob','nLoans100k_per_totaljob','amtLoans100k_per_totaljob']\n",
    "new_columns += ['nLoans250k_per_totaljob','amtLoans250k_per_totaljob','nLoansToSmallest_per_totaljob','amtLoansToSmallest_per_totaljob']\n",
    "new_columns += ['nLoans_per_totaljob','amtLoans_per_totaljob','nWorkingLoans_per_totaljob','amtWorkingLoans_per_totaljob']\n",
    "new_columns += ['nLoans1_per_sbjob','amtLoans1_per_sbjob','nLoans100k_per_sbjob','amtLoans100k_per_sbjob','nLoans250k_per_sbjob','amtLoans250k_per_sbjob']\n",
    "new_columns += ['nLoansToSmallest_per_sbjob','amtLoansToSmallest_per_sbjob','nLoans_per_sbjob','amtLoans_per_sbjob','nWorkingLoans_per_sbjob','amtWorkingLoans_per_sbjob']\n",
    "\n",
    "# census count estimates\n",
    "new_columns += ['pop_total','pop_by_race_total','pop_white','pop_black']\n",
    "new_columns += ['total_households','white_householders','black_householders','female_householders']\n",
    "new_columns += ['total_housing_units','occupied_housing_units','vacant_housing_units']\n",
    "new_columns += ['educated','unemployed_16plus','poverty_past_12_months','mfi'] # high school graduates (ages 25+), unemployed, in poverty, MFI, travel time to work\n",
    "new_columns += ['owner_occ_housing_units','median_home_value','year_built', 'median_year_built'] # home\n",
    "new_columns += ['pop_25plus','labor_force_16plus','poverty_status_known']\n",
    "new_columns += ['household_type_family', 'household_type_nonfamily','household_type_total', 'tenure_rent', 'tenure_total']\n",
    "\n",
    "\n",
    "# census count errors\n",
    "new_columns += ['pop_total_err','pop_by_race_total_err','pop_white_err','pop_black_err']\n",
    "new_columns += ['total_households_err','white_householders_err','black_householders_err','female_householders_err']\n",
    "new_columns += ['total_housing_units_err','occupied_housing_units_err','vacant_housing_units_err']\n",
    "new_columns += ['educated_serr','unemployed_16plus_err','poverty_past_12_months_err','mfi_err'] # high school graduates (ages 25+), unemployed, in poverty, MFI, travel time to work\n",
    "new_columns += ['owner_occ_housing_units_err','median_home_value_err', 'median_year_built_err'] # home\n",
    "new_columns += ['pop_25plus_err','labor_force_16plus_err','poverty_status_known_err']\n",
    "new_columns += ['household_type_family_err', 'household_type_nonfamily_err', 'household_type_total_err', 'tenure_rent_err', 'tenure_total_err']\n",
    "\n",
    "# census percentages and their errors\n",
    "new_columns += ['pct_white','pct_black','pct_white_householders','pct_black_householders','pct_female_householders']\n",
    "new_columns += ['peducated','pct_unemployed','pct_poverty','pct_vacant', 'pct_household_family', 'pct_rent']\n",
    "new_columns += ['pct_white_serr','pct_black_serr','pct_white_householders_serr','pct_black_householders_serr','pct_female_householders_serr']\n",
    "new_columns += ['peducated_serr','pct_unemployed_serr','pct_poverty_serr','pct_vacant_serr', 'pct_household_family_serr', 'pct_rent_serr']\n",
    "\n",
    "# # inflation-adjusted money values\n",
    "new_columns += ['amtLoans1_adj','amtLoans100k_adj','amtLoans250k_adj','amtLoansToSmallest_adj']\n",
    "new_columns += ['avgSmallLoan_adj','amtLoans_adj','amtWorkingLoans_adj']\n",
    "new_columns += ['amtLoans1_per_totaljob_adj','amtLoans100k_per_totaljob_adj','amtLoans250k_per_totaljob_adj','amtLoansToSmallest_per_totaljob_adj']\n",
    "new_columns += ['amtLoans_per_totaljob_adj','amtWorkingLoans_per_totaljob_adj']\n",
    "new_columns += ['amtLoans1_per_sbjob_adj','amtLoans100k_per_sbjob_adj','amtLoans250k_per_sbjob_adj','amtLoansToSmallest_per_sbjob_adj','amtLoans_per_sbjob_adj','amtWorkingLoans_per_sbjob_adj']\n",
    "new_columns += ['mfi_adj','median_home_value_adj','mfi_err_adj','median_home_value_err_adj']\n",
    "\n",
    "## geographic information about the census tract\n",
    "new_columns += ['NAME', 'NAMELSAD', 'MTFCC', 'FUNCSTAT', 'ALAND', 'AWATER', 'INTPTLAT', 'INTPTLON']\n",
    "\n",
    "## more info from loans data\n",
    "new_columns += ['split_county_indicator', 'population_classification', 'split_county_indicator', 'population_classification']#, 'geography_name_loans_jobs']\n",
    "\n",
    "## ACS5 raw census data under original variable names\n",
    "new_columns += [x for x in data_gdf.columns if (x not in new_columns) and (x[0] == 'B')]\n",
    "## ACS5 tract metadata\n",
    "#data_gdf = data_gdf.rename(columns =  {'Geography Name_acs5': 'geography_name_acs5'})\n",
    "#new_columns += ['FILEID', 'FILETYPE', 'STUSAB', 'CHARITER', 'SEQUENCE', 'LOGRECNO']#, 'state_acs5',  'Name']#, 'geography_name_acs5']\n",
    "\n",
    "# # if (city == 'Baltimore'):\n",
    "# #     new_columns += ['GEOID10', 'TRACTCE10', 'NAME10', 'CSA2010']\n",
    "\n",
    "# as traditional, geometry goes last\n",
    "#new_columns += ['geometry']\n",
    "\n",
    "# # # check that we got them all\n",
    "data_gdf = data_gdf[new_columns]\n",
    "# # #pprint([x for x in data_gdf.columns if x not in new_columns])\n",
    "\n",
    "\n",
    "e = time.time()\n",
    "g = g + (e-s)\n",
    "print('Got {0:,.0f} columns for {1:,.0f} tract-years in {2:,.0f} minutes {3:,.0f} seconds!'.format(data_gdf.shape[1], data_gdf.shape[0], np.floor((e-s)/60), np.floor((e-s)%60)))\n",
    "# Got 365 columns for 1,600 tract-years in 0 minutes 0 seconds!\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 26,
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "outputting data...\n",
      "\tWrote 1,480 tract-years by 369 columns (without geometry data) to:\n",
      "\t\t/home/idies/workspace/21cc/raddick/community_reinvestment_act/city_data/pittsburgh_2012_2022_year_by_year.csv\n",
      "\n",
      "\n",
      "Saved 365 variables for 1,480 tract-years in 6 minutes 24 seconds!\n"
     ]
    }
   ],
   "source": [
    "print('outputting data...')\n",
    "\n",
    "# print('outputting without shapefile data...')\n",
    "# s = time.time()\n",
    "city_name_for_output = city.replace(\".\",\"\").replace(\" \",\"_\").lower()\n",
    "outfilename = '{0:}_{1:.0f}_{2:.0f}_year_by_year.csv'.format(city_name_for_output, earliest_year, latest_year)\n",
    "\n",
    "\n",
    "\n",
    "data_gdf[new_columns].to_csv(output_data_dir+outfilename)\n",
    "\n",
    "print('\\tWrote {0:,.0f} tract-years by {1:,.0f} columns (without geometry data) to:\\n\\t\\t{2:}'.format(\n",
    "    data_gdf[new_columns].shape[0], data_gdf[new_columns].shape[1], output_data_dir+outfilename\n",
    "     ))\n",
    "print('\\n')\n",
    "\n",
    "# print('outputting shapefile data separately...')\n",
    "# outfilename_shapefiles_geo = '{0:}_{1:.0f}_{2:.0f}_shapefiles.shp'.format(city_name_for_output, earliest_year, latest_year)\n",
    "# shapefile_writer_gdf = geopandas.GeoDataFrame(data=data_gdf['geometry'], columns=['geometry'], crs=data_gdf.crs, geometry='geometry')\n",
    "# shapefile_writer_gdf.to_file(output_data_dir+outfilename_shapefiles_geo)\n",
    "\n",
    "# print('\\tWrote {0:,.0f} tract-years by {1:,.0f} columns (SHAPEFILE DATA ONLY) to \\n\\t{2:}'.format(shapefile_writer_gdf.shape[0], shapefile_writer_gdf.shape[1], output_data_dir+outfilename_shapefiles_geo))\n",
    "\n",
    "e = time.time()\n",
    "g = g + (e-s)\n",
    "print('Saved {0:,.0f} variables for {1:,.0f} tract-years in {2:,.0f} minutes {3:,.0f} seconds!'.format(data_gdf.shape[1], data_gdf.shape[0], np.floor((g)/60), np.floor((g)%60)))\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 27,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>GEOID</th>\n",
       "      <th>year</th>\n",
       "      <th>state</th>\n",
       "      <th>county</th>\n",
       "      <th>census_tract</th>\n",
       "      <th>city_name</th>\n",
       "      <th>msa</th>\n",
       "      <th>loan_type</th>\n",
       "      <th>action_taken_type</th>\n",
       "      <th>income_group_code</th>\n",
       "      <th>...</th>\n",
       "      <th>B19013_001_err</th>\n",
       "      <th>B19013A_001_err</th>\n",
       "      <th>B19013B_001_err</th>\n",
       "      <th>B25003A_001_err</th>\n",
       "      <th>B25003A_002_err</th>\n",
       "      <th>B25003A_003_err</th>\n",
       "      <th>B25003B_001_err</th>\n",
       "      <th>B25003B_002_err</th>\n",
       "      <th>B25003B_003_err</th>\n",
       "      <th>B25034_001_err</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>14000US42003010300</td>\n",
       "      <td>2012</td>\n",
       "      <td>42</td>\n",
       "      <td>3</td>\n",
       "      <td>103.0</td>\n",
       "      <td>Pittsburgh</td>\n",
       "      <td>38300.0</td>\n",
       "      <td>4.0</td>\n",
       "      <td>1.0</td>\n",
       "      <td>3.0</td>\n",
       "      <td>...</td>\n",
       "      <td>27810.0</td>\n",
       "      <td>31273.0</td>\n",
       "      <td>7064.0</td>\n",
       "      <td>68</td>\n",
       "      <td>59</td>\n",
       "      <td>38</td>\n",
       "      <td>30</td>\n",
       "      <td>19</td>\n",
       "      <td>21</td>\n",
       "      <td>82</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>14000US42003010300</td>\n",
       "      <td>2013</td>\n",
       "      <td>42</td>\n",
       "      <td>3</td>\n",
       "      <td>103.0</td>\n",
       "      <td>Pittsburgh</td>\n",
       "      <td>38300.0</td>\n",
       "      <td>4.0</td>\n",
       "      <td>1.0</td>\n",
       "      <td>3.0</td>\n",
       "      <td>...</td>\n",
       "      <td>21260.0</td>\n",
       "      <td>20932.0</td>\n",
       "      <td>15990.0</td>\n",
       "      <td>73</td>\n",
       "      <td>69</td>\n",
       "      <td>40</td>\n",
       "      <td>24</td>\n",
       "      <td>16</td>\n",
       "      <td>19</td>\n",
       "      <td>60</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>14000US42003010300</td>\n",
       "      <td>2014</td>\n",
       "      <td>42</td>\n",
       "      <td>3</td>\n",
       "      <td>103.0</td>\n",
       "      <td>Pittsburgh</td>\n",
       "      <td>38300.0</td>\n",
       "      <td>4.0</td>\n",
       "      <td>1.0</td>\n",
       "      <td>3.0</td>\n",
       "      <td>...</td>\n",
       "      <td>8976.0</td>\n",
       "      <td>9315.0</td>\n",
       "      <td>35325.0</td>\n",
       "      <td>54</td>\n",
       "      <td>51</td>\n",
       "      <td>36</td>\n",
       "      <td>26</td>\n",
       "      <td>11</td>\n",
       "      <td>21</td>\n",
       "      <td>49</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>14000US42003010300</td>\n",
       "      <td>2015</td>\n",
       "      <td>42</td>\n",
       "      <td>3</td>\n",
       "      <td>103.0</td>\n",
       "      <td>Pittsburgh</td>\n",
       "      <td>38300.0</td>\n",
       "      <td>4.0</td>\n",
       "      <td>1.0</td>\n",
       "      <td>3.0</td>\n",
       "      <td>...</td>\n",
       "      <td>13372.0</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>39</td>\n",
       "      <td>28</td>\n",
       "      <td>32</td>\n",
       "      <td>35</td>\n",
       "      <td>14</td>\n",
       "      <td>28</td>\n",
       "      <td>31</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>14000US42003010300</td>\n",
       "      <td>2016</td>\n",
       "      <td>42</td>\n",
       "      <td>3</td>\n",
       "      <td>103.0</td>\n",
       "      <td>Pittsburgh</td>\n",
       "      <td>38300.0</td>\n",
       "      <td>4.0</td>\n",
       "      <td>1.0</td>\n",
       "      <td>3.0</td>\n",
       "      <td>...</td>\n",
       "      <td>16306.0</td>\n",
       "      <td>11616.0</td>\n",
       "      <td>NaN</td>\n",
       "      <td>39</td>\n",
       "      <td>29</td>\n",
       "      <td>39</td>\n",
       "      <td>27</td>\n",
       "      <td>12</td>\n",
       "      <td>30</td>\n",
       "      <td>24</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>...</th>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1475</th>\n",
       "      <td>14000US42003982200</td>\n",
       "      <td>2018</td>\n",
       "      <td>42</td>\n",
       "      <td>3</td>\n",
       "      <td>9822.0</td>\n",
       "      <td>Pittsburgh</td>\n",
       "      <td>38300.0</td>\n",
       "      <td>4.0</td>\n",
       "      <td>1.0</td>\n",
       "      <td>14.0</td>\n",
       "      <td>...</td>\n",
       "      <td>111211.0</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>6</td>\n",
       "      <td>6</td>\n",
       "      <td>10</td>\n",
       "      <td>10</td>\n",
       "      <td>10</td>\n",
       "      <td>10</td>\n",
       "      <td>23</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1476</th>\n",
       "      <td>14000US42003982200</td>\n",
       "      <td>2019</td>\n",
       "      <td>42</td>\n",
       "      <td>3</td>\n",
       "      <td>9822.0</td>\n",
       "      <td>Pittsburgh</td>\n",
       "      <td>38300.0</td>\n",
       "      <td>4.0</td>\n",
       "      <td>1.0</td>\n",
       "      <td>14.0</td>\n",
       "      <td>...</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>10</td>\n",
       "      <td>6</td>\n",
       "      <td>6</td>\n",
       "      <td>10</td>\n",
       "      <td>10</td>\n",
       "      <td>10</td>\n",
       "      <td>10</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1477</th>\n",
       "      <td>14000US42003982200</td>\n",
       "      <td>2020</td>\n",
       "      <td>42</td>\n",
       "      <td>3</td>\n",
       "      <td>9822.0</td>\n",
       "      <td>Pittsburgh</td>\n",
       "      <td>38300.0</td>\n",
       "      <td>4.0</td>\n",
       "      <td>1.0</td>\n",
       "      <td>14.0</td>\n",
       "      <td>...</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>8</td>\n",
       "      <td>6</td>\n",
       "      <td>8</td>\n",
       "      <td>11</td>\n",
       "      <td>11</td>\n",
       "      <td>11</td>\n",
       "      <td>9</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1478</th>\n",
       "      <td>14000US42003982200</td>\n",
       "      <td>2021</td>\n",
       "      <td>42</td>\n",
       "      <td>3</td>\n",
       "      <td>9822.0</td>\n",
       "      <td>Pittsburgh</td>\n",
       "      <td>38300.0</td>\n",
       "      <td>4.0</td>\n",
       "      <td>1.0</td>\n",
       "      <td>14.0</td>\n",
       "      <td>...</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>9</td>\n",
       "      <td>5</td>\n",
       "      <td>7</td>\n",
       "      <td>11</td>\n",
       "      <td>11</td>\n",
       "      <td>11</td>\n",
       "      <td>9</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1479</th>\n",
       "      <td>14000US42003982200</td>\n",
       "      <td>2022</td>\n",
       "      <td>42</td>\n",
       "      <td>3</td>\n",
       "      <td>9822.0</td>\n",
       "      <td>Pittsburgh</td>\n",
       "      <td>38300.0</td>\n",
       "      <td>4.0</td>\n",
       "      <td>1.0</td>\n",
       "      <td>14.0</td>\n",
       "      <td>...</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>11</td>\n",
       "      <td>7</td>\n",
       "      <td>7</td>\n",
       "      <td>16</td>\n",
       "      <td>16</td>\n",
       "      <td>16</td>\n",
       "      <td>11</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "<p>1480 rows × 371 columns</p>\n",
       "</div>"
      ],
      "text/plain": [
       "                   GEOID  year  state  county  census_tract   city_name  \\\n",
       "0     14000US42003010300  2012     42       3         103.0  Pittsburgh   \n",
       "1     14000US42003010300  2013     42       3         103.0  Pittsburgh   \n",
       "2     14000US42003010300  2014     42       3         103.0  Pittsburgh   \n",
       "3     14000US42003010300  2015     42       3         103.0  Pittsburgh   \n",
       "4     14000US42003010300  2016     42       3         103.0  Pittsburgh   \n",
       "...                  ...   ...    ...     ...           ...         ...   \n",
       "1475  14000US42003982200  2018     42       3        9822.0  Pittsburgh   \n",
       "1476  14000US42003982200  2019     42       3        9822.0  Pittsburgh   \n",
       "1477  14000US42003982200  2020     42       3        9822.0  Pittsburgh   \n",
       "1478  14000US42003982200  2021     42       3        9822.0  Pittsburgh   \n",
       "1479  14000US42003982200  2022     42       3        9822.0  Pittsburgh   \n",
       "\n",
       "          msa  loan_type  action_taken_type  income_group_code  ...  \\\n",
       "0     38300.0        4.0                1.0                3.0  ...   \n",
       "1     38300.0        4.0                1.0                3.0  ...   \n",
       "2     38300.0        4.0                1.0                3.0  ...   \n",
       "3     38300.0        4.0                1.0                3.0  ...   \n",
       "4     38300.0        4.0                1.0                3.0  ...   \n",
       "...       ...        ...                ...                ...  ...   \n",
       "1475  38300.0        4.0                1.0               14.0  ...   \n",
       "1476  38300.0        4.0                1.0               14.0  ...   \n",
       "1477  38300.0        4.0                1.0               14.0  ...   \n",
       "1478  38300.0        4.0                1.0               14.0  ...   \n",
       "1479  38300.0        4.0                1.0               14.0  ...   \n",
       "\n",
       "     B19013_001_err B19013A_001_err  B19013B_001_err  B25003A_001_err  \\\n",
       "0           27810.0         31273.0           7064.0               68   \n",
       "1           21260.0         20932.0          15990.0               73   \n",
       "2            8976.0          9315.0          35325.0               54   \n",
       "3           13372.0             NaN              NaN               39   \n",
       "4           16306.0         11616.0              NaN               39   \n",
       "...             ...             ...              ...              ...   \n",
       "1475       111211.0             NaN              NaN                6   \n",
       "1476            NaN             NaN              NaN               10   \n",
       "1477            NaN             NaN              NaN                8   \n",
       "1478            NaN             NaN              NaN                9   \n",
       "1479            NaN             NaN              NaN               11   \n",
       "\n",
       "      B25003A_002_err  B25003A_003_err  B25003B_001_err  B25003B_002_err  \\\n",
       "0                  59               38               30               19   \n",
       "1                  69               40               24               16   \n",
       "2                  51               36               26               11   \n",
       "3                  28               32               35               14   \n",
       "4                  29               39               27               12   \n",
       "...               ...              ...              ...              ...   \n",
       "1475                6               10               10               10   \n",
       "1476                6                6               10               10   \n",
       "1477                6                8               11               11   \n",
       "1478                5                7               11               11   \n",
       "1479                7                7               16               16   \n",
       "\n",
       "      B25003B_003_err  B25034_001_err  \n",
       "0                  21              82  \n",
       "1                  19              60  \n",
       "2                  21              49  \n",
       "3                  28              31  \n",
       "4                  30              24  \n",
       "...               ...             ...  \n",
       "1475               10              23  \n",
       "1476               10              10  \n",
       "1477               11               9  \n",
       "1478               11               9  \n",
       "1479               16              11  \n",
       "\n",
       "[1480 rows x 371 columns]"
      ]
     },
     "execution_count": 27,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "df = pandas.read_csv(output_data_dir+outfilename)\n",
    "df\n",
    "\n",
    "# Baltimore: Saved 365 variables for 2,194 tract-years in 5 minutes 1 seconds!\n",
    "# Washington DC: Saved 365 variables for 2,050 tract-years in 4 minutes 51 seconds!\n",
    "# Detroit: Saved 365 variables for 3,195 tract-years in 6 minutes 28 seconds!\n",
    "# Cleveland: Saved 365 variables for 1,893 tract-years in 6 minutes 33 seconds!\n",
    "# Newark: Saved 365 variables for 960 tract-years in 5 minutes 10 seconds!\n",
    "# St. Louis: Saved 365 variables for 1,160 tract-years in 5 minutes 1 seconds!\n",
    "# Richmond: Saved 365 variables for 753 tract-years in 5 minutes 10 seconds!\n",
    "# San Francisco: Saved 365 variables for 2,308 tract-years in 9 minutes 15 seconds!\n",
    "# Philadelphia: Saved 365 variables for 4,296 tract-years in 7 minutes 31 seconds!\n",
    "# Pittsburgh: Saved 365 variables for 1,480 tract-years in 6 minutes 24 seconds!"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3.8 (py38)",
   "language": "python",
   "name": "py38"
  },
  "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.8.5"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}