10  Model Building for Supervised Machine Learning

Click here to download this chapter as a Jupyter (.ipynb) file.

This chapter explains and demonstrates how the following tasks can be structured into processes to build a supervised machine learning model, using scikit-learn tools:

• Preprocessing - Getting the data ready for the specific machine-learning algorithms that will be tried
• Estimating generalization performance - Estimating how a machine-learning algorithm will perform on new data that was not used to train it. This is also called performance on out-of-sample data
• Model tuning - Determining the best parameter setting(s) for a machine-learning algorithm
• Algorithm selection - Determining which predictive algorithm to use for our final deployed model

In particular, this chapter describes and shows examples of two methods of model-building:

• The One-Outer-Split Model-Building Process
• Nested Cross Validation

10.1 Module and Function Imports

import time
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt

from sklearn.base import clone

from sklearn.model_selection import train_test_split
from sklearn.model_selection import cross_validate
from sklearn.model_selection import StratifiedKFold, RepeatedStratifiedKFold
from sklearn.model_selection import GridSearchCV

from sklearn.preprocessing import MinMaxScaler
from sklearn.preprocessing import StandardScaler
from sklearn.preprocessing import OneHotEncoder
from sklearn.preprocessing import PolynomialFeatures

from sklearn.compose import ColumnTransformer

from sklearn.pipeline import Pipeline

from sklearn.neighbors import KNeighborsClassifier
from sklearn.linear_model import LogisticRegression

10.2 Load the Urban Delivery Delays Data

This synthetic dataset contains information about 5,000 deliveries made by a fictional delivery company operating in a metropolitan area. The goal is to predict whether a delivery will arrive late based on its location, expected rainfall, and service type.

The dataset contains four predictors:

• East_km: East-west coordinate of the delivery location.
• North_km: North-south coordinate of the delivery location.
• Rain_mm: Expected rainfall in millimeters.
• Service_Type: Standard or Express delivery.

The target variable, Late_Delivery, equals 1 if the delivery misses its promised delivery window and 0 otherwise. Approximately 40% of deliveries are late. The data include several localized areas with high congestion, creating nonlinear relationships between location and delivery delays.

Business question: Can we build a model to identify deliveries at risk of arriving late, allowing the company to take preventive action?

deliveries = pd.read_csv("https://neuronjolt.com/data/delivery_delays.csv")

deliveries.info()

<class 'pandas.DataFrame'>
RangeIndex: 5000 entries, 0 to 4999
Data columns (total 5 columns):
 #   Column         Non-Null Count  Dtype  
---  ------         --------------  -----  
 0   East_km        5000 non-null   float64
 1   North_km       5000 non-null   float64
 2   Rain_mm        5000 non-null   float64
 3   Service_Type   5000 non-null   str    
 4   Late_Delivery  5000 non-null   int64  
dtypes: float64(3), int64(1), str(1)
memory usage: 195.4 KB

10.2.1 Assign the features to X and the target to y

X = deliveries.drop(columns="Late_Delivery")
y = deliveries["Late_Delivery"]

10.2.2 Define lists for the categorical and numeric features

We can then use these lists in our column transformer, rather than having to list out the specific features for each step.

num_features = ["East_km", "North_km", "Rain_mm"]
cat_features = ["Service_Type"]

10.3 Cross Validation

Cross-validation is a method for evaluating a model or model-building procedure by repeatedly dividing the available data into training and validation subsets. We will use cross-validation to compare hyperparameter settings and select algorithms. We begin by explaining how K-fold cross-validation works.”

10.3.1 K-Fold Cross Validation

When fitting a machine-learning model, we need a way to determine which hyperparameter settings produce good predictions on observations that were not used to fit the model.

K-fold cross-validation provides a systematic way to do this. We divide the observations into K approximately equal-sized groups called folds. Often, we shuffle the observations before creating the folds. We then repeat the following process K times:

1. Fit the model using K − 1 folds, called the training folds.
2. Evaluate the model on the remaining fold, called the validation fold.

Each fold serves as the validation fold exactly once. We average the K validation fold scores to summarize the performance of that algorithm and hyperparameter setting.

We can then compare these average validation scores across different hyperparameter settings. This is the basic idea behind grid-search cross-validation, which we will discuss later in this chapter.

Although cross-validation scores provide useful information about predictive performance, we’ll see later why we also reserve an independent test set when building and evaluating a final model.

five_fold_cv_diagram_2.png

One advantage of K-fold cross-validation over repeated random splitting is that every observation is used for validation exactly once and for training K − 1 times. With repeated random splitting, some observations may appear in validation sets several times, while others may never be included in a validation set.

10.3.2 Define a Splitter to use other types of cross validation

Splitter objects define the splitting strategy for cross validation. They are created with constructor functions similar to how estimators or transformers are created. Examples of the splitters available in scikit-learn include the following:

• K-fold cross validation (create with KFold() function) - divide the data into k folds; each fold takes a turn as the validation fold, with the other folds used for model-training. The KFold() splitter function has a shuffle parameter that can be set to True to shuffle the data before the folds are defined, which can help prevent problems caused by data that is sorted in a systematic order.
• Stratified K-fold cross validation (create with StratifiedKFold() function) - like k-fold cross validation, but the splits maintain the relative proportions of the target class values in each fold.
• Repeated K-fold cross validation (create with RepeatedKFold() function) - perform k-fold more than once, shuffling the data in between
• Repeated Stratified K-fold cross validation (create with RepeatedStratifiedKFold() function) - like repeated K-fold, but stratified
• Shuffle-split cross validation (create with ShuffleSplit() function) - like random train-test splitting, specifying either number of instances in train and test, or percentage of dataset in train and test.
• Stratified shuffle-split cross validation (create with StratifiedShuffleSplit() function) - shuffle split cross validation with stratification added.

Note that instances of these Splitter classes may be passed to the cross_validate(), or GridSearchCV() functions via the cv parameter.

10.4 Model Tuning with Grid Search

Model Tuning is the general term for the process of finding the parameter settings for a prediction algorithm that result in the best performance on a particular dataset. Other names for this process are Parameter Tuning and Hyperparameter Tuning. Grid Search is a common technique used for model tuning. Grid search refers to the process of estimating out-of-sample performance for an algorithm for each of a set of parameter settings in order to find the parameter setting at which the model performs best on the data. Grid search uses cross validation internally. At each split of the cross validation, models using each of the parameter settings are built on the training data and scored on the validation data. In this way, each parameter combination will have several validation fold scores which may then be averaged to obtain an estimate of the model’s performance with that particular parameter setting.

We actually performed a grid search in Chapter 7, when we generated the data for a validation curve graph for the K-nn model. We used nested looping code, with data splitting in the outer loop and training and validation for all values of \(k\) in the inner loop. In this chapter we will use the scikit-learn GridSearchCV() function to perform grid search rather than constructing our own nested looping code.

The diagram below illustrates a grid search that uses 5-fold cross validation and tries 4 different settings for k in a K-nn model. At each split, models for each k parameter setting are trained on the training folds and evaluated on the validation fold. When this process completes there are 5 validation scores for each value of the k parameter. In its final step the grid search fits a model to all 5 folds, using the k parameter value that had the highest mean validation fold score.

grid_search_diagram_3.png

10.5 The One-Outer Split Model-Building Process

The purpose of a grid search is to find the best parameter setting(s) for an algorithm. After the cross validation within the grid search completes, the mean validation fold performance across the splits for each parameter setting is calculated. The setting with the highest mean validation fold performance is the best parameter setting. The grid search will then train the algorithm again, at its best parameter setting, on all the model-training data.

It is important to note, however, that the highest mean validation score from a grid search may be an optimistic estimate of generalization performance because the same validation scores were used to select the winning hyperparameter setting. Although cross-validation scores provide estimates of predictive performance, some hyperparameter combinations may receive unusually high scores simply by chance. Because grid search selects the combination with the highest observed score, the winning score tends to be somewhat optimistic as an estimate of performance on new data. Therefore, we use these scores for model selection and reserve the independent test set for our final estimate of generalization performance.

In the one-outer-split model-building strategy the data is first split into training data and test data. The training data is used for model tuning and selection and the test data is used to estimate generalization performance. As part of the model tuning process the training data from the first level of splitting gets repeatedly split into model-training data and validation data by the grid searches that test each algorithm at each parameter setting.

So, when comparing multiple algorithms, we perform a separate grid search for each one using the training data. We select the algorithm whose best hyperparameter setting produces the highest mean cross-validation score. Only after selecting the winning algorithm do we evaluate it on the independent test set.

one_outer_split_model_building_process.png

10.6 One-Outer-Split Model Building with the Deliveries Data

10.6.1 Choose between K-nn and logistic regression for the final model

In this section we will practice the one-outer-split model-building process. We will use it to build a model to make predictions for the deliveries data, with logistic regression and K-nn being the candidate prediction algorithms. We will conduct a grid search on the training data with a logistic regression algorithm, and then view the mean validation fold score at the best parameter setting. Then we will conduct a similar grid search for K-nn and view the mean validation fold score at its best parameter setting. The algorithm with the best mean performance on the validation folds at its best setting is the one we will choose for our final model. We score the winning algorithm’s fitted grid search object on the test data to estimate its generalization performance on this data.

The final step is to then re-run the grid search for the best-performing algorithm on ALL the data to tune the model and prepare it to make predictions on future data.

10.6.2 Split the deliveries data into training and test sets

Split the deliveries data into training and test sets. Use a stratified split, since the target variable is categorical. We set the random_state parameter so that the random split is repeatable.

X_train, X_test, y_train, y_test = train_test_split(X, 
                                                    y, 
                                                    stratify = y, 
                                                    random_state = 102)

10.6.3 Logistic Regression

10.6.3.1 Combine the pre-processing and estimation steps in a pipeline

Logistic regression is a linear classification algorithm that models the probability of a categorical outcome. Like linear regression, it uses a linear combination of the predictors, but it transforms that combination into a predicted probability. It has a parameter, C, that controls the model’s flexibility to fit the training data. The higher C is, the more freedom the algorithm has to fit training data. It is recommended to standard scale numeric features and one-hot encode categorical features before running logistic regression.

