gen_business_density.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Business Density"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Generating measures of the number of businesses linked to a hex grid. Data for businesses comes from Statistics Canada's business register division\n",
    "\n",
    "https://dataverse.scholarsportal.info/dataset.xhtml?persistentId=doi:10.5683/SP/FLLHOV&version=2.0"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "First we want to clean the business data and query it by NAICS code. This essentially removes factories / heavy industry, farms, wholesale trade, etc. - businesses that don't serve the public directly"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [
    {
     "ename": "FileNotFoundError",
     "evalue": "[Errno 2] No such file or directory: 'tabular_data/business_data_dec2016/CBP/Ontario-Employees.csv'",
     "output_type": "error",
     "traceback": [
      "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m",
      "\u001b[0;31mFileNotFoundError\u001b[0m                         Traceback (most recent call last)",
      "\u001b[0;32m<ipython-input-1-5a7dda17c603>\u001b[0m in \u001b[0;36m<module>\u001b[0;34m\u001b[0m\n\u001b[1;32m      9\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m     10\u001b[0m \u001b[0;31m# open up the csv\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m---> 11\u001b[0;31m \u001b[0;32mwith\u001b[0m \u001b[0mopen\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m\"tabular_data/business_data_dec2016/CBP/Ontario-Employees.csv\"\u001b[0m\u001b[0;34m)\u001b[0m \u001b[0;32mas\u001b[0m \u001b[0mcsvfile\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m     12\u001b[0m     \u001b[0mreader\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mcsv\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mreader\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mcsvfile\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m     13\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;31mFileNotFoundError\u001b[0m: [Errno 2] No such file or directory: 'tabular_data/business_data_dec2016/CBP/Ontario-Employees.csv'"
     ]
    }
   ],
   "source": [
    "import pandas as pd\n",
    "import csv\n",
    "\n",
    "# which 2-digit naics code to remove from the data\n",
    "no_naics_2digit = [11,21,22,23,31,32,33,41,48,49]\n",
    "\n",
    "# output array\n",
    "out_data = []\n",
    "\n",
    "# open up the csv\n",
    "with open(\"input_data/tabular_data/business_data_dec2016/CBP/Ontario-Employees.csv\") as csvfile:\n",
    "    reader = csv.reader(csvfile)\n",
    "   \n",
    "    # loop over the rows\n",
    "    r = 0\n",
    "    for row in reader:\n",
    "       \n",
    "        # for the first row, with the NAICS code\n",
    "        if r == 0:\n",
    "            \n",
    "            out_header = []\n",
    "            remove_by_index = []\n",
    "            c = 0\n",
    "            \n",
    "            # selecting the columns to keep with if statements\n",
    "            # putting \n",
    "            for cell in row:\n",
    "                if cell == \"Unclassified\":\n",
    "                    outcell = \"000000\"\n",
    "                    out_header.append(outcell)\n",
    "                elif cell == \"\":\n",
    "                    outcell = \"dauid\"\n",
    "                    out_header.append(outcell)\n",
    "                elif cell == \"Total\" or cell == \"Sub-total, classified\":\n",
    "                    remove_by_index.append(c)\n",
    "                else:\n",
    "                    naics2 = (int(cell[0:2]))\n",
    "                    if naics2 in no_naics_2digit:\n",
    "                        remove_by_index.append(c)\n",
    "                    else:\n",
    "                        outcell = cell\n",
    "                        out_header.append(cell[0:6])\n",
    "                c += 1\n",
    "            \n",
    "\n",
    "        # the second row pertains to if there is employment / no employment - ignore for now\n",
    "        elif r == 1:\n",
    "            None\n",
    "        \n",
    "        # removing columns by the index in each row\n",
    "        else:\n",
    "            for index in sorted(remove_by_index, reverse=True):\n",
    "                del row[index]\n",
    "            \n",
    "            int_row = []\n",
    "            for cell in row:\n",
    "                int_row.append(int(cell))\n",
    "                \n",
    "            out_data.append(int_row)\n",
    "            \n",
    "        r += 1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# put this array into pandas data frame\n",
    "\n",
    "df = pd.DataFrame(out_data, columns=out_header)\n",
    "\n",
    "df = df.set_index(\"dauid\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# sum all the rows\n",
    "\n",
    "dfs = pd.DataFrame(df.sum(axis=1))\n",
    "\n",
    "# write to csv\n",
    "\n",
    "dfs.to_csv(\"input_data/tabular_data/DA_2016_business_store_subset.csv\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We then use areal interpolation to link the business data from DAs ($B_d$) to the hex grid ($B_h$)\n",
    "\n",
    "### $B_h = \\sum B_{d} {w_{d,h}}$\n",
    "\n",
    "### $B_h = \\sum B_{d} \\frac{A_{d \\cap h}}{A_b}$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "With postgreSQL and PostGIS ..\n",
    "\n",
    "```sql\n",
    "-- creating a spatial index for the hex grid to speed up the intersection\n",
    "\n",
    "DROP INDEX IF EXISTS hex_grid_200m_gix;\n",
    "CREATE INDEX hex_grid_200m_gix ON hex_grid_200m USING GIST (geom);\n",
    "\n",
    "-- creating an area column $A_b$ for the dissemination area boundaries\n",
    "\n",
    "ALTER TABLE zones_DA11 ADD COLUMN area double precision;\n",
    "UPDATE zones_DA11 SET area = ST_AREA(geom);\n",
    "\n",
    "-- Join the business count data to DAs (2011 geometry since this is how the data was coded)\n",
    "\n",
    "DROP TABLE IF EXISTS temp_DA_business;\n",
    "CREATE TABLE temp_DA_business AS ( \n",
    "    SELECT\n",
    "    zones_DA11.geom AS geom,\n",
    "    zones_DA11.dauid AS dauid,\n",
    "    zones_DA11.area AS area,\n",
    "    coalesce(table_da_business_2016_in_11da.business2016,0) AS business2016 -- nulls to 0\n",
    "    FROM zones_DA11\n",
    "    LEFT OUTER JOIN \n",
    "    table_da_business_2016_in_11da ON table_da_business_2016_in_11da.dauid = zones_DA11.dauid\n",
    ");   \n",
    "\n",
    "-- create a spatial index for this dauid polygon file\n",
    "\n",
    "DROP INDEX IF EXISTS temp_DA_business_gix;\n",
    "CREATE INDEX temp_DA_business_gix ON temp_DA_business USING GIST (geom); \n",
    "\n",
    "-- intersecting with the hex grid, then grouping to get the a set of weights of block to hex grid\n",
    "\n",
    "CREATE TABLE temp_int_DA11_hex AS\n",
    "(\n",
    "SELECT\n",
    "hex_grid_200m.id AS hexid,\n",
    "temp_DA_business.dauid AS dauid,\n",
    "temp_DA_business.area AS area_full,\n",
    "ST_Intersection(hex_grid_200m.geom,temp_DA_business.geom) AS geom\n",
    "FROM\n",
    "hex_grid_200m INNER JOIN temp_DA_business ON ST_Intersects(hex_grid_200m.geom,temp_DA_business.geom)\n",
    ");\n",
    "\n",
    "-- update area of intersected geoms\n",
    "\n",
    "ALTER TABLE temp_int_DA11_hex ADD COLUMN area_int double precision;\n",
    "UPDATE temp_int_DA11_hex SET area_int = ST_AREA(geom);\n",
    "\n",
    "ALTER TABLE temp_int_DA11_hex ADD COLUMN area_ratio double precision;\n",
    "UPDATE temp_int_DA11_hex SET area_ratio = area_int / area_full;\n",
    "\n",
    "-- grouping by hex and DA unique IDs - i.e. this is a weights table that can be used for apportioning data\n",
    "\n",
    "CREATE TABLE weights_da11_hex AS\n",
    "(\n",
    "SELECT\n",
    "hexid,\n",
    "dauid,\n",
    "sum(area_ratio) AS weight\n",
    "FROM\n",
    "temp_int_DA11_hex\n",
    "GROUP BY dauid, hexid\n",
    "ORDER BY dauid, hexid\n",
    ");\n",
    "\n",
    "\n",
    "```"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We now use the weights table to apportion the business data to the hex grid\n",
    "\n",
    "```sql\n",
    "-- joining the businesses to the weights table\n",
    "-- multiplying the businesses by the weight and aggregating by hex grid\n",
    "\n",
    "DROP TABLE IF EXISTS out_data_hex_business2016;\n",
    "CREATE TABLE out_data_hex_business2016 AS (\n",
    "    WITH temp_weight_join AS (\n",
    "        SELECT\n",
    "        weights_da11_hex.dauid AS dauid,\n",
    "        weights_da11_hex.hexid AS hexid,\n",
    "        weights_da11_hex.weight AS weight,\n",
    "        temp_DA_business.business2016 AS business2016\n",
    "        FROM weights_da11_hex\n",
    "        INNER JOIN temp_DA_business ON temp_DA_business.dauid = weights_da11_hex.dauid\n",
    "    ) \n",
    "    SELECT\n",
    "    temp_weight_join.hexid AS hexid,\n",
    "    SUM(temp_weight_join.weight * temp_weight_join.business2016) AS business2016\n",
    "    FROM temp_weight_join GROUP BY temp_weight_join.hexid\n",
    ");\n",
    "\n",
    "-- Create a density column                                                                 \n",
    "                                                                  \n",
    "ALTER TABLE out_data_hex_business2016 ADD COLUMN businessdensity2016 double precision;\n",
    "UPDATE out_data_hex_business2016 SET businessdensity2016 = business2016 / (34641.0161513719 / (1000 * 1000));\n",
    "                                                                                               \n",
    "                                                                                               -- lets write that to a file\n",
    "\n",
    "\\COPY out_data_hex_business2016 TO 'out_data_hex_business2016.csv' WITH (FORMAT CSV, HEADER);\n",
    "\n",
    "-- we can delete the temp tables if we want, but they take a little while to compute so it may be nice to keep them if space isnt an issue\n",
    "\n",
    "DROP TABLE temp_DA_business;\n",
    "DROP TABLE temp_int_DA11_hex;\n",
    "\n",
    "```"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.6.6"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}