gen_population_density.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Population Density"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Joining population data from census blocks to a hex grid via areal interpolation\n",
    "\n",
    "Areal interpolation is the process of transfering data from one set of areal units to another set of areal units\n",
    "\n",
    "In this case, we are transfering population data from census blocks $b$ to a hex cell $h$ based on the ratio of the area of the block that intersects the hex cell $A_{b \\cap h}$ by the total area of the block $A_b$\n",
    "\n",
    "### $P_h = \\sum P_{b} {w_{b,h}}$\n",
    "\n",
    "### $P_h = \\sum P_{b} \\frac{A_{b \\cap h}}{A_b}$\n"
   ]
  },
  {
   "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",
    "\n",
    "-- creating an area column $A_b$ for the block data\n",
    "\n",
    "ALTER TABLE zones_db ADD COLUMN area double precision;\n",
    "UPDATE zones_db SET area = ST_AREA(geom);\n",
    "\n",
    "-- Then we join the block populationo data to the spatial boundaries, $P_b$\n",
    "\n",
    "CREATE TABLE temp_DB_pop AS ( \n",
    "    SELECT\n",
    "    zones_DB.geom AS geom,\n",
    "    zones_DB.dbuid AS dbuid,\n",
    "    zones_DB.area AS area,\n",
    "    table_DB_pop_2016.pop2016int AS pop2016int\n",
    "    FROM zones_DB\n",
    "    INNER JOIN table_DB_pop_2016 ON table_DB_pop_2016.dbuid = zones_DB.dbuid\n",
    ");    \n",
    "\n",
    "-- create a spatial index for this block polygon file\n",
    "\n",
    "DROP INDEX IF EXISTS temp_DB_pop_gix;\n",
    "CREATE INDEX temp_DB_pop_gix ON temp_DB_pop 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_DB_hex AS\n",
    "(\n",
    "SELECT\n",
    "hex_grid_200m.id AS hexid,\n",
    "temp_DB_pop.dbuid AS dbuid,\n",
    "temp_DB_pop.area AS area_full,\n",
    "ST_Intersection(hex_grid_200m.geom,temp_DB_pop.geom) AS geom\n",
    "FROM\n",
    "hex_grid_200m INNER JOIN temp_DB_pop ON ST_Intersects(hex_grid_200m.geom,temp_DB_pop.geom)\n",
    ");\n",
    "\n",
    "-- update area of intersected geoms\n",
    "\n",
    "ALTER TABLE temp_int_DB_hex ADD COLUMN area_int double precision;\n",
    "UPDATE temp_int_DB_hex SET area_int = ST_AREA(geom);\n",
    "\n",
    "ALTER TABLE temp_int_DB_hex ADD COLUMN area_ratio double precision;\n",
    "UPDATE temp_int_DB_hex SET area_ratio = area_int / area_full;\n",
    "\n",
    "-- grouping by hex and DB unique IDs - i.e. this is a weights table that can be used for apportioning data\n",
    "\n",
    "CREATE TABLE weights_db_hex AS\n",
    "(\n",
    "SELECT\n",
    "hexid,\n",
    "dbuid,\n",
    "sum(area_ratio) AS weight\n",
    "FROM\n",
    "temp_int_DB_hex\n",
    "GROUP BY dbuid, hexid\n",
    "ORDER BY dbuid, hexid\n",
    ");\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_DB_pop;\n",
    "DROP TABLE temp_int_DB_hex;\n",
    "\n",
    "```\n",
    "\n",
    "We now use the weights table to apportion the population data to the hex grid\n",
    "\n",
    "```sql\n",
    "-- joining the population to the weights table\n",
    "-- multiplying the pop by the weight and aggregating by hex grid - also normalize as persons/km2\n",
    "\n",
    "DROP TABLE IF EXISTS out_data_hex_pop2016;\n",
    "CREATE TABLE out_data_hex_pop2016 AS (\n",
    "    WITH temp_weight_join AS (\n",
    "        SELECT\n",
    "        weights_db_hex.dbuid AS dbuid,\n",
    "        weights_db_hex.hexid AS hexid,\n",
    "        weights_db_hex.weight AS weight,\n",
    "        table_DB_pop_2016.pop2016int AS pop2016int\n",
    "        FROM weights_db_hex\n",
    "        INNER JOIN table_DB_pop_2016 ON table_DB_pop_2016.dbuid = weights_db_hex.dbuid\n",
    "    ) \n",
    "    SELECT\n",
    "    temp_weight_join.hexid AS hexid,\n",
    "    SUM(temp_weight_join.weight * temp_weight_join.pop2016int) AS pop2016\n",
    "    FROM temp_weight_join GROUP BY temp_weight_join.hexid\n",
    ");\n",
    "  \n",
    "-- Create a pop density column                                                                 \n",
    "                                                                  \n",
    "ALTER TABLE out_data_hex_pop2016 ADD COLUMN popdensity2016 double precision;\n",
    "UPDATE out_data_hex_pop2016 SET popdensity2016 = pop2016 / (34641.0161513719 / (1000 * 1000));\n",
    "                                                                  \n",
    "-- lets write that to a file\n",
    "\n",
    "\\COPY out_data_hex_pop2016 TO 'output_hex_data/out_data_hex_pop2016.csv' WITH (FORMAT CSV, HEADER);\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
}