In the code cells below we create a column transformer for the deliveries data and then package the column transformer together with the logistic regression estimator in a pipeline. In the column transformer we will one-hot encode the categorical features. We will try doing two things to the numeric features. First we will apply the PolynomialFeatures() function out to degree 3 to add squared features, cubed features, and interactions, and then we will standard-scale the resulting numeric features. Because we are doing two things to the numeric features - adding polynomial features and then standard scaling - we will package those two things into a pipeline.

Adding polynomial and interaction features gives logistic regression the ability to model more complex relationships. This may reduce underfitting, but also increases the potential for overfitting. The C parameter gives us a tool to use to make fine-grained reductions in the flexibility level of the algorithm, and the grid search will help us to set the C parameter to the right level. This is a common pre-processing approach that we will learn more about in the next chapter.

ct_lreg = ColumnTransformer([
    ('poly_std', Pipeline([('poly', PolynomialFeatures(degree = 3, include_bias = False)),
                            ('std', StandardScaler())]), num_features),
    ('one_hot', OneHotEncoder(handle_unknown = 'ignore'), cat_features)
])

pipe_lreg = Pipeline([
    ('pre', ct_lreg),
    ('lreg', LogisticRegression(max_iter = 10_000))
])

10.6.3.2 Specify the parameter values to try

The pipeline we created can serve as the estimator in the GridSearchCV() function because its last step is an estimator. When we specify the values of parameter C that we want to try we need to specify them in a way that provides the information needed to pass the parameters to the correct step of the pipeline. To do so, we need to specify the parameters in a dictionary and use the following convention for the key value in the dictionary:

the step in the pipeline, followed by two underscores, followed by the parameter name.

The parameter C is the parameter for logistic regression that influences the model’s flexibility to fit the training data. The higher C is the more flexible the model is to fit the training data. It should be varied on a logarithmic scale.

Below we set up a parameter grid to test various settings of the C parameter, generating 16 values on a logarithmic scale between \(10^{-3}\) and \(10^3\), or between 0.001 and 1,000.

How do we know what range to use for the parameter setting? Although experience can eventually provide you with good starting points for various parameter settings, the short answer is that we don’t know before running the grid search what the best parameter setting will be, so we typically set the parameter to a fairly wide range and try the grid search. Then we examine the results of the grid search to see if we need to adjust the parameter setting values. The best way to do so is to make a validation curve plot. A validation curve shows how model performance changes as a hyperparameter varies. It can help us identify promising parameter ranges and determine whether additional values should be investigated. Sometimes validation performance reaches a peak and then declines, while other times it levels off as model flexibility increases.

pgrid_lreg = {
    'lreg__C': np.logspace(-3, 3, 16)
}

10.6.3.3 Set up Cross Validation (Splitting) Strategy

Let’s use 2 x stratified 8-fold cross validation, which results in 16 splits. Since this is a classification task we will use stratified splits. We will also set the random_state parameter so that the shuffling and splitting are repeatable. Then, when we set up the grid search for K-nn we can use the same splitter to ensure that we are comparing the two algorithms over the same splits.

rskf8x2 = RepeatedStratifiedKFold(n_splits = 8, 
                                  n_repeats = 2, 
                                  random_state = 102)

10.6.3.4 Instantiate an object of the GridSearchCV class

Specify the estimator, the parameter grid to search, and the cross validation strategy. The cross validation will be performed for every parameter combination from the parameter grid.

Note that we also retain the scores on the training data so that we can make a validation curve plot, and we set the scoring metric to ‘accuracy’. This is the default performance metric, but in a later chapter we will learn about other metrics that can be used instead.

gs_lreg = GridSearchCV(estimator = pipe_lreg,
                           param_grid = pgrid_lreg,
                           cv = rskf8x2,
                           scoring = 'accuracy', 
                           return_train_score = True)

10.6.3.5 Do an initial fit of the grid search to the training data

The fit() method of GridSearchCV finds the best parameters through a cross validation process in which at each split a model is fitted and evaluated for every value of the parameter setting. The best parameter setting is the one that has the highest mean performance on the validation folds. After the best parameter setting is identified the grid search automatically fits a model with the best parameter settings to the whole training set. The fitted GridSearchCV object is thus itself an estimator, with fit(), predict(), and score() methods.

gs_lreg.fit(X_train, y_train)

