<SYSTEM_TASK:>
Given the following text description, write Python code to implement the functionality described below step by step
<END_TASK>
<USER_TASK:>
Description:
Step1:  Chicago taxi fare training experience
This experiment using Scikit-learn Random Forest to train a ML model on Chicago taxi dataset to estimate taxi trip fare in a given time and start-end locations. Selected approach, feature engineering is based on https
Step2:  Query dataset
Step3:  Column info
Watch amount of null values in 'Non-Null Count column'
Step4:  Raw descriptive statistics
Step5:  Feature engineering
Step6:  Remaining null values per column after feature engineering
Step7:  Data profiling
(executing the next cell takes long time)
Step8:  Visual dropoff locations
Step9:  Location histograms
Step10:  Time based explorations
Trip start distribution
Step11:  Trip loginess
Step12:  Fare by trip start hour
Step13:  Split dataframe to examples and output
Step14:  Training pipeline
Step15:  Option 1
Step16:  Option 2
Step17:  Prediction test
Step18:  Cross validation score to test set
<ASSISTANT_TASK:>
Python Code:
import numpy as np
import pandas as pd
from pandas_profiling import ProfileReport
from scipy import stats
from sklearn.ensemble import RandomForestRegressor
from sklearn.compose import ColumnTransformer
from sklearn.model_selection import train_test_split, cross_val_score, GridSearchCV
from sklearn.pipeline import Pipeline
from sklearn.preprocessing import OneHotEncoder, StandardScaler
# MLflow
import mlflow
import mlflow.sklearn
# plotting libraries:
import matplotlib.pyplot as plt
import matplotlib as mpl
import seaborn as sns
# Google clients
import google.auth
from google.cloud import bigquery
from google.cloud import bigquery_storage
# Set default appearance
# - overide maplot libs ugly colours.
# - default figure size
sns.set(color_codes=True)
mpl.rcParams['figure.figsize'] = [13, 8]
%matplotlib inline
BQ_DATASET = 'chicago_taxi_trips'
BQ_TABLE = 'taxi_trips'
BQ_QUERY = 
with tmp_table as (
SELECT trip_seconds, trip_miles, fare, tolls, 
    company, pickup_latitude, pickup_longitude, dropoff_latitude, dropoff_longitude,
    DATETIME(trip_start_timestamp, 'America/Chicago') trip_start_timestamp,
    DATETIME(trip_end_timestamp, 'America/Chicago') trip_end_timestamp,
    CASE WHEN (pickup_community_area IN (56, 64, 76)) OR (dropoff_community_area IN (56, 64, 76)) THEN 1 else 0 END is_airport,
FROM `bigquery-public-data.chicago_taxi_trips.taxi_trips`
WHERE
  dropoff_latitude IS NOT NULL and
  dropoff_longitude IS NOT NULL and
  pickup_latitude IS NOT NULL and
  pickup_longitude IS NOT NULL and
  fare > 0 and 
  trip_miles > 0 and
  MOD(ABS(FARM_FINGERPRINT(unique_key)), 100) {}
ORDER BY RAND()
LIMIT 20000)
SELECT *,
    EXTRACT(YEAR FROM trip_start_timestamp) trip_start_year,
    EXTRACT(MONTH FROM trip_start_timestamp) trip_start_month,
    EXTRACT(DAY FROM trip_start_timestamp) trip_start_day,
    EXTRACT(HOUR FROM trip_start_timestamp) trip_start_hour,
    FORMAT_DATE('%a', DATE(trip_start_timestamp)) trip_start_day_of_week
