# ---
# jupyter:
#   jupytext:
#     text_representation:
#       extension: .py
#       format_name: light
#       format_version: '1.5'
#       jupytext_version: 1.14.4
#   kernelspec:
#     display_name: Python 3
#     language: python
#     name: python3
# ---

# %load_ext autoreload
# %autoreload 2

# import libraries
import pandas as pd
pd.set_option('display.max_columns', None)
import matplotlib.pyplot as plt
# %matplotlib inline
import seaborn as sns
import numpy as np
from datetime import timedelta, datetime
from termcolor import colored

# aims:
# 1. extract features from test failing vehicles
# 2. build classification model for predicting outcome of emission test for vehicles
#
# additionals:
#
# 1. new customer - prediction
# 2. return customer - prediction
# 3. predict test cost

# # import data
# [here](https://opendata.cabq.gov/dataset/vehicle-emissions/resource/2663c138-549e-4d55-932f-6aa549c7d158) is the MetaData
#

# +
df_all = pd.read_csv('../data/sample201320.csv', low_memory=False)

print("data size:", df_all.shape)
df_all.head(2)
# -

# # target
#
# 1. keep all the records with OVERALL_RESULT = P or F
# 2. change to numeric
# - P    173762
# - F     13741

import seaborn as sns
y = df_all['OVERALL_RESULT']
print('Num of rows with missing values:', y.isnull().sum())
tmp = y.value_counts()*100.0/y.shape[0]
sns.barplot(x=tmp.index, y=tmp);

y.value_counts()

df_all.loc[(df_all.OVERALL_RESULT=='A') & (df_all.ABORT_CODE_OTHER_DESC.notnull()), 'ABORT_CODE_OTHER_DESC'].head(30)

# drop the rows where the test was aborted
df_all = df_all[~df_all.OVERALL_RESULT.isin(['A','O'])]
df_all.OVERALL_RESULT.value_counts()

# P=1, F=1
df_all['RESULT'] = df_all['OVERALL_RESULT'].map({'P':1, 'F':0})
df_all.RESULT.value_counts()

# # vehicle features
#

# find the unique id column
print(df_all.shape[0])
print(df_all[['RecordID', 'RECORD_NUM', 'VIR_CERT_NUMBER']].nunique())

# set RecordID as index
df_all.set_index('RecordID', inplace=True)
print("RecordID has been set as index")

# ## quick exploration for selecting features

# + [markdown] heading_collapsed=true
# ### vehicle features and missing values
# - VID_TEST_TYPE:    80462
# - GVW_TYPE:            71
# - VEH_LICENSE:         13
# - MODEL:                2
#

# + hidden=true
cols = ['RESULT',
        'TEST_TYPE', 'VID_TEST_TYPE', # which
        'TEST_SDATE', 'TEST_EDATE', 
        'VIN', 
        'ZIP_CODE',
        'VEHICLE_TYPE',
        'MODEL_YEAR', 'MAKE', 'MODEL',
        'GVW_TYPE', 'GVWR', # which
        'CYL', 
        'ENGINE_SIZE', 
        'TRANS_TYPE',
        'DUAL_EXHAUST', 
        'ODOMETER', 
        'FUEL_TYPE']
df = df_all[cols]
df.head()

# + hidden=true
df.isnull().sum().sort_values(ascending=False)[:5]

# + hidden=true
# how many unique vehicles
df.VIN.nunique()

# + [markdown] heading_collapsed=true
# ### TEST_TYPE VS VID_TEST_TYPE
# - drop VID_TEST_TYPE due to the missing values

# + hidden=true
tmp = df[['RESULT','TEST_TYPE', 'VID_TEST_TYPE']]
fig, axs = plt.subplots(1, 2, figsize=(20, 4))
sns.countplot(data=tmp, x='RESULT', hue='TEST_TYPE', ax=axs[0])
sns.countplot(data=tmp, x='RESULT', hue='VID_TEST_TYPE', ax=axs[1])

# + hidden=true
df = df.drop(columns=['VID_TEST_TYPE'])
# -

# ### GVW_TYPE vs GVWR
# - drop GVW_TYPE
#

tmp = df[['RESULT','GVW_TYPE', 'GVWR']]
sns.countplot(data=tmp, x='RESULT', hue='GVW_TYPE')