GridSearchCV(cv=RepeatedStratifiedKFold(n_repeats=2, n_splits=8, random_state=102),
             estimator=Pipeline(steps=[('pre',
                                        ColumnTransformer(transformers=[('poly_std',
                                                                         Pipeline(steps=[('poly',
                                                                                          PolynomialFeatures(degree=3,
                                                                                                             include_bias=False)),
                                                                                         ('std',
                                                                                          StandardScaler())]),
                                                                         ['East_km',
                                                                          'North_km',
                                                                          'Rain_mm']),
                                                                        ('one_hot',
                                                                         OneHotEncoder(handle_unknown='ignore'),
                                                                         ['Servic...
                                        LogisticRegression(max_iter=10000))]),
             param_grid={'lreg__C': array([1.00000000e-03, 2.51188643e-03, 6.30957344e-03, 1.58489319e-02,
       3.98107171e-02, 1.00000000e-01, 2.51188643e-01, 6.30957344e-01,
       1.58489319e+00, 3.98107171e+00, 1.00000000e+01, 2.51188643e+01,
       6.30957344e+01, 1.58489319e+02, 3.98107171e+02, 1.00000000e+03])},
             return_train_score=True, scoring='accuracy')

In a Jupyter environment, please rerun this cell to show the HTML representation or trust the notebook.
On GitHub, the HTML representation is unable to render, please try loading this page with nbviewer.org.

10.6.3.6 View the pattern of validation performance vs parameter setting

It can be useful to check the best parameter setting after the grid search and compare it to the range of values you specified for the parameter. If the best parameter setting is at or very near the ends of the ranges of values you specified in your parameter grid for the grid search you may want to widen the parameter grid and re-run the grid search, depending on the prediction algorithm.

Conversely, if the best parameter setting is in the middle of the range of values you specified for your parameter grid you might want to narrow the range and re-run the grid search.

Below we will view the best parameter setting and also make a validation curve plot to see the pattern of validation fold performance vs the parameter settings.

# Best parameter setting from the grid search
gs_lreg.best_params_

{'lreg__C': np.float64(0.25118864315095824)}

To make the validation curve plot we first need to put the relevant grid search results into a DataFrame. The grid search results are stored in a dictionary accessed with the grid search’s cv_results_ property. This results dictionary can be easily converted to a pandas DataFrame for viewing and analysis.

results_lreg = pd.DataFrame(gs_lreg.cv_results_)

# Create the validation curve plot
plt.plot(results_lreg['param_lreg__C'], results_lreg['mean_train_score'])
plt.plot(results_lreg['param_lreg__C'], results_lreg['mean_test_score'])
plt.axhline(y = max(results_lreg['mean_test_score']),
           color = 'gray',
           linestyle = ':')
plt.axvline(x = results_lreg['param_lreg__C'][results_lreg['mean_test_score'].idxmax()],
           color = 'gray',
           linestyle = ':')
plt.title("Validation Curve Plot for Logistic Regression")
plt.xscale('log')
plt.xlabel("Regularization Parameter (C)")
plt.ylabel("Accuracy")
plt.legend(["Mod-Train Accuracy", "Validation Accuracy", "Max Validation Accuracy"], 
           loc="best");

10.6.3.7 Redefine the parameter grid (if desired) and then re-run the grid search

Since the best parameter setting is between \(10^{-2}\) and \(10^1\) let’s redefine the parameter grid to that range and include more parameter candidates so that they are spaced closer together, then re-run the grid search. Narrowing the parameter range and trying more closely spaced values may help us identify a parameter setting with better validation performance. However, small differences in validation scores may simply reflect random variation in the data.

For linear models such as logistic regression (and ridge and lasso, which we will see in the next chapter), the validation curve sometimes rises and then levels off rather than declining as model flexibility increases. In these cases, additional flexibility provides little or no improvement in validation performance. Although grid search may select a highly flexible setting with a slightly higher validation score, we may prefer a less flexible setting because it produces a simpler model. This reflects a general principle in machine learning: when two models perform similarly, we often prefer the simpler one.

# re-define the parameter settings, based on what we saw in the 
# validation curve plot
pgrid_lreg_adj = {
    'lreg__C': np.logspace(-2, 1, 15)
}

# re-define the grid search
gs_lreg_adj = GridSearchCV(estimator = pipe_lreg,
                           param_grid = pgrid_lreg_adj,
                           cv = rskf8x2,
                           scoring = 'accuracy', 
                           return_train_score = True)

gs_lreg_adj.fit(X_train, y_train)

GridSearchCV(cv=RepeatedStratifiedKFold(n_repeats=2, n_splits=8, random_state=102),
             estimator=Pipeline(steps=[('pre',
                                        ColumnTransformer(transformers=[('poly_std',
                                                                         Pipeline(steps=[('poly',
                                                                                          PolynomialFeatures(degree=3,
                                                                                                             include_bias=False)),
                                                                                         ('std',
                                                                                          StandardScaler())]),
                                                                         ['East_km',
                                                                          'North_km',
                                                                          'Rain_mm']),
                                                                        ('one_hot',
                                                                         OneHotEncoder(handle_unknown='ignore'),
                                                                         ['Service_Type'])])),
                                       ('lreg',
                                        LogisticRegression(max_iter=10000))]),
             param_grid={'lreg__C': array([ 0.01      ,  0.01637894,  0.02682696,  0.04393971,  0.07196857,
        0.11787686,  0.19306977,  0.31622777,  0.51794747,  0.8483429 ,
        1.38949549,  2.27584593,  3.72759372,  6.1054023 , 10.        ])},
             return_train_score=True, scoring='accuracy')

In a Jupyter environment, please rerun this cell to show the HTML representation or trust the notebook.
On GitHub, the HTML representation is unable to render, please try loading this page with nbviewer.org.

# View the best parameter setting from the grid search
gs_lreg_adj.best_params_

{'lreg__C': np.float64(0.517947467923121)}

10.6.3.8 Obtain the mean CV validation score at the best parameter setting

The grid search did 16 splits of the training data into model-training and validation data. At each split it built and scored a logistic regression model at each of the C parameter settings. So, we now have 16 validation scores for each parameter setting. The grid search then identified the parameter setting with the highest mean validation score and fit a logistic regression model at that setting to all the training data (model-training data + validation data). You can extract the best parameter setting with the grid search’s best_params_ property, as we have seen. You can also extract the mean validation score at that setting with the best_score_ property. This is the score that we can use to compare logistic regression to K-nn on this data.

gs_lreg_adj.best_score_

np.float64(0.6018661841889454)

10.6.4 K-nn

Now we will repeat the process with K-nn to find its best_score_. The algorithm with the highest best_score_ is selected for future predictions.

# Set up column transformer for K-nn
ct_knn = ColumnTransformer([
    ('mms', MinMaxScaler(), num_features),
    ('ohe', OneHotEncoder(handle_unknown = 'ignore'), cat_features)
])

# Set up the pipeline for K-nn
pipe_knn = Pipeline([
    ('prep', ct_knn),
    ('knn', KNeighborsClassifier())
])

# Set up the parameter grid for K-nn
pgrid_knn = {'knn__n_neighbors': range(1, 100, 4)}

# Set up the grid search for K-nn. We will use the same splitter we
# defined for logistic regression so that we compare them over the same splits.
gs_knn = GridSearchCV(estimator = pipe_knn,
                      param_grid = pgrid_knn,
                      cv = rskf8x2,
                      n_jobs = -1,
                      return_train_score = True)

# Fit the grid search to the training data
gs_knn.fit(X_train, y_train)

GridSearchCV(cv=RepeatedStratifiedKFold(n_repeats=2, n_splits=8, random_state=102),
             estimator=Pipeline(steps=[('prep',
                                        ColumnTransformer(transformers=[('mms',
                                                                         MinMaxScaler(),
                                                                         ['East_km',
                                                                          'North_km',
                                                                          'Rain_mm']),
                                                                        ('ohe',
                                                                         OneHotEncoder(handle_unknown='ignore'),
                                                                         ['Service_Type'])])),
                                       ('knn', KNeighborsClassifier())]),
             n_jobs=-1, param_grid={'knn__n_neighbors': range(1, 100, 4)},
             return_train_score=True)

In a Jupyter environment, please rerun this cell to show the HTML representation or trust the notebook.
On GitHub, the HTML representation is unable to render, please try loading this page with nbviewer.org.

results_knn = pd.DataFrame(gs_knn.cv_results_)

# Create the validation curve plot
plt.plot(results_knn['param_knn__n_neighbors'], results_knn['mean_train_score'])
plt.plot(results_knn['param_knn__n_neighbors'], results_knn['mean_test_score'])
plt.axhline(y = max(results_knn['mean_test_score']),
           color = 'gray',
           linestyle = ':')
plt.axvline(x = results_knn['param_knn__n_neighbors'][results_knn['mean_test_score'].idxmax()],
           color = 'gray',
           linestyle = ':')
plt.gca().invert_xaxis()
plt.title("Validation Curve Plot for K-nn")
plt.xlabel("Number of Neighbors")
plt.ylabel("Accuracy")
plt.legend(["Mod-Train Accuracy", "Validation Accuracy", "Max Validation Accuracy"], 
           loc="best");

# Adjust parameter grid, re-define and re-fit grid search
pgrid_knn_adj = {'knn__n_neighbors': range(10, 61)}

gs_knn_adj = GridSearchCV(estimator = pipe_knn,
                          param_grid = pgrid_knn_adj,
                          cv = rskf8x2,
                          n_jobs = -1,
                          return_train_score = True)

gs_knn_adj.fit(X_train, y_train)

GridSearchCV(cv=RepeatedStratifiedKFold(n_repeats=2, n_splits=8, random_state=102),
             estimator=Pipeline(steps=[('prep',
                                        ColumnTransformer(transformers=[('mms',
                                                                         MinMaxScaler(),
                                                                         ['East_km',
                                                                          'North_km',
                                                                          'Rain_mm']),
                                                                        ('ohe',
                                                                         OneHotEncoder(handle_unknown='ignore'),
                                                                         ['Service_Type'])])),
                                       ('knn', KNeighborsClassifier())]),
             n_jobs=-1, param_grid={'knn__n_neighbors': range(10, 61)},
             return_train_score=True)

In a Jupyter environment, please rerun this cell to show the HTML representation or trust the notebook.
On GitHub, the HTML representation is unable to render, please try loading this page with nbviewer.org.

# Get mean validation fold score for K-nn
gs_knn_adj.best_score_

np.float64(0.704927457492756)

10.6.5 Identify the winning algorithm

print("Logistic regression:", gs_lreg_adj.best_score_)
print("KNN:", gs_knn_adj.best_score_)

Logistic regression: 0.6018661841889454
KNN: 0.704927457492756

10.6.6 Score the grid search of the winning algorithm on the test data to estimate generalization performance

The training data from the original split was used in the grid search to find the best algorithm and parameter setting. K-nn performed substantially better than logistic regression. Now we can use the test data to estimate that K-nn model’s generalization performance, since the test data wasn’t used to set the parameters or choose among the algorithms.

gs_knn_adj.score(X_test, y_test)

0.6792

We have now selected an algorithm for future predictions. Our estimate of its accuracy on new deliveries is 67.92%. Next we will prepare the model to make future predictions.

10.6.7 Prepare best algorithm for future predictions

To prepare the K-nn algorithm for future predictions we perform a final tuning of the model by re-running the grid search on all the data. First, however, we clone it so that we can preserve the original grid search object and all of its results. After the cloned grid search object is fit to all the data it may be used as an estimator to make predictions for future data, but here we extract best_estimator_, which is the complete fitted pipeline containing the preprocessing steps and the prediction algorithm with its selected hyperparameters.

final_search = clone(gs_knn_adj)

final_search.fit(X, y)

final_model = final_search.best_estimator_

10.6.8 Summary of the One-Outer-Split Model-Building Process

We used cross-validation on the training data to tune logistic regression and K-nn, then compared their best mean validation scores to select K-nn. Next, we evaluated the selected model on the independent test set to estimate its performance on new observations. Finally, we repeated grid-search cross-validation for K-nn using all available data to produce our deployment model.

The key principle is that validation scores guide model selection, while the independent test set provides our final evaluation of the selected model.

10.7 Model Building with Multiple Outer Splits: Nested Cross Validation

In the model-building process described and implemented above we used only one outer split into training and test sets, so our estimate of generalization performance is based on just one hold-out dataset. The one-outer-split approach is especially attractive when the dataset is large enough to provide substantial training and test sets. With smaller datasets, performance estimates from a single test set may be unstable, making repeated evaluation methods such as nested cross validation useful. Nested cross validation is another model-building approach that makes its estimates of generalization performance based on multiple train-test splits. To do so, the single outer split of the one-outer-split model-building process may be replaced by a cross validation process. The resulting process is called nested cross validation because both the outer split and the inner model tuning process utilize cross validation. Nested cross validation is implemented by using the cross_validate() function with a grid search as the estimator.

Nested cross validation works in the following way:

• Outer cross validation repeatedly splits the data into training and test sets.
• For each split of the outer cross validation, a grid search is fit to the training set. The cross validation in the grid search is the inner cross validation. It repeatedly splits the training data from the outer cross validation into model-training and validation sets to find the best parameter setting and then fits the model to all the training data from the outer cross validation, using the best parameter setting.
• The fitted grid search at the best parameter setting is then scored on the test set for that iteration of the outer loop.

The result of the nested cross validation procedure is a list of test set scores that give us an idea how the model-building procedure generalizes. Since each outer loop split can result in different best parameters on the inner cross validation this procedure is NOT used to find best parameters, or to find a model to use on new data. Rather, it is used to investigate which predictive algorithm performs best on a particular dataset. After the best-performing algorithm is determined you can then perform a grid search on that algorithm, using all the data, to find the best parameter settings for that algorithm and fit a model at the best parameter setting to all the data. This fitted model may then be used to make predictions on future data.

Nested cross validation is implemented in scikit-learn by running cross_validate() with an instance of GridSearchCV() as the estimator.

10.7.1 Compare generalization performance of Logistic Regression and K-nn

Below we perform a separate nested cross validation for two different predictive algorithms, K-nn and logistic regression, to compare their generalization performance on the deliveries data. Then, we will select the algorithm with the highest mean outer-CV score on all the data and run its grid search to create our final model to use to make predictions on future data.

Because K-nn and logistic regression require different preprocessing and have different parameters to tune, we set up separate column transformers, pipelines, parameter grids, and grid searches for each algorithm. We use splitter objects to set up the cross validation strategies. Specifically, we will use repeated stratified 8-fold cross-validation, with two repetitions, for the inner grid searches. For the outer cross-validation, we will use stratified 8-fold cross-validation without repetition. Finally, we also specify the scoring metric. The metric we are specifying, accuracy, is the default for classification tasks, so we don’t really need to set it here. We do so, however, to learn how to set the performance metric to prepare to use more sophisticated performance metrics in the future.

10.7.2 Logistic Regression

The code cell below sets up the column transformer, pipeline, parameter grid, cross validation splitter, and grid search.

# Make column transformer
ct_lreg = ColumnTransformer([
    ('poly_std', Pipeline([('poly', PolynomialFeatures(degree = 3, include_bias = False)),
                            ('std', StandardScaler())]), num_features),
    ('one_hot', OneHotEncoder(handle_unknown = 'ignore'), cat_features)
])

# Make pipeline to combine pre-processing with the estimator
pipe_lreg = Pipeline([
    ('pre', ct_lreg),
    ('lreg', LogisticRegression(max_iter = 2_000))
])

# Define initial parameter grid
pgrid_lreg = {'lreg__C': np.logspace(-3, 3, 15)}

# Set up a splitter for the inner splits within the grid search
inner_cv = RepeatedStratifiedKFold(
    n_splits=8,
    n_repeats=2,
    random_state=102
)


# Define grid search
gs_lreg = GridSearchCV(estimator = pipe_lreg, 
                       param_grid = pgrid_lreg, 
                       cv = inner_cv, 
                       n_jobs = -1, 
                       scoring = 'accuracy',
                       return_train_score = True)

10.7.3 Run the grid search once to check if the parameter range is appropriate

Below we run the grid search once to check to see if the best parameter setting is within the range we chose when we set up the grid search.

Important limitation: We are inspecting grid-search results using the full dataset to help choose parameter ranges before running nested cross-validation. This means the eventual outer test folds have indirectly influenced model development. A fully independent nested-CV evaluation would specify the parameter ranges without using those outer test observations.

gs_lreg.fit(X, y)

GridSearchCV(cv=RepeatedStratifiedKFold(n_repeats=2, n_splits=8, random_state=102),
             estimator=Pipeline(steps=[('pre',
                                        ColumnTransformer(transformers=[('poly_std',
                                                                         Pipeline(steps=[('poly',
                                                                                          PolynomialFeatures(degree=3,
                                                                                                             include_bias=False)),
                                                                                         ('std',
                                                                                          StandardScaler())]),
                                                                         ['East_km',
                                                                          'North_km',
                                                                          'Rain_mm']),
                                                                        ('one_hot',
                                                                         OneHotEncoder(handle_unknown='ignore'),
                                                                         ['Servic...
                                       ('lreg',
                                        LogisticRegression(max_iter=2000))]),
             n_jobs=-1,
             param_grid={'lreg__C': array([1.00000000e-03, 2.68269580e-03, 7.19685673e-03, 1.93069773e-02,
       5.17947468e-02, 1.38949549e-01, 3.72759372e-01, 1.00000000e+00,
       2.68269580e+00, 7.19685673e+00, 1.93069773e+01, 5.17947468e+01,
       1.38949549e+02, 3.72759372e+02, 1.00000000e+03])},
             return_train_score=True, scoring='accuracy')

In a Jupyter environment, please rerun this cell to show the HTML representation or trust the notebook.
On GitHub, the HTML representation is unable to render, please try loading this page with nbviewer.org.

gs_lreg.best_params_

{'lreg__C': np.float64(0.0517947467923121)}

This time we will use matplotlib code to create the validation curve plot. It is a bit easier than plotnine, because we don’t have to create a subset of the results DataFrame and melt it to long format before making the plot.

# Extract the grid search results to a DataFrame
results_lreg = pd.DataFrame(gs_lreg.cv_results_)

# Create the validation curve plot with matplotlib
plt.plot(results_lreg['param_lreg__C'], 
         results_lreg['mean_train_score'])
plt.plot(results_lreg['param_lreg__C'], 
         results_lreg['mean_test_score'])
plt.axhline(y = max(results_lreg['mean_test_score']),
           color = 'gray',
           linestyle = ':')
plt.axvline(x = results_lreg['param_lreg__C'][results_lreg['mean_test_score'].idxmax()],
           color = 'gray',
           linestyle = ':')
plt.title("Validation Curve Graph for Logistic Regression")
plt.xscale('log')
plt.xlabel("C Parameter")
plt.ylabel("Accuracy")
plt.legend(["Train Accuracy", "Validation Accuracy", "Max Validation Accuracy"], 
           loc="best");

The initial parameter range we chose for the logistic regression grid searches is fine, as the best parameter setting is near the middle of the range. It looks like the peak performance is between \(10^{-2}\) and \(10^1\), so we will focus the search on that range of C values by re-defining the parameter grid and grid search before running the nested cross validation.

# Re-define the parameter grid
pgrid_lreg_adj = {'lreg__C': np.logspace(-2, 1, 15)}

# Re-define the grid search
gs_lreg_adj = GridSearchCV(estimator = pipe_lreg, 
                           param_grid = pgrid_lreg_adj, 
                           cv = inner_cv, 
                           n_jobs = -1, 
                           scoring = 'accuracy')

Now we can run the nested cross validation. A nested cross validation is implemented with a cross_validate() function in which the estimator is a grid search. A timer has been added to the code, using the perf_counter() function from the time module.

# Start the timer
start = time.perf_counter()

# Set up the splitting strategy for the outer splits in the nested cv
outer_cv = StratifiedKFold(
    n_splits=8,
    random_state=204,
    shuffle = True
)

# Perform the nested cross validation for logistic regression
lreg_gen_scores = cross_validate(estimator = gs_lreg_adj, 
                                 X = X, 
                                 y = y,
                                 cv = outer_cv,
                                 scoring = 'accuracy')
# Stop the timer
end = time.perf_counter()

# Print out message with the elapsed time
print(f"Elapsed time is {end - start:0.2f} seconds")

Elapsed time is 6.00 seconds

Next, we calculate the mean and standard deviation of test-fold performance for logistic regression. These represent our estimate of generalization performance for the logistic regression algorithm on this data. After we go through a similar process for K-nn we can compare the estimated generalization performance of logistic regression and K-nn on the deliveries data.

print(f"lreg mean test score: {np.mean(lreg_gen_scores['test_score'])}")
print(f"lreg std test score: {np.std(lreg_gen_scores['test_score'])}")

lreg mean test score: 0.6012
lreg std test score: 0.01470238075959127

10.7.4 K-Nearest Neighbors

Now we will conduct a nested cross validation to estimate generalization performance for the K-Nearest Neighbors model on the deliveries data.

The code cell below sets up the column transformer, pipeline, parameter grid, cross validation splitter, and grid search.

# Make column transformer to do preprocessing
ct_knn = ColumnTransformer([
    ('mm_scale', MinMaxScaler(), num_features),
    ('one_hot', OneHotEncoder(handle_unknown = 'ignore'), cat_features)
])

# Make pipeline to combine preprocessing with estimator
pipe_knn = Pipeline([
    ('pre', ct_knn),
    ('knn', KNeighborsClassifier())
])

# Define parameter grid
pgrid_knn = {'knn__n_neighbors': range(1, 100, 4)}

# Set up a splitter for the inner splits within the grid search
inner_cv = RepeatedStratifiedKFold(
    n_splits=8,
    n_repeats=2,
    random_state=102
)

# Define grid search
gs_knn = GridSearchCV(estimator = pipe_knn, 
                      param_grid = pgrid_knn, 
                      cv = inner_cv, 
                      n_jobs = -1, 
                      scoring = 'accuracy',
                      return_train_score = True)

10.7.5 Run the grid search once to check if the parameter range is appropriate

Below we run the grid search once to check to see if the best parameter setting is within the range we chose when we set up the grid search.

gs_knn.fit(X, y)

GridSearchCV(cv=RepeatedStratifiedKFold(n_repeats=2, n_splits=8, random_state=102),
             estimator=Pipeline(steps=[('pre',
                                        ColumnTransformer(transformers=[('mm_scale',
                                                                         MinMaxScaler(),
                                                                         ['East_km',
                                                                          'North_km',
                                                                          'Rain_mm']),
                                                                        ('one_hot',
                                                                         OneHotEncoder(handle_unknown='ignore'),
                                                                         ['Service_Type'])])),
                                       ('knn', KNeighborsClassifier())]),
             n_jobs=-1, param_grid={'knn__n_neighbors': range(1, 100, 4)},
             return_train_score=True, scoring='accuracy')

In a Jupyter environment, please rerun this cell to show the HTML representation or trust the notebook.
On GitHub, the HTML representation is unable to render, please try loading this page with nbviewer.org.

gs_knn.best_params_

{'knn__n_neighbors': 37}

# Extract the grid search results to a DataFrame
results_knn = pd.DataFrame(gs_knn.cv_results_)

# Create the validation curve plot with matplotlib
plt.plot(results_knn['param_knn__n_neighbors'], 
         results_knn['mean_train_score'])
plt.plot(results_knn['param_knn__n_neighbors'], 
         results_knn['mean_test_score'])
plt.axhline(y = max(results_knn['mean_test_score']),
           color = 'gray',
           linestyle = ':')
plt.axvline(x = results_knn['param_knn__n_neighbors'][results_knn['mean_test_score'].idxmax()],
           color = 'gray',
           linestyle = ':')
# Flip the x-axis limits so that flexibility increases from left to right
plt.xlim(max(results_knn['param_knn__n_neighbors']), 
         min(results_knn['param_knn__n_neighbors']))
plt.title("Validation Curve Graph for K-nn")
plt.xlabel("Number of Neighbors")
plt.ylabel("Accuracy")
plt.legend(["Train Accuracy", "Validation Accuracy", "Max Validation Accuracy"], 
           loc="best");

The initial parameter range we chose for the grid search looks a little wide. We will use a narrower range with more closely spaced values to examine the region of stronger validation performance.

# Re-define parameter grid
pgrid_knn_adj = {'knn__n_neighbors': range(10, 66)}

# Re-define grid search
gs_knn_adj = GridSearchCV(estimator = pipe_knn, 
                      param_grid = pgrid_knn_adj, 
                      cv = inner_cv, 
                      n_jobs = -1, 
                      scoring = 'accuracy')

Now we can run the nested cross validation. A nested cross validation is implemented with a cross_validate() function in which the estimator is a grid search. A timer has been added to the code, using the perf_counter() function from the time module.

# Start the timer
start = time.perf_counter()

# Set up a splitter for the outer splits in the nested cv
outer_cv = StratifiedKFold(
    n_splits=8,
    random_state=204,
    shuffle = True
)

# Perform the nested cross validation for K-nn
knn_gen_scores = cross_validate(estimator = gs_knn_adj, 
                                X = X, 
                                y = y,
                                cv = outer_cv,
                                scoring = 'accuracy')

# Stop the timer
end = time.perf_counter()

# Print message with the elapsed time
print(f"Elapsed time is {end - start:0.2f} seconds")

Elapsed time is 15.26 seconds

Next, we calculate the mean and standard deviation of test-fold performance for K-nn. These represent our estimate of generalization performance for the K-nn algorithm on this data. We can compare them to the estimated generalization performance of logistic regression on the deliveries data.

print(f"knn mean test score: {np.mean(knn_gen_scores['test_score'])}")
print(f"knn std test score: {np.std(knn_gen_scores['test_score'])}")

knn mean test score: 0.6968
knn std test score: 0.011027239001672186

10.7.6 Perform grid search on model with the highest mean outer-CV score to set parameters

As we saw in the one-outer split model-building process, the K-nn model performs better on this data. Let’s go ahead and move forward with the K-nn model.

The final step in our model building process is to run a new grid search for the best model using all our labeled data. When it finishes the grid search will retrain the model with the best parameters on all the labeled data, and it will then function as an estimator that we can use to make predictions for other delivery data.

final_search = clone(gs_knn_adj)

final_search.fit(X, y)

final_model = final_search.best_estimator_

Now we can use final_model to make predictions on new data