FROM tmp_table
# Create BigQuery client
credentials, your_project_id = google.auth.default(
    scopes=['https://www.googleapis.com/auth/cloud-platform']
)
bqclient = bigquery.Client(credentials=credentials, project=your_project_id,)
bqstorageclient = bigquery_storage.BigQueryReadClient(credentials=credentials)
Explanation: Chicago taxi fare training experience
This experiment using Scikit-learn Random Forest to train a ML model on Chicago taxi dataset to estimate taxi trip fare in a given time and start-end locations. Selected approach, feature engineering is based on https://github.com/v-loves-avocados/chicago-taxi data exploration and analysis by Aradhana Chaturvedi.
End of explanation
df = (
    bqclient.query(BQ_QUERY.format('between 0 and 99'))
    .result()
    .to_dataframe(bqstorage_client=bqstorageclient)
)
Explanation: Query dataset
End of explanation
display(df.info())
Explanation: Column info
Watch amount of null values in 'Non-Null Count column'
End of explanation
display(df.describe())
Explanation: Raw descriptive statistics
End of explanation
def feature_engineering(data):
    # Add 'N/A' for missing 'Company'
    data.fillna(value={'company':'N/A','tolls':0}, inplace=True)
    # Drop rows contains null data.
    data.dropna(how='any', axis='rows', inplace=True)
    # Pickup and dropoff locations distance
    data['abs_distance'] = (np.hypot(data['dropoff_latitude']-data['pickup_latitude'], data['dropoff_longitude']-data['pickup_longitude']))*100
    # Remove extremes, outliers
    possible_outliers_cols = ['trip_seconds', 'trip_miles', 'fare', 'abs_distance']
    data=data[(np.abs(stats.zscore(data[possible_outliers_cols])) < 3).all(axis=1)].copy()
    # Reduce location accuracy
    data=data.round({'pickup_latitude': 3, 'pickup_longitude': 3, 'dropoff_latitude':3, 'dropoff_longitude':3})
    return data
df=feature_engineering(df)
display(df.describe())
Explanation: Feature engineering
End of explanation
print(df.isnull().sum().sort_values(ascending=False))
Explanation: Remaining null values per column after feature engineering
End of explanation
ProfileReport(df, title='Chicago taxi dataset profiling Report').to_notebook_iframe()
Explanation: Data profiling
(executing the next cell takes long time)
End of explanation
sc = plt.scatter(df.dropoff_longitude, df.dropoff_latitude, c = df['fare'], cmap = 'summer')
plt.colorbar(sc)
Explanation: Visual dropoff locations
End of explanation
fig, axs = plt.subplots(2)
fig.suptitle('Pickup location histograms')
df.hist('pickup_longitude', bins=100, ax=axs[0])
df.hist('pickup_latitude', bins=100, ax=axs[1])
plt.show()
fig, axs = plt.subplots(2)
fig.suptitle('Dropoff location histograms')
df.hist('dropoff_longitude', bins=100, ax=axs[0])
df.hist('dropoff_latitude', bins=100, ax=axs[1])
plt.show()
Explanation: Location histograms
End of explanation
fig, axs = plt.subplots(4)
fig.suptitle('Trip start histograms')
fig.set_size_inches(18, 12, forward=True)
df.hist('trip_start_year', bins=8, ax=axs[0], )
df.hist('trip_start_month', bins=12, ax=axs[1])
df.hist('trip_start_day', bins=31, ax=axs[2])
df.hist('trip_start_hour', bins=24, ax=axs[3])
plt.show()
Explanation: Time based explorations
Trip start distribution
End of explanation
fig, axs = plt.subplots(2)
fig.set_size_inches(18, 8, forward=True)
df.hist('trip_miles', bins=50, ax=axs[0])
df.hist('trip_seconds', bins=50, ax=axs[1])
plt.show()
Explanation: Trip loginess
End of explanation
display(df.groupby('trip_start_hour')['fare'].mean().plot())
Explanation: Fare by trip start hour
End of explanation
# Drop complex fields and split dataframe to examples and output
mlflow.log_param('training_shape', f'{df.shape}')
X=df.drop(['trip_start_timestamp'],axis=1)
y=df['fare']
Explanation: Split dataframe to examples and output
End of explanation
# global variables
experiment_name = 'chicago-taxi-1'
ct_pipe = ColumnTransformer(transformers=[
    ('hourly_cat', OneHotEncoder(categories=[range(0,24)], sparse = False), ['trip_start_hour']),
    ('dow', OneHotEncoder(categories=[['Mon', 'Tue', 'Sun', 'Wed', 'Sat', 'Fri', 'Thu']], sparse = False), ['trip_start_day_of_week']),
    ('std_scaler', StandardScaler(), [
        'trip_start_year',
        'abs_distance',
        'pickup_longitude',
        'pickup_latitude',
        'dropoff_longitude',
        'dropoff_latitude',
        'trip_miles',
        'trip_seconds'])
])
X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.25, random_state=123)
X_train=X_train.drop('fare', axis=1)
# for more details: https://scikit-learn.org/stable/modules/generated/sklearn.ensemble.RandomForestRegressor.html
rfr_pipe = Pipeline([
    ('ct', ct_pipe),
    ('forest_reg', RandomForestRegressor(
        n_estimators = 20,
        max_features = 'auto',
        n_jobs = -1,
        random_state = 3,
        max_depth=None,
        max_leaf_nodes=None,
    ))
])
rfr_score = cross_val_score(rfr_pipe, X_train, y_train, scoring = 'neg_mean_squared_error', cv = 5)
rfr_rmse = np.sqrt(-rfr_score)
rfr_rmse.mean()
mlflow.log_metric('train_cross_valid_score_rmse_mean', np.sqrt(-rfr_score).mean())
mlflow.log_param('number_of_estimators', 20)
Explanation: Training pipeline
End of explanation
# To see all RandomForestRegressor hyper parameters:
# estimator=RandomForestRegressor()
# display(estimator.get_params())
# Train model
mlflow.set_experiment('chicago-taxi-0')
# mlflow.sklearn.autolog()
with mlflow.start_run(nested=True) as mlflow_run:
    final_model=rfr_pipe.fit(X_train, y_train)
    mlflow.sklearn.log_model(final_model, 'chicago_rnd_forest')