sns.pairplot(tmp)

# drop GVW
df = df.drop(columns=['GVW_TYPE'])

# + [markdown] heading_collapsed=true
# ### categorical cols distribution
# - drop DUAL_EXHAUST

# + hidden=true
df.nunique().sort_values()

# + hidden=true
df.nunique().sort_values().index

# + hidden=true
categorical_cols = ['TEST_TYPE', 
                    'TRANS_TYPE', 
                    'VEHICLE_TYPE', 
                    'FUEL_TYPE', 
                    'CYL']

# distribution of categorical values in Pass and Fail data sets
df['count'] = 1
n = len(categorical_cols[1:])
fig, axs = plt.subplots(n, 1, figsize=(5, 6*n))
for col, ax in zip(categorical_cols[1:], axs.flat):
    tmp = df.groupby(['RESULT', col])['count'].count().groupby('RESULT').apply(lambda x: 100*x/x.sum()).unstack(level=1)
    tmp.plot(kind='bar', stacked=True, ax=ax)
    ax.set_title(col)
df = df.drop(columns=['count'])

# + hidden=true
# let's look closer to DUAL_EXHAUST
df.groupby(['RESULT'])['DUAL_EXHAUST'].value_counts()/df.shape[0]

# + hidden=true
# let's drop dual_exhaust
df = df.drop(columns=['DUAL_EXHAUST'])
# -

# ## cleaning and feature engineering
#
#
# **set RecordID as index**
#
#     - total records: 187503
#     - P: 173762
#     - F: 13741
#
# **string columns : VIN, MAKE, MODEL, ZIP_CODE, MODEL_YEAR**
#
#
#     - missing values in MODEL: 2
#         - using the historical records of that VIN fill in the missing MODEL values
#
#
#     - engineer new feature: MAKE + MODEL + MODEL_YEAR
#         - unique MAKE_MODEL_YEAR: 11684
#         - total records: 187503
#         
#         
#     - drop: MODEL, MAKE, MODEL_YEAR
#     
#     
#     - extract features from ZIP_CODE: pause!!!
#         - typos
#         - zip codes of ABQ are not structured
#     
#     
#     - drop: ZIP_CODE
#     
#     
# final columns: VIN, MAKE_MODEL_YEAR
#     
#
# **datetime columns: TEST_SDATE, TEST_EDATE**
#
#     - if a vehicle has multiple records from the same date, keep first record
#     
#         - records where vehicles got tested more than 1 time a day: 4046
#         
#         
#         - after removing repeated tests from same day
#         
#         
#             - total number of records: 187503 > 183457 
#             - P: 173762 > 170120
#             - F: 13741 > 13337
#             - data range: 2013-01-01/2020-12-31 
#             
#             
#     - engineer feature: VEHICLE_AGE   
#         
#     - drop TEST_SDATE
#     
# final columns: VEHICLE_AGE
#         
# **categorical columns: TEST_TYPE, TRANS_TYPE, VEHICLE_TYPE, FUEL_TYPE, CYL**
#
#     - map to numeric: TEST_TYPE:{'I':0, 'A':1}, TRANS_TYPE:{'A':0, 'M':1}
#
#     - drop VEHICLE_TYPE
#         - near 1000 MAKE_MODEL_YEARs had more than 1 VEHICLE_TYPE 
#         - it's not clear how vehicle types were defined     
#         
#     - drop CYL
#         - share of fail is higher when CYL is 9, 12 and R
#         - however only 82 records whose CYL is 9, 12 and R
#
#     - drop FUEL_TYPE (FUEL_TYPE L only has 4 records and all are Pass)
#         - more than 95% of the tests have FUEL_TYPE G
#
# final columns: TEST_TYPE, TRANS_TYPE
#
#
# **numeric columns: GVWR, ENGINE_SIZE, ODOMETER**
#     
#     - 0 in GVWR: 11577 (P: 11090, F: 487)
#         
#         - took the median GVWR values for each VIN
#             
#         - drop records where GVWR values are missing: 7942 (P: 7560, F: 382) 
#         
#         - range of GVWR: 847 - 10000
#         
#     - drop records where ODOMETER = 0 : 875 (P: 746, F: 129) 
#     
#     - possible outliers in ODOMETER
#         - records where ODOMETER > 1 million miles: 229
#         - vehicles with ODOMETER > 1 million miles: 218
#         
# **final columns: RESULT, TEST_TYPE, ENGINE_SIZE, TRANS_TYPE, ODOMETER, VEHICLE_AGE, GVWR**
#
#
# **final data size: 174651, 8**

