taxi_clean.ipynb

{
 "cells": [
  {
   "cell_type": "raw",
   "metadata": {},
   "source": [
    "!mkdir -p taxi_raw\n",
    "!mkdir -p data/taxi_clean"
   ]
  },
  {
   "cell_type": "raw",
   "metadata": {},
   "source": [
    "!wget -nc -O data/taxi_raw/yellow_201601.csv https://s3.amazonaws.com/nyc-tlc/trip+data/yellow_tripdata_2016-01.csv\n",
    "!wget -nc -O data/taxi_raw/yellow_201602.csv https://s3.amazonaws.com/nyc-tlc/trip+data/yellow_tripdata_2016-02.csv\n",
    "!wget -nc -O data/taxi_raw/yellow_201603.csv https://s3.amazonaws.com/nyc-tlc/trip+data/yellow_tripdata_2016-03.csv\n",
    "!wget -nc -O data/taxi_raw/yellow_201604.csv https://s3.amazonaws.com/nyc-tlc/trip+data/yellow_tripdata_2016-04.csv\n",
    "!wget -nc -O data/taxi_raw/yellow_201605.csv https://s3.amazonaws.com/nyc-tlc/trip+data/yellow_tripdata_2016-05.csv"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import os\n",
    "\n",
    "import geopandas as gpd\n",
    "import pandas as pd\n",
    "pd.options.display.max_columns = None\n",
    "import matplotlib\n",
    "%matplotlib inline\n",
    "import shapely.geometry\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "I am densified (external_values, 1 elements)\n"
     ]
    }
   ],
   "source": [
    "cols = [\"tpep_pickup_datetime\", \"tpep_dropoff_datetime\", \n",
    "        \"pickup_longitude\", \"pickup_latitude\", \n",
    "        \"dropoff_longitude\", \"dropoff_latitude\",\n",
    "        \"trip_distance\"]\n",
    "boroughs = gpd.read_file(\"data/nybb_18a/nybb.shp\")\n",
    "manhattan = boroughs[boroughs.BoroName == \"Manhattan\"].to_crs(epsg=4326).reset_index(drop=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "def process_month(month):\n",
    "    print(\"processing month {}\".format(month))\n",
    "    raw_fname = \"data/taxi_raw/yellow_2016{}.csv\".format(\"0\" + str(month))\n",
    "    df = pd.read_csv(raw_fname)[cols]\n",
    "    print(\"raw df has {} rows\".format(len(df)))\n",
    "    before_5am = df[df[\"tpep_pickup_datetime\"].str.split(\" \", 1).apply(lambda x: int(x[1].split(\":\")[0])) < 5]\n",
    "    before_5am[\"pickup_point\"] = before_5am[[\"pickup_longitude\", \"pickup_latitude\"]].apply(shapely.geometry.Point, axis=1)\n",
    "    before_5am[\"dropoff_point\"] = before_5am[[\"dropoff_longitude\", \"dropoff_latitude\"]].apply(shapely.geometry.Point, axis=1)\n",
    "    before_5am = before_5am.reset_index(drop=True)\n",
    "    print(\"before 5am has {} rows\".format(len(before_5am)))\n",
    "    before_5am[\"geometry\"] = before_5am[\"pickup_point\"]\n",
    "    before_5am = gpd.GeoDataFrame(before_5am)\n",
    "    #before_5am = before_5am.set_geometry(before_5am[\"pickup_point\"], crs={\"epsg\": 4326})\n",
    "    before_filter_pickup = gpd.sjoin(before_5am, manhattan, how='inner', op='within')\n",
    "    del before_filter_pickup[\"index_right\"]\n",
    "    before_filter_pickup[\"geometry\"] = before_filter_pickup[\"dropoff_point\"]\n",
    "    before_filter_pickup.set_geometry(\"geometry\", \n",
    "                                      crs={\"epsg\": 4326}, inplace=True)\n",
    "    before_filter_dropoff = gpd.sjoin(before_filter_pickup, manhattan, how='inner', op='within')\n",
    "    before_filter_dropoff = before_filter_dropoff[cols]\n",
    "    before_filter_dropoff.to_pickle(\"data/taxi_clean/2016{}_filtered_cython.pkl\".format(\"0\" + str(month)))\n",
    "    print(\"filtered df has {} rows\".format(len(before_filter_dropoff)))\n",
    "    del before_filter_dropoff\n",
    "    del before_filter_pickup\n",
    "    del before_5am\n",
    "    del df\n",
    "    print(\"======\")\n",
    "    return \"OK\""
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "processing month 2\n",
      "processing month 1\n",
      "raw df has 10906858 rows\n",
      "raw df has 11382049 rows\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/home/christian/.miniconda2/envs/new_gpd_3/lib/python3.6/site-packages/ipykernel_launcher.py:7: SettingWithCopyWarning: \n",
      "A value is trying to be set on a copy of a slice from a DataFrame.\n",
      "Try using .loc[row_indexer,col_indexer] = value instead\n",
      "\n",
      "See the caveats in the documentation: http://pandas.pydata.org/pandas-docs/stable/indexing.html#indexing-view-versus-copy\n",
      "  import sys\n",
      "/home/christian/.miniconda2/envs/new_gpd_3/lib/python3.6/site-packages/ipykernel_launcher.py:7: SettingWithCopyWarning: \n",
      "A value is trying to be set on a copy of a slice from a DataFrame.\n",
      "Try using .loc[row_indexer,col_indexer] = value instead\n",
      "\n",
      "See the caveats in the documentation: http://pandas.pydata.org/pandas-docs/stable/indexing.html#indexing-view-versus-copy\n",
      "  import sys\n",
      "/home/christian/.miniconda2/envs/new_gpd_3/lib/python3.6/site-packages/ipykernel_launcher.py:8: SettingWithCopyWarning: \n",
      "A value is trying to be set on a copy of a slice from a DataFrame.\n",
      "Try using .loc[row_indexer,col_indexer] = value instead\n",
      "\n",
      "See the caveats in the documentation: http://pandas.pydata.org/pandas-docs/stable/indexing.html#indexing-view-versus-copy\n",
      "  \n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "before 5am has 1158880 rows\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/home/christian/.miniconda2/envs/new_gpd_3/lib/python3.6/site-packages/geopandas/tools/sjoin.py:53: UserWarning: CRS of frames being joined does not match!\n",
      "  warn('CRS of frames being joined does not match!')\n",
      "/home/christian/.miniconda2/envs/new_gpd_3/lib/python3.6/site-packages/ipykernel_launcher.py:8: SettingWithCopyWarning: \n",
      "A value is trying to be set on a copy of a slice from a DataFrame.\n",
      "Try using .loc[row_indexer,col_indexer] = value instead\n",
      "\n",
      "See the caveats in the documentation: http://pandas.pydata.org/pandas-docs/stable/indexing.html#indexing-view-versus-copy\n",
      "  \n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "before 5am has 1220310 rows\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/home/christian/.miniconda2/envs/new_gpd_3/lib/python3.6/site-packages/geopandas/tools/sjoin.py:53: UserWarning: CRS of frames being joined does not match!\n",
      "  warn('CRS of frames being joined does not match!')\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "filtered df has 821880 rows\n",
      "======\n",
      "processing month 3\n",
      "filtered df has 867801 rows\n",
      "======\n",
      "processing month 4\n",
      "raw df has 12210952 rows\n",
      "raw df has 11934338 rows\n",
      "before 5am has 1245046 rows\n",
      "filtered df has 871515 rows\n",
      "before 5am has 1294082 rows\n",
      "======\n",
      "processing month 5\n",
      "filtered df has 922057 rows\n",
      "======\n",
      "raw df has 11836853 rows\n",
      "before 5am has 1292115 rows\n",
      "filtered df has 913739 rows\n",
      "======\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "['OK', 'OK', 'OK', 'OK', 'OK']"
      ]
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "import joblib as jl\n",
    "p = jl.Parallel(n_jobs=2)\n",
    "gen = (jl.delayed(process_month)(month) for month in range(1, 6))\n",
    "p(gen)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Part 2: restart notebook before running this! \n",
    "\n",
    "## WARNING: Requires 13GB memory"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "ename": "NameError",
     "evalue": "name 'nx' is not defined",
     "output_type": "error",
     "traceback": [
      "\u001b[0;31m-------------------------------------------------------\u001b[0m",
      "\u001b[0;31mNameError\u001b[0m             Traceback (most recent call last)",
      "\u001b[0;32m<ipython-input-2-1ecd7d5f7d4b>\u001b[0m in \u001b[0;36m<module>\u001b[0;34m()\u001b[0m\n\u001b[0;32m----> 1\u001b[0;31m \u001b[0mnx\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0m__version__\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m",
      "\u001b[0;31mNameError\u001b[0m: name 'nx' is not defined"
     ]
    }
   ],
   "source": [
    "nx.__version__"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "import pandas as pd\n",
    "import geopandas as gpd\n",
    "import shapely.geometry\n",
    "from shapely.ops import nearest_points\n",
    "import networkx as nx\n",
    "\n",
    "import glob\n",
    "import numpy as np\n",
    "\n",
    "cols = [\"tpep_pickup_datetime\", \"tpep_dropoff_datetime\", \n",
    "        \"pickup_longitude\", \"pickup_latitude\", \n",
    "        \"dropoff_longitude\", \"dropoff_latitude\",\n",
    "        \"trip_distance\"]\n",
    "\n",
    "node_cols = [\"NODEID\", \"geometry\"]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "I am densified (external_values, 131335 elements)\n",
      "131335\n",
      "12263\n"
     ]
    }
   ],
   "source": [
    "files = glob.glob(\"data/taxi_clean/*_filtered_cython.pkl\")\n",
    "nodes = gpd.read_file(\"data/lion/lion.shp/node.shp\").to_crs(epsg=4326)[node_cols]\n",
    "def uniform_str(x):\n",
    "    strd = str(x)\n",
    "    while len(strd) < 7:\n",
    "        strd = '0' + strd\n",
    "    return strd\n",
    "\n",
    "nodes[\"NODEID_STR\"] = nodes[\"NODEID\"].apply(uniform_str)\n",
    "nodes.set_index(\"NODEID_STR\", inplace=True)\n",
    "\n",
    "g = set(nx.read_gpickle(\"data/final_graph_1st_pass_nx_2.1.pkl\").nodes())\n",
    "for x in g:\n",
    "    if not x in nodes.index:\n",
    "        print(x)\n",
    "print(len(nodes))\n",
    "nodes = nodes.loc[g].reset_index()[node_cols]\n",
    "print(len(nodes))\n",
    "del g"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [],
   "source": [
    "i = 0\n",
    "def closest_p(g):\n",
    "    global i\n",
    "    i = i + 1\n",
    "    if i % 100000 == 0:\n",
    "        print(i)\n",
    "    if len(g) == 1 or i > 10:\n",
    "        return g.iloc[0]\n",
    "  \n",
    "    s = g.pickup_point.iloc[0]\n",
    "    nod = nodes.loc[g.index_right].geometry\n",
    "    distances = [s.distance(s2) for s2 in nod]\n",
    "    return g.iloc[np.argmin(distances)]\n",
    "\n",
    "def closest_d(g):\n",
    "    global i\n",
    "    i = i + 1\n",
    "    if i % 100000 == 0:\n",
    "        print(i)\n",
    "    if len(g) == 1 or i > 10:\n",
    "        return g.iloc[0]\n",
    "  \n",
    "    s = g.dropoff_point.iloc[0]\n",
    "    nod = nodes.loc[g.index_right].geometry\n",
    "    distances = [s.distance(s2) for s2 in nod]\n",
    "    return g.iloc[np.argmin(distances)]\n",
    "\n",
    "def handle_file(f):\n",
    "    global i \n",
    "    i = 0\n",
    "    out = f.split(\".\")[0] + \"_od_v2.pkl\"\n",
    "    p = pd.read_pickle(f)\n",
    "    orig = len(p)\n",
    "    p[\"pickup_point\"] = p[[\"pickup_longitude\", \"pickup_latitude\"]].apply(shapely.geometry.Point, axis=1)\n",
    "    p[\"geometry\"] = p[\"pickup_point\"].apply(lambda x: x.buffer(.0005))\n",
    "    p = gpd.GeoDataFrame(p)\n",
    "    joined = gpd.sjoin(p, nodes, how='inner', op='contains')\n",
    "    del p\n",
    "    closest_only = joined.groupby(joined.index).apply(closest_p)\n",
    "    closest_only[\"NODEID_O\"] = closest_only[\"NODEID\"]\n",
    "\n",
    "    del joined\n",
    "    del closest_only[\"NODEID\"]\n",
    "    del closest_only[\"pickup_point\"]\n",
    "    closest_only = closest_only[cols + [\"NODEID_O\"]]\n",
    "    closest_only[\"dropoff_point\"] = closest_only[[\"dropoff_longitude\", \"dropoff_latitude\"]].apply(shapely.geometry.Point, axis=1)\n",
    "    closest_only[\"geometry\"] = closest_only[\"dropoff_point\"].apply(lambda x: x.buffer(.0005))\n",
    "    closest_only = gpd.GeoDataFrame(closest_only)\n",
    "    joined = gpd.sjoin(closest_only, nodes, how='inner', op='contains')\n",
    "    \n",
    "    i = 0\n",
    "    final = joined.groupby(joined.index).apply(closest_d)\n",
    "    del joined\n",
    "    del closest_only\n",
    "    final[\"NODEID_D\"] = final[\"NODEID\"]\n",
    "    final = final[cols + [\"NODEID_O\", \"NODEID_D\"]]\n",
    "    final.to_pickle(out)\n",
    "    \n",
    "    print(\"orig had {}\".format(orig))\n",
    "    print(\"final had {}\".format(len(final)))\n",
    "    print(\"{} %\".format(orig / float(len(final))))\n",
    "    return final"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "handling data/taxi_clean/201601_filtered_cython.pkl\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/home/christian/.miniconda2/envs/new_gpd_3/lib/python3.6/site-packages/geopandas/tools/sjoin.py:53: UserWarning: CRS of frames being joined does not match!\n",
      "  warn('CRS of frames being joined does not match!')\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "100000\n",
      "200000\n",
      "300000\n",
      "400000\n",
      "500000\n",
      "600000\n",
      "700000\n",
      "100000\n",
      "200000\n",
      "300000\n",
      "400000\n",
      "500000\n",
      "600000\n",
      "orig had 867801\n",
      "final had 604965\n",
      "1.4344648037489771 %\n",
      "===\n",
      "handling data/taxi_clean/201603_filtered_cython.pkl\n",
      "100000\n",
      "200000\n",
      "300000\n",
      "400000\n",
      "500000\n",
      "600000\n",
      "700000\n",
      "100000\n",
      "200000\n",
      "300000\n",
      "400000\n",
      "500000\n",
      "600000\n",
      "orig had 871515\n",
      "final had 604488\n",
      "1.44174077897328 %\n",
      "===\n",
      "handling data/taxi_clean/201605_filtered_cython.pkl\n",
      "100000\n",
      "200000\n",
      "300000\n",
      "400000\n",
      "500000\n",
      "600000\n",
      "700000\n",
      "100000\n",
      "200000\n",
      "300000\n",
      "400000\n",
      "500000\n",
      "600000\n",
      "orig had 913739\n",
      "final had 639095\n",
      "1.4297389277024541 %\n",
      "===\n",
      "handling data/taxi_clean/201604_filtered_cython.pkl\n",
      "100000\n",
      "200000\n",
      "300000\n",
      "400000\n",
      "500000\n",
      "600000\n",
      "700000\n",
      "100000\n",
      "200000\n",
      "300000\n",
      "400000\n",
      "500000\n",
      "600000\n",
      "orig had 922057\n",
      "final had 638255\n",
      "1.4446529991931125 %\n",
      "===\n",
      "handling data/taxi_clean/201602_filtered_cython.pkl\n",
      "100000\n",
      "200000\n",
      "300000\n",
      "400000\n",
      "500000\n",
      "600000\n",
      "700000\n",
      "100000\n",
      "200000\n",
      "300000\n",
      "400000\n",
      "500000\n",
      "orig had 821880\n",
      "final had 566836\n",
      "1.4499431934457232 %\n",
      "===\n"
     ]
    }
   ],
   "source": [
    "for f in files:\n",
    "    print(\"handling {}\".format(f))\n",
    "    out_thing = handle_file(f)\n",
    "    print(\"===\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "!ls -lh data/taxi_clean/*_od.pkl"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "joined = gpd.sjoin(p, nodes, how='inner', op='contains')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "del joined[\"VIntersect\"]\n",
    "del joined[\"index_righ\"]\n",
    "del joined[\"BoroCode\"]\n",
    "del joined[\"BoroName\"]\n",
    "del joined[\"Shape_Leng\"]\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "closest_only = joined.groupby(joined.index).apply(closest)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print len(closest_only) / float(len(p))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "closest_only.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# pts = nodes.geometry.unary_union\n",
    "# i = 0\n",
    "# def near(row, pts=pts):\n",
    "#     global i\n",
    "#     i = i + 1\n",
    "#     if i % 100 == 0:\n",
    "#         print \"{} out of {}\".format(i, len(p))\n",
    "#     point = row.geometry\n",
    "#     n = nearest_points(point, pts)\n",
    "#     nearest = nodes.geometry == n[0]\n",
    "#     return nodes[nearest].NODEID.get_values()[0]\n",
    "# p[\"nearest_origin\"] = p.apply(lambda row: near(row), axis=1)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Part 3: Time partitioning"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "import glob\n",
    "import pandas as pd\n",
    "final_files = glob.glob(\"data/taxi_clean/*od_v2.pkl\")\n",
    "all_dfs = pd.concat(pd.read_pickle(f) for f in final_files)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {},
   "outputs": [],
   "source": [
    "all_dfs[\"pickup_hr\"] = all_dfs[\"tpep_pickup_datetime\"].str.split(\" \").str[-1].str.split(\":\").str[0].apply(int)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 26,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0, 1062644\n",
      "1, 769259\n",
      "2, 556066\n",
      "3, 396180\n",
      "4, 269490\n"
     ]
    }
   ],
   "source": [
    "for name, grouped in all_dfs.groupby(\"pickup_hr\"):\n",
    "    print(\"{}, {}\".format(name, len(grouped)))\n",
    "    grouped.to_pickle(\"data/taxi_clean/FINAL_HR_{}.pkl\".format(name))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.6.5"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}