Explanation: Option 1: Simple training
(~fast)
End of explanation
param_grid = {'forest_reg__n_estimators': [5, 250], 'forest_reg__max_features': [6, 16, 'auto']}
forest_gs = GridSearchCV(rfr_pipe, param_grid, cv = 5, scoring = 'neg_mean_squared_error', n_jobs = -1)
forest_gs.fit(X_train, y_train)
print(f'Best parameters: {forest_gs.best_params_}')
print(f'Best score: {np.sqrt(-forest_gs.best_score_)}')
print(f'(All scores: {np.sqrt(-forest_gs.cv_results_['mean_test_score'])})')
final_model=forest_gs.best_estimator_
Explanation: Option 2: Parameter search + training
(time consuming)
End of explanation
X_pred = pd.DataFrame(X_test, columns=X_test.columns)
X_pred['fare_pred'] = final_model.predict(X_test.drop('fare',axis=1))
X_pred.head(5)
Explanation: Prediction test
End of explanation
rfr_score = cross_val_score(final_model, X_test, y_test, scoring='neg_mean_squared_error', cv = 5)
rfr_rmse = np.sqrt(-rfr_score)
rfr_rmse.mean()
mlflow.log_metric('eval_cross_valid_score_rmse_mean', np.sqrt(-rfr_score).mean())
# Comparer test
def model_comparer(job_name, **kwargs):
    print(f'Model blessing: "{job_name}"')
    experiment = mlflow.get_experiment_by_name(experiment_name)
    filter_string = f"tags.job_name ILIKE '{job_name}_%'"
    df = mlflow.search_runs([experiment.experiment_id], filter_string=filter_string)
    display(df)
    # Compare
    # Available columns:
    # run_id	experiment_id	status	artifact_uri	start_time	end_time	metrics.train_cross_valid_score_rmse_mean	params.number_of_estimators	tags.job_name	tags.mlflow.source.name	tags.mlflow.user	tags.mlflow.source.type	tags.version
    eval_max = df.loc[df['metrics.eval_cross_valid_score_rmse_mean'].idxmax()]
    train_max= df.loc[df['metrics.train_cross_valid_score_rmse_mean'].idxmax()]
    display(eval_max)
    return eval_max
# You need to set a previous training job name manually. Which is following this naming pattern: training_job_...time stamp...
best_run = model_comparer('training_job_20210119T220534')
client = mlflow.tracking.MlflowClient()
def register_model(run_id, model_name):
    model_uri = f'runs:/{run_id}/{model_name}'
    registered_model = mlflow.register_model(model_uri, model_name)
    print(registered_model)
registered_models=client.search_registered_models(filter_string=f"name='{experiment_name}'", max_results=1, order_by=['timestamp DESC'])
if len(registered_models) ==0:
    register_model(best_run.run_id, experiment_name)
else:
    last_version = registered_models[0].latest_versions[0]
    run = client.get_run(last_version.run_id)
    print(run)
    if not run:
        print(f'Registered version run missing!')            
        
    last_eval_metric=run.data.metrics['eval_cross_valid_score_rmse_mean']
    best_run_metric=best_run['metrics.eval_cross_valid_score_rmse_mean']
    if last_eval_metric<best_run_metric:
        print(f'Register better version with metric: {best_run_metric}')
        register_model(best_run.run_id, experiment_name)
    else:
        print(f'Registered version still better. Metric: {last_eval_metric}')
Explanation: Cross validation score to test set
End of explanation
<END_TASK>