# + [markdown] heading_collapsed=true
# ### missing values 
#

# + hidden=true
df.head(2)

# + hidden=true
df['RESULT'].value_counts()

# + hidden=true
df.isnull().sum().sort_values(ascending=False)[:5]

# + hidden=true
numeric_cols = ['MODEL_YEAR', 'GVWR', 'ENGINE_SIZE', 'ODOMETER']

fig, axs = plt.subplots(4, 1, figsize=(20, 12))
for col, ax in zip(numeric_cols, axs.flat):
    print(col, df[col].sort_values().values[0], df[col].sort_values().values[-1])
    sns.boxplot(x=df[col], ax=ax)

# + hidden=true
df.head()

# + [markdown] heading_collapsed=true
# ### strings

# + hidden=true
# select string columns

string_cols = ['VIN', 'ZIP_CODE', 'MODEL_YEAR', 'MAKE', 'MODEL']

# tranform string columns

def to_string(df, string_cols):
    for col in string_cols:
        df[col] = df[col].astype('string').str.strip().str.lower()
        df
    return df

df = to_string(df, string_cols)
print(df[string_cols].dtypes)
df[string_cols].head()

# + hidden=true
df[string_cols].isnull().sum()

# + hidden=true
# using the historical values for that VIN fill in the missing MODEL values
vin = df[df.MODEL.isnull()].VIN.values[0]
model = df[df.VIN==vin].MODEL.values[-1]
df.loc[df.MODEL.isnull(), 'MODEL'] = model
df.MODEL.isnull().sum()

# + hidden=true
# make + model + year
df['MAKE_MODEL_YEAR'] = df.MAKE + '/' + df.MODEL + '/' + df.MODEL_YEAR
print("unique MAKE_MODEL_YEAR:", df.MAKE_MODEL_YEAR.nunique())

# drop MAKE and MODEL
df = df.drop(columns=['MAKE', 'MODEL'])

print("total records:", df.shape[0])

# + hidden=true
ABQ_zipcodes = [87101, 87102, 87103, 87104, 
                87105, 87106, 87107, 87108, 87109,
                87110, 87111, 87112, 87113, 87114,
                87115, 87116, 87119, 87120, 87121,
                87122, 87123, 87125, 87131, 87151,
                87153, 87154, 87158, 87176, 87181, 
                87184, 87185, 87187, 87190, 87191,
                87192, 87193, 87194, 87195, 87196, 
                87197, 87198, 87199]

# extract ZIP_CODE
def extract_zipcode(df, col):
    # whether from the city
    df['ABQ'] = 0
    df.loc[df[col].str[:3]=='871', 'ABQ'] = 1
    # which part of the city
    df['DISTRICT'] = df[col].str[3:].str.strip('-')
    return df

# drop ZIP_CODE
df = df.drop(columns=['ZIP_CODE'])

# + hidden=true
df.select_dtypes(include='string').head()

# + hidden=true
# vehicles from a specific brand/model/year fails the test?
tmp = df.groupby('MAKE_MODEL_YEAR')['RESULT'].value_counts().unstack(level=1)
tmp[(tmp[0] > 2) & tmp[1].isnull()]
# -

# ### dates
#

# +
# select date columns
date_cols = ['TEST_SDATE', 'TEST_EDATE']

# transform date columns
def to_datetime(df, date_cols):
    for col in date_cols:
        df[col] = pd.to_datetime(df[col])
    return df

df = to_datetime(df, date_cols)
print(df[date_cols].dtypes)
df[date_cols].head()
# -

# **⚠️ if a vehicle has multiple records from the same date, keep earliest record**
#

# the gap between test ending time and test starting time
df['gap'] = (df['TEST_EDATE'] - df['TEST_SDATE'])/timedelta(hours=1)
# distribution of gap
plt.figure(figsize=(20, 4))
sns.boxplot(x=df.gap);
# drop gap column
df = df.drop(columns=['gap'])

