Regularization and tuning

At the end of Bending the line, the degree-4 model sent a 1978 Buick to 477 mpg, and I promised to keep the columns and put a price on large coefficients instead of picking a degree by hand.

A price on large coefficients

Least squares picks the coefficients with the smallest sum of squared errors. Ridge regression adds the squared coefficients to the number it minimizes:

The first sum is the familiar sum of squared errors. The second grows with every coefficient, so a model that wants a coefficient of 1669 has to earn it with a big drop in the errors. alpha decides how much a coefficient costs compared to an error. At 0, you get plain least squares back. With a large enough alpha, every coefficient shrinks toward 0, and the model predicts the average car for everything. The penalty doesn’t touch the intercept.

The penalty depends on the size of a coefficient, and the size depends on the column’s units. Weight in pounds needs a tiny coefficient, weight in tons a 2,000 times larger one. So without scaling, Ridge would punish some columns for their units. The scaler goes inside the pipeline, and with curves, it needs a second one after PolynomialFeatures, because squares and products of standardized columns aren’t standardized anymore.

from sklearn.linear_model import Ridge

def penalized(degree, model):
    return make_pipeline(
        ColumnTransformer([
            ('numeric', StandardScaler(), numeric),
            ('origin', OneHotEncoder(drop=['USA'], sparse_output=False), ['origin']),
        ], verbose_feature_names_out=False),
        PolynomialFeatures(degree, include_bias=False),
        StandardScaler(),
        model,
    )

Ridge on degree 4

Here is the degree-4 model, first without a penalty and then with growing values of alpha, next to the largest coefficient it ends up with.

ridge_scores = {}
for alpha in [0, 0.1, 1, 10, 100, 1000]:
    model = penalized(4, Ridge(alpha=alpha) if alpha else LinearRegression())
    ridge_scores[alpha] = {
        'fold RMSE': fold_rmse(model, X_named).mean(),
        'largest coefficient': np.abs(model.fit(X_named, y_train)[-1].coef_).max(),
    }

pd.DataFrame(ridge_scores).T.rename_axis('alpha').to_markdown()

Without a penalty, this pipeline scores 80 mpg on the folds. In Bending the line, it was 70. The only thing I changed was the second scaler, and I didn’t expect it to matter at all. It wouldn’t, if the data pinned the coefficients down. With 329 columns and at most 278 cars in a fold, many sets of coefficients fit the training cars equally well, and as When the answer isn’t unique showed, scikit-learn returns the one with the smallest coefficients. “Smallest” depends on the units of each column, so rescaling the columns picks a different solution, and that one misses the validation years by 10 mpg more. Then a penalty of 0.1 takes the RMSE from 80 to 4.7, and around 10, it’s within 0.06 mpg of degree 2 without a penalty (3.38). At 1000, it’s too stiff again, and we’re back on the other side of the bias-variance trade-off.

no_penalty_errors = fold_residuals(penalized(4, LinearRegression()), X_named)
ridge_errors = fold_residuals(penalized(4, Ridge(alpha=10)), X_named)

pd.DataFrame({
    'car name': X_train.loc[worst, 'car name'],
    'actual': y_train[worst],
    'no penalty': y_train[worst] - no_penalty_errors[worst],
    'Ridge, alpha=10': y_train[worst] - ridge_errors[worst],
}).to_markdown()

The Buick comes back from 478 mpg to 18.2, half an mpg from the truth. The Mercedes-Benz 300D is still 8 mpg too pessimistic, but at least it no longer produces fuel. It’s a diesel, and the dataset still has no column for that.

The 280S also shows what the second scaler did. In Bending the line, it was predicted at +311 mpg. Here, with only the second scaler changed, the model picked a different solution and predicts −359.

Lasso: zeros instead of small numbers

Lasso penalizes the absolute values of the coefficients instead of their squares. Scikit-learn also divides the errors by twice the number of rows, which is why Lasso’s useful values of alpha are much smaller than Ridge’s:

Under a squared penalty, shrinking a coefficient from 0.1 to 0 saves almost nothing (0.01), so Ridge leaves small coefficients small. Under an absolute penalty, the same step saves the full 0.1, so Lasso pushes every coefficient that doesn’t reduce the errors enough all the way to 0. The column drops out of the model, and Lasso has picked the features for us.

In Choosing what to keep, I picked horsepower, weight, and acceleration by hand. Let’s see what Lasso picks from the same columns, without curves, as alpha grows.

from sklearn.linear_model import Lasso

lasso_path = {}
for alpha in [0.01, 0.1, 0.5, 1, 2]:
    # A high max_iter, because the solver needs more rounds at small alphas
    model = penalized(1, Lasso(alpha=alpha, max_iter=100_000)).set_output(transform='pandas')
    coefficients = pd.Series(model.fit(X_named, y_train)[-1].coef_, index=model[-1].feature_names_in_)
    lasso_path[alpha] = {
        'fold RMSE': fold_rmse(model, X_named).mean(),
        'kept columns': ', '.join(coefficients[coefficients != 0].index),
    }

pd.DataFrame(lasso_path).T.rename_axis('alpha').to_markdown()

At alpha=1, Lasso keeps three columns, and they aren’t mine. It traded acceleration for displacement, one of the columns I dropped (its coefficient flipped sign between folds). Lasso doesn’t know about that. All it sees is that displacement, horsepower, and weight correlate at 0.84 to 0.95, so whichever of them fits the training rows slightly better stays in. The path already shows the trade: at 0.01, Lasso keeps acceleration and drops displacement, and at 0.1, it’s the other way around. Here are the columns it keeps at alpha=1 when it trains on each fold’s years:

lasso_1 = penalized(1, Lasso(alpha=1, max_iter=100_000)).set_output(transform='pandas')
kept_per_fold = {}
for (train_rows, _), year in zip(folds, range(75, 80)):
    fitted = clone(lasso_1).fit(X_named.iloc[train_rows], y_train.iloc[train_rows])
    kept_per_fold[year] = ', '.join(fitted[-1].feature_names_in_[fitted[-1].coef_ != 0])

pd.Series(kept_per_fold, name='kept columns').rename_axis('validation year').to_markdown()

Trained on the years up to 1974, 1975, 1976, or 1977, Lasso keeps cylinders, displacement, and weight. Once the 1978 cars join, it swaps cylinders for horsepower. That’s Rule 15 again, with Lasso doing the choosing.

The 4.06 at alpha=1 also includes coefficients pulled toward 0, so it can’t tell us whether Lasso’s three columns are worse than mine. Without any penalty, the two sets score:

pd.Series({
    "Lasso's three": fold_rmse(LinearRegression(), X_train[['displacement', 'horsepower', 'weight']]).mean(),
    'my three': fold_rmse(LinearRegression(), X_train[['horsepower', 'weight', 'acceleration']]).mean(),
}, name='fold RMSE').to_markdown()

0.01 mpg apart. Lasso and I picked different columns from the same correlated block, and the folds can’t tell the two choices apart.

On the degree-4 columns, Lasso does the pruning on a bigger scale:

lasso_4 = penalized(4, Lasso(alpha=0.05, max_iter=100_000))
fold_rmse(lasso_4, X_named).mean(), (lasso_4.fit(X_named, y_train)[-1].coef_ != 0).sum()

(np.float64(3.364243058328075), np.int64(36))

Lasso set 293 of the 329 coefficients to exactly 0, and the model that sent the Buick to 478 mpg now scores 3.36 on the folds. Apparently, it only needed 36 of them.

If you can’t decide between the two penalties, ElasticNet uses both.

Choosing alpha

So far, I picked the values of alpha by trying a few. Scikit-learn can do the trying. GridSearchCV takes a pipeline and a grid of settings, scores every combination with cross-validation, and keeps the best one. These settings are hyperparameters: values we choose before fitting, unlike the coefficients, which fit() chooses. Choosing them is step 8 of the machine learning process.

We pass our year-by-year folds as cv, so every combination is judged the same way as everything else in this book. The degree goes into the grid too.

from sklearn.model_selection import GridSearchCV

search = GridSearchCV(
    penalized(2, Ridge()),
    {'polynomialfeatures__degree': [1, 2, 3, 4], 'ridge__alpha': np.logspace(-2, 3, 11)},
    cv=folds,
    scoring='neg_root_mean_squared_error',
).fit(X_named, y_train)

search.best_params_, -search.best_score_

({'polynomialfeatures__degree': 3, 'ridge__alpha': np.float64(3.1622776601683795)}, np.float64(3.292990721299124))

For Ridge and Lasso alone, RidgeCV and LassoCV would be the more proper tools. LassoCV computes the whole range of alpha in one pass, and RidgeCV has a fast shortcut when you let it use its own cross-validation. Inside our pipeline, though, they would tune only the last step, after the scalers had already learned from the validation years. I use GridSearchCV because it refits the whole pipeline in every fold and works the same way for every model and every setting, including the degree.

The search tried 44 combinations (four degrees times 11 values of alpha) on five folds each. The winner is degree 3 with alpha at 3.16: 3.29 on the folds, against 3.33 for the pipeline. Before you celebrate, look at who came second.

results = pd.DataFrame(search.cv_results_).sort_values('rank_test_score')
pd.DataFrame({
    'degree': results['param_polynomialfeatures__degree'],
    'alpha': results['param_ridge__alpha'],
    'fold RMSE': -results['mean_test_score'],
}).head(5).to_markdown(index=False)

The top five sit within 0.06 mpg of each other. An average over five folds hides how the two models trade places, so here is the winner against the pipeline from Improving the model, year by year:

pd.DataFrame({
    'tuned Ridge': fold_rmse(search.best_estimator_, X_named),
    'final pipeline': fold_rmse(final, X_train),
}).assign(difference=lambda folds_table: folds_table['final pipeline'] - folds_table['tuned Ridge']).to_markdown()

The tuned Ridge wins four years by 0.14 to 0.33 mpg and loses 1979 by 0.74, and 1979 is the year closest to the cars we will have to predict. The lead changes more from one year to the next than its average size of 0.04. Compare that with the reciprocal target in Improving the model, which won four years out of five and lost the fifth by 0.07. The best of 44 tries also gets more help from the winner’s curse than the best of the six degrees in Bending the line. (The proper measurement is nested cross-validation: run the whole search again inside each fold, then score the winner on the year the search never saw. I’m skipping it, because its answer wouldn’t change what we do next.) I’m keeping the pipeline.

Rule 19: The best score of a search is the best of many tries, so it flatters the winner. Compare the winner with the model you already have, fold by fold, and if its lead ever changes between two neighboring folds by more than its average size, call it a tie. The tie goes to the model you can explain.

So the penalty fixes the blowups without me picking a degree by hand, and it still lands at 3.3 on the folds, like everything since Improving the model. Every one of these models adds up columns multiplied by coefficients, and none of them has the year. The next chapter stops adding things up and predicts from the cars that look most similar.