4k/2014 Rate Calculation using CART.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Cross-Classification Refinement \n",
    "For 2014 Survey data for 4k\n",
    "\n",
    "Using a decision (classification) tree approach to update trip rates\n",
    "Based on guidance from this report: http://citeseerx.ist.psu.edu/viewdoc/download?doi=10.1.1.539.3755&rep=rep1&type=pdf\n",
    "\n",
    "Using python tools to build a decision tree for generating classifications\n",
    "This process helps combine classes with small samle sizes to match existing crossclass definitions"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "import pandas as pd\n",
    "import numpy as np\n",
    "from sklearn.tree import DecisionTreeClassifier, export_graphviz\n",
    "from sklearn.model_selection import train_test_split"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# Load household records to prepare data and generate cross-classification\n",
    "hh = pd.read_excel(r'J:\\Projects\\Surveys\\HHTravel\\Survey2014\\Data\\Final database\\Release 4\\2014-pr3-M-hhsurvey-households.xlsx',\n",
    "                  sheetname='Data')\n",
    "\n",
    "# Reclassify household columns for cross-classification\n",
    "\n",
    "# Number of Workers 0 -> 3+\n",
    "hh['numworkers_crossclass'] = hh['numworkers']\n",
    "hh.ix[hh['numworkers'] >= 3, 'numworkers_crossclass'] = 3\n",
    "hh['numworkers_crossclass'] = hh['numworkers_crossclass'].astype('int')\n",
    "\n",
    "# Household size\n",
    "hh['hhsize_crossclass'] = hh['hhsize']\n",
    "hh.ix[hh['hhsize'] >= 4, 'hhsize_crossclass'] = 4\n",
    "hh['hhsize_crossclass'] = hh['hhsize_crossclass'].astype('int')\n",
    "\n",
    "# Household income\n",
    "# Note that exact ranges from 2014 do not match 2006\n",
    "# old ranges: 30, 60, 90+\n",
    "# new ranges: 35, 75, 100+\n",
    "hh['income_crossclass'] = hh['hh_income_detailed_imp']\n",
    "hh.ix[hh['hh_income_detailed_imp'] <= 3, 'income_crossclass'] = 1    # $35k  \n",
    "hh.ix[(hh['hh_income_detailed_imp'] > 3) & (hh['hh_income_detailed_imp'] <= 5),    # $35-75k\n",
    "           'income_crossclass'] = 2\n",
    "hh.ix[(hh['hh_income_detailed_imp'] > 5) & (hh['hh_income_detailed_imp'] <= 6),    # $75-100k\n",
    "           'income_crossclass'] = 3\n",
    "hh.ix[hh['hh_income_detailed_imp'] >= 7, 'income_crossclass'] = 4    # >$100k"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# Data needs to be in format of each household on a row and columns for trips by purpose\n",
    "# Load trip data and compute 4k purposes\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# Data needs to be in format of each household on a row and columns for trips by purpose\n",
    "\n",
    "# Load trip data and compute 4k purposes\n",
    "\n",
    "######\n",
    "# NOTE: different bins will be created based on which trip file is used\n",
    "# For now, we are using the GPS trip weights to create the bins\n",
    "# and using this one throughout\n",
    "######\n",
    "# 2014 File with adjustments from 2017 smartphone survey\n",
    "trip = pd.read_csv(r'J:\\Projects\\Surveys\\HHTravel\\Survey2014\\Data\\Final database\\Release 4\\Adjusted\\trips_2014_adjusted_wt_daysim.csv')\n",
    "\n",
    "# GPS Weighted trips (in Daysim format)\n",
    "# trip = pd.read_csv(r'R:\\SoundCastDocuments\\2014Estimation\\Files_From_Mark_2014\\gps_weights_11_3_16\\formatted\\skims_attached\\tripP14_w.dat')\n",
    "\n",
    "# Non GPS weighted\n",
    "# trip = pd.read_csv(r'R:\\SoundCast\\estimation\\2014\\Surveys\\2014 Survey\\P5\\tripP5.dat', delim_whitespace=True)\n",
    "\n",
    "# Separate college student trips from regional survey trips\n",
    "college_trips = trip[trip['hhno'] < 14000000]\n",
    "trip = trip[trip['hhno'] >= 14000000]\n",
    "\n",
    "# Home-Based Work (HBW) Trips, directly from home-to-work and work-to-home\n",
    "trip.ix[(trip['opurp'] == 0) & (trip['dpurp'] == 1), '4k_purp'] = 'HBW'\n",
    "trip.ix[(trip['opurp'] == 1) & (trip['dpurp'] == 0), '4k_purp'] = 'HBW'\n",
    "\n",
    "# Home-Based shopping\n",
    "trip.ix[(trip['opurp'] == 0) & (trip['dpurp'] == 5), '4k_purp'] = 'Home-Based Shopping'\n",
    "trip.ix[(trip['opurp'] == 5) & (trip['dpurp'] == 0), '4k_purp'] = 'Home-Based Shopping'\n",
    "\n",
    "# Home-Based School\n",
    "trip.ix[(trip['opurp'] == 0) & (trip['dpurp'] == 2), '4k_purp'] = 'School'\n",
    "trip.ix[(trip['opurp'] == 2) & (trip['dpurp'] == 0), '4k_purp'] = 'School'\n",
    "\n",
    "# Home-Based Other\n",
    "trip.ix[(trip['opurp'] == 0) & (trip['dpurp'].isin([3,4,6,7,8,9,10])), '4k_purp'] = 'Home-Based Other'\n",
    "trip.ix[(trip['opurp'].isin([3,4,6,7,8,9,10]) & (trip['dpurp'] == 0)), '4k_purp'] = 'Home-Based Other'\n",
    "\n",
    "# NHB Work-to-Other\n",
    "trip.ix[(trip['opurp'] == 1) & (trip['dpurp'] != 0) & (trip['4k_purp'].isnull()), '4k_purp'] = 'NHB WtO'\n",
    "trip.ix[(trip['opurp'] != 0) & (trip['dpurp'] == 1) & (trip['4k_purp'].isnull()), '4k_purp'] = 'NHB WtO'\n",
    "\n",
    "# # NHB Other-to-Other (Destination and Origins are neither work nor home)\n",
    "trip.ix[(trip['opurp'] != 1) & (trip['dpurp'] != 0) & (trip['opurp'] != 0) \n",
    "        & (trip['dpurp'] != 1) & (trip['4k_purp'].isnull()), '4k_purp'] = 'NHB OtO'"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# Export GPS Trips with 4k_purp\n",
    "# trip.to_excel(r'R:\\4K\\2014\\Trip Generation\\Trip Rates\\trip_2014_survey_adjusted.xlsx', index=False)\n",
    "\n",
    "# Add household fields to data\n",
    "# hh_fields = ['vehicle_count','hhsize','numadults','numchildren','numworkers','lifecycle',\n",
    "#              'rent_own','res_type','income_crossclass','h_city','h_address','h_county_name','h_district_name',\n",
    "#              'h_taz2010','hhid']\n",
    "# export_trip_hh = pd.merge(trip, hh[hh_fields], left_on='hhno', right_on='hhid', how='left')\n",
    "# export_trip_hh.to_csv(r'R:\\4K\\2014\\Trip Generation\\Trip Rates\\trip_2014_survey_gps.csv', index=False)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# Apply a decision tree to figure out which rates should be averaged and combined\n",
    "\n",
    "min_samples_size = 50\n",
    "seed_value = 99\n",
    "int_conversion = 1000\n",
    "\n",
    "# Reorganize the trip record data\n",
    "# Group by household id\n",
    "_df = trip.groupby(['hhno','4k_purp']).sum()[['trexpfac']]\n",
    "_df = _df.unstack()\n",
    "_df.fillna(0,inplace=True)\n",
    "_df = _df.reset_index()\n",
    "purp_cols = ['HBW','Home-Based Other','Home-Based Shopping','NHB OtO','NHB WtO','Home-Based School']\n",
    "_df.columns = ['hhno']+purp_cols\n",
    "\n",
    "# # Join hh information \n",
    "df = pd.merge(_df, hh[['hhid','hhsize_crossclass','numworkers_crossclass', 'income_crossclass']], \n",
    "         left_on='hhno', right_on='hhid',how='left')\n",
    "\n",
    "# Get total trips per household\n",
    "df['tot_trips'] = df['HBW']+df['Home-Based Other']+df['Home-Based Shopping']+df['NHB OtO']+df['NHB WtO']+df['Home-Based School']\n",
    "\n",
    "# For decision tree analysis, need a target (y) and features (x)\n",
    "# We are modeling trips as a function of HH size, number of workers, and income classes\n",
    "\n",
    "# Features for classification\n",
    "features = [\"hhsize_crossclass\",\"numworkers_crossclass\",\"income_crossclass\"]\n",
    "X = df[features]\n",
    "\n",
    "# Decision trees are used to predict total trips\n",
    "# since we are classifying based on household characteristics and not trips\n",
    "# we can use total trips as dependent variable to generate bins\n",
    "\n",
    "# Must be in integer form for prediction\n",
    "# Convert to int by multiplying by factor\n",
    "df['tot_trips_int'] = (df['tot_trips']*int_conversion).astype('int')\n",
    "\n",
    "# Define the target y value as trips by purpsoe\n",
    "y = df['tot_trips_int']\n",
    "\n",
    "# Create training and testing sample (maybe not necessary...)\n",
    "X_train, X_test, y_train, y_test = train_test_split(X, y, random_state=0)\n",
    "\n",
    "# Create decision tree, set min sample size and random seed\n",
    "# This part generates the optimal category separations\n",
    "estimator = DecisionTreeClassifier(min_samples_leaf=min_samples_size, random_state=seed_value)\n",
    "\n",
    "# Export for visualization\n",
    "# visualize_tree(estimator, features)\n",
    "\n",
    "# Use the new estimator object to assign categories to original data\n",
    "estimator.fit(X_train, y_train)\n",
    "new_categories = estimator.apply(X)\n",
    "df['optimum_grouping'] = new_categories"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 29,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# Group by old and new categories to apply new cross class blending\n",
    "# For rows that have the same value for 'new_grouping' col, results for these should be averaged\n",
    "df_cat = df.groupby(['optimum_grouping','hhsize_crossclass','income_crossclass','numworkers_crossclass']).count()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "df_cat = df_cat.reset_index()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "hh = pd.merge(hh, df_cat[['hhsize_crossclass','income_crossclass','numworkers_crossclass','optimum_grouping']], \n",
    "         how='left', on=['hhsize_crossclass','income_crossclass','numworkers_crossclass'])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# Total households by cross-classification group (hhsize, numworkers, and income)\n",
    "total_hh = hh.groupby(['optimum_grouping']).sum()[['expwt_2']]\n",
    "total_hh = total_hh.reset_index()\n",
    "\n",
    "# Join the category values\n",
    "total_hh = pd.merge(total_hh, df_cat[['hhsize_crossclass','numworkers_crossclass','income_crossclass','optimum_grouping']], \n",
    "                    how='left', on='optimum_grouping')\n",
    "total_hh.to_clipboard()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 30,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "hh.groupby(['optimum_grouping']).count()[['expwt_2']].to_clipboard()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Calculate trips rates\n",
    "\n",
    "# Load trip files\n",
    "# GPS Weighted trips (in Daysim format)\n",
    "\n",
    "\n",
    "# Merge household grouping to trip\n",
    "trip_hh = pd.merge(trip, hh, left_on='hhno', right_on='hhid', how='left')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "hhno\n",
      "pno\n",
      "day\n",
      "tour\n",
      "half\n",
      "tseg\n",
      "tsvid\n",
      "opurp\n",
      "dpurp\n",
      "oadtyp\n",
      "dadtyp\n",
      "opcl\n",
      "otaz\n",
      "dpcl\n",
      "dtaz\n",
      "mode\n",
      "pathtype\n",
      "dorp\n",
      "deptm\n",
      "arrtm\n",
      "endacttm\n",
      "trexpfac\n",
      "id\n",
      "travcost\n",
      "travdist\n",
      "travtime\n",
      "tripID\n",
      "4k_purp\n"
     ]
    }
   ],
   "source": [
    "for i in trip.columns:\n",
    "    print i"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 28,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Trips by Purpose\n",
    "# purp = 'HBW'\n",
    "# purp = 'Home-Based Shopping'\n",
    "# purp = 'Home-Based Other'\n",
    "# purp = 'School'\n",
    "purp = 'NHB OtO'\n",
    "# purp = 'NHB WtO'\n",
    "_df = trip_hh[trip_hh['4k_purp'] == purp].groupby(['optimum_grouping']).sum()[['trexpfac']]\n",
    "_df = _df.reset_index()\n",
    "_df = pd.merge(_df, df_cat[['hhsize_crossclass','numworkers_crossclass','income_crossclass','optimum_grouping']], \n",
    "                    how='left', on='optimum_grouping')\n",
    "_df.to_clipboard()\n",
    "\n",
    "# Get trip totals by purpose\n",
    "# trip_hh[trip_hh['4k_purp'] == purp]['trexpfac'].sum()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 39,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "4k_purp\n",
       "HBW                    2.339437e+06\n",
       "Home-Based Other       3.944246e+06\n",
       "Home-Based Shopping    1.383019e+06\n",
       "NHB OtO                2.504022e+06\n",
       "NHB WtO                1.270762e+06\n",
       "School                 1.291760e+06\n",
       "Name: trip_wt_final, dtype: float64"
      ]
     },
     "execution_count": 39,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Get trip totals\n",
    "trip_hh.groupby('4k_purp').sum()['trexpfac']"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 65,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "12960491.149999935"
      ]
     },
     "execution_count": 65,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "trip_hh.groupby('4k_purp').sum()['trexpfac'].sum()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 2",
   "language": "python",
   "name": "python2"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 2
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython2",
   "version": "2.7.14"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}