# +
# add a helper column which has date but not hours of the day
df['TEST_DATE'] = df['TEST_SDATE'].dt.date

# keep the first record from the same date for each vehicle
df1 = df.loc[df.groupby(['VIN', 'TEST_DATE'])['TEST_SDATE'].idxmin(),:].copy()

# select the records whose vehicles got tested more than 1 time a day, can capture some extra information here
df2 = df[~df.index.isin(df1.index)]
print("records where vehicles got tested more than 1 time a day:", df2.shape[0])
print("records after removing repeated tests from same day:", df1.shape[0])

# date range
print('data range:', df1.TEST_DATE.min(), df1.TEST_DATE.max())

# drop the some columns
df = df1.drop(columns=['TEST_DATE', 'TEST_EDATE'])
# -

# check the distribution of target
df.RESULT.value_counts()

# +
# engineering VEHICLE_AGE
df['VEHICLE_AGE'] = df.TEST_SDATE.dt.year.astype('int') - df.MODEL_YEAR.astype('int') + 2

# drop TEST_SDATE
# df = df.drop(columns=['TEST_SDATE'])

# check out the VEHICLE_AGE distribution
sns.histplot(data=df, x='VEHICLE_AGE');
# -

print("Age distribution of vehicles which PASSED in emission test")
plt.figure(figsize=(20,4))
sns.countplot(data=df[df.RESULT==1], x='VEHICLE_AGE');

print("Age distribution of vehicles which FAILED in emission test")
plt.figure(figsize=(20,4))
sns.countplot(data=df[df.RESULT==0], x='VEHICLE_AGE');


# +
def extract_time(df, col):
    df['month'] = df[col].dt.month
    df['weekday'] = df[col].dt.weekday
    df['hour'] = df[col].dt.hour
    return df

# df = extract_time(df, 'TEST_SDATE')


# +
# sns.pairplot(data=df[['RESULT', 'month', 'weekday', 'hour']], corner=True)

# + [markdown] heading_collapsed=true
# ### categorical columns 

# + hidden=true
categorical_cols = ['TEST_TYPE', 
                    'TRANS_TYPE', 
                    'VEHICLE_TYPE', 
                    'FUEL_TYPE', 
                    'CYL']

df[categorical_cols].isnull().sum()

# + hidden=true
# distribution of categorical values in Pass and Fail data sets
df['count'] = 1
n = len(categorical_cols)
fig, axs = plt.subplots(n, 1, figsize=(5, 6*n))
for col, ax in zip(categorical_cols, axs.flat):
    tmp = df.groupby(['RESULT', col])['count'].count().groupby('RESULT').apply(lambda x: 100*x/x.sum()).unstack(level=1)
    tmp.plot(kind='bar', stacked=True, ax=ax)
    ax.set_title(col)
df = df.drop(columns=['count'])

# + hidden=true
# TEST_TYPE
df = df.replace({'TEST_TYPE':{'I':0, 'A':1}})
df.TEST_TYPE.value_counts()

# + hidden=true
# TRANS_TYPE
df = df.replace({'TRANS_TYPE':{'A':0, 'M':1}})
df.TRANS_TYPE.value_counts()

# + hidden=true
# VEHICLE_TYPE
# for each vehicle tpye, compare pass and fail
df.groupby(['VEHICLE_TYPE', 'RESULT']).size()\
    .groupby('VEHICLE_TYPE').apply(lambda x: 100*x/x.sum())\
    .unstack(level=1).plot(kind='bar', stacked=True);

# check how many records have more than 1 VEHICLE_TYPE were assigned to the same MAKE_MODEL_YEAR
tmp = df[['VEHICLE_TYPE', 'MAKE_MODEL_YEAR']].groupby('MAKE_MODEL_YEAR').nunique().sort_values('VEHICLE_TYPE')
print('MAKE_MODEL_YEARs with more than 1 VEHICLE_TYPE assigned to them:', tmp[tmp.VEHICLE_TYPE > 1].shape[0])

# + hidden=true
# CYL
# for each CYL, compare pass and fail
df.groupby(['CYL', 'RESULT']).size()\
    .groupby('CYL').apply(lambda x: 100*x/x.sum())\
    .unstack(level=1).plot(kind='bar', stacked=True);

df['CYL_912R'] = 0
df.loc[df.CYL.isin(['9', '12', 'R']), 'CYL_912R'] = 1
print(df.CYL_912R.value_counts())

# + hidden=true
# drop CYL
df = df.drop(columns=['CYL', 'CYL_912R'])

# + hidden=true
# FUEL TYPE 
# for each FUEL_TYPE, compare pass and fail
df.groupby(['FUEL_TYPE', 'RESULT']).size()\
    .groupby('FUEL_TYPE').apply(lambda x: 100*x/x.sum())\
    .unstack(level=1).plot(kind='bar', stacked=True);

# check FUEL_TYPE L
print(df[df.FUEL_TYPE=='L'].groupby('RESULT')['FUEL_TYPE'].value_counts())

# drop FUEL_TYPE
df = df.drop(columns=['FUEL_TYPE'])

# + hidden=true
df.head()
# -

# ### numeric columns
#
#

# +
numeric_cols = ['GVWR', 'ENGINE_SIZE', 'ODOMETER']

fig, axs = plt.subplots(3, 1, figsize=(20, 12))
for col, ax in zip(numeric_cols, axs.flat):
    sns.boxplot(x=df[col], ax=ax)
# -

# statistics
df[numeric_cols].describe()

# +
# let's clean GVWR a little bit

# get median GVWR of each VIN
tmp = df[['VIN', 'GVWR']].groupby('VIN').agg({'GVWR':'median'})
tmp.reset_index(inplace=True)

# merge tmp with df
df.reset_index(inplace=True)
df = df.merge(tmp, how='left', on='VIN', suffixes=('_0',''))
df.set_index('RecordID', inplace=True)

# missing values in GVWR and GVWR_0
print('GVWR_0 = 0:', df[df.GVWR_0==0].shape[0])
print('GVWR = 0:', df[df.GVWR==0].shape[0])

# outliers in GVWR and GVWR_0
print('GVWR_0 > 9000:', df[df.GVWR_0 > 9000].shape[0])
print('GVWR > 9000:', df[df.GVWR > 9000].shape[0])

# distribution of GVWR and GVWR_0
cols = ['GVWR_0', 'GVWR']

fig, axs = plt.subplots(2, 1, figsize=(20, 6))
for col, ax in zip(cols, axs.flat):
    sns.boxplot(x=df[col], ax=ax)

# replace 0 with np.nan
df.loc[df.GVWR==0, 'GVWR'] = np.nan
df.loc[df.GVWR_0==0, 'GVWR_0'] = np.nan

# using GVWR_0 fill missing values in GVWR
df['GVWR'] = df.GVWR.fillna(df.GVWR_0)

# keep GVW and drop GVWR_0
df = df.drop(columns=['GVWR_0'])

# check share of P and F in missing values
df.loc[df.GVWR.isnull(), 'RESULT'].value_counts()
# -

# 0 in ODOMETER
df.loc[df.ODOMETER==0, 'ODOMETER'] = np.nan
print('number of ODOMETER=0:', df.ODOMETER.isnull().sum())
df.loc[df.ODOMETER.isnull(), 'RESULT'].value_counts()

# +
# drop ODOMETER = 0
df = df[df.ODOMETER!=0]

# engineer MILE_YEAR from ODOMETER
df['MILE_YEAR'] = df['ODOMETER']/df['VEHICLE_AGE']

# Americans drive an average of 14,300 miles per year
# outliers in MILE_YEAR
print('records where MILE_YEAR > 100,000 miles:', df[df.MILE_YEAR > 100000].shape[0])
print('vehicles with ODOMETER > 100,000 miles:', df[df.MILE_YEAR > 100000].VIN.nunique())

fig, axs = plt.subplots(2, 1, figsize=(20, 8))
sns.boxplot(x=df['MILE_YEAR'], ax=axs[0])
sns.boxplot(x=df[df.MILE_YEAR < 100000]['MILE_YEAR'], ax=axs[1])
# -

df['MILE_YEAR_CAT'] = df.MILE_YEAR.apply(lambda x: 1 if x > 40000 else 0)
df.groupby(['MILE_YEAR_CAT', 'RESULT']).size()\
    .groupby('MILE_YEAR_CAT').apply(lambda x: 100*x/x.sum())\
    .unstack(level=1).plot(kind='bar', stacked=True);
df = df.drop(columns=['MILE_YEAR_CAT'])

# +
# drop useless columns
tmp = df.drop(columns=['TEST_SDATE', 'VIN', 'MODEL_YEAR'])

# drop na
tmp = tmp.dropna()

# data size
print(tmp.shape)

tmp.head()
# -

# ## vehicles with historical records
#

# how many vehicle were tested how many times
tmp = df.VIN.value_counts().sort_values()
ax = sns.countplot(x=tmp)
ax.set(ylabel='number of vehicles', xlabel='number of tests');

df.VIN.value_counts().value_counts()

# let's check a vehicle that was tested 10 times
df[df.VIN == df.VIN.value_counts().sort_values(ascending=False).index[5]]

# + [markdown] heading_collapsed=true
# ### current records + last test records
#

# + hidden=true
# all the records where vehicles didn't receive test before
df1 = df.loc[df.groupby('VIN')['TEST_SDATE'].idxmin(),:]

# records whose vehicle has historical data
tmp = df.VIN.value_counts()
df2 = df[df.VIN.isin(tmp[tmp>1].index)]

print("Records where vehicles didn't receive test before", df1.shape)
print("Records where vehicles have historical data", df2.shape)

# + hidden=true
# merge historical data once on VIN
df3 = df2.reset_index()
df4 = df2.drop(columns=['MODEL_YEAR', 'MAKE_MODEL_YEAR']).reset_index()
merge = pd.merge(df3, df4, how='left', on='VIN', suffixes=(None, '_1'))


# add a helper column
merge['DELTA'] = merge['TEST_SDATE']-merge['TEST_SDATE_1']

# keep the entries where TEST_DATE_1 is earlier than TEST_DATE
merge = merge[merge.DELTA > timedelta(days=0)]
merge.reset_index(drop=True, inplace=True)


# keep the rows whose DELTA is > 90, meaning vehicle's last check up was at least 90 days earlier than this time

# print(merge.shape)
# merge = merge[merge.DELTA > timedelta(days=90)]
# merge.reset_index(inplace=True)
# print(merge.shape)

# for each RecordID only keep most recent record as previous testing record
df_hist = merge.loc[merge.groupby('RecordID')['DELTA'].idxmin(), :]
df_hist.set_index('RecordID', inplace=True)

print("records in historical dataframe", df_hist.shape[0])
print("unique records in historical dataframe", df_hist.index.nunique())

print("\nStatistic summary of time gaps between current and last tests\n")
df_hist.DELTA.describe()

# + hidden=true
# stack df1 and df_hist
updated_df = pd.concat([df1, df_hist])
print("Shape of the new dataframe", updated_df.shape)

updated_df.sort_index().tail()

# + [markdown] heading_collapsed=true
# ### one entry for one vehicle
#

# + hidden=true
# keep most recent record for each vehicle

df1 = df.loc[df.groupby('VIN')['TEST_SDATE'].idxmax(),:]

# let's take the historical data out

df2 = df[~df.index.isin(df1.index)]

print("number records that don't belong to historical records:", df1.shape[0])
print('number of unique vehicles in df1:', df1.VIN.nunique())    
print('number records that belong to historical records:', df2.shape[0])

# + hidden=true
# let's peel the historical df2 like an onion into multiple dataframes according to timeline

older_records = [] # a list to collect dataframes

updated_df2 = df2

for i in range(10):
    print(f'getting layer {i}')
    latest = updated_df2.loc[updated_df2.groupby('VIN')['TEST_SDATE'].idxmax(), :]
    older_records.append(latest.reset_index())
    tmp = updated_df2[~updated_df2.index.isin(latest.index)]
    updated_df2 = tmp.copy()

# + hidden=true
# merge the dataframes in older_records with df1

updated_df2 = df1
for index, record in enumerate(older_records):
    merged_df = pd.merge(updated_df2, record, how='left', on='VIN', suffixes=(None, '_'+str(index+1)))
    updated_df2 = merged_df.copy()  

print('The size of one-entry-one-vehicle dataframe:', updated_df2.shape)

# let's check out the vehicles that were tested for 11 times in last 7 years
updated_df2[~updated_df2.RecordID_10.isnull()]

# + hidden=true


# + hidden=true


StarcoderdataJupyterScriptsDedupFiltered