In this unit, we will cover the idea of subset selection.

- Know what subset selection is and when to use it.

In our previous discussion of model performance, you were introduced to the problem of **overfitting.** Here is a quick repeat. The idea is that you generally want to make inferences or predictions **not for the data you used to fit your model, but for new/different data.**

If the goal was to fit the data you have as well as possible, you could always achieve that with a complex model. If you have N observations, you can find a polynomial of degree N-1 that hits all of your data points. However, as discussed when we talked about the **bias-variance trade-off**, a model that gets too close to the data used for fitting will generally not perform well when applied to new data.

You want a model that gets the trade-off between high bias (model not flexible enough to capture all the patterns seen in the data) and high variance (model too flexible, will capture noise as well) right.

I showed the figure below previously. I’m repeating it here as a refresher. The yellow model (straight line), in this example, has too much bias. The model is too simple to capture all the important patterns, and it **underfits** the data. The performance metric (mean square error here) is high for both the training and test/CV data (yellow squares on the left of the right panel.) The green model (very wiggly line) has too much variance. It goes too close to the data used for fitting and doesn’t generalize to new data, it **overfits** the data. This reflects itself in good performance (low MSE) for the train data, but high MSE for the test/CV data (green squares on left of the left panel) In this example, the blue model has the optimal model complexity. It captures the overall pattern in the data, without trying to fit the noise/variability. This manifests itself with the best performance on the test/CV data (blue squares).

We talked about cross-validation and how it can be used to help you determine if your model overfits (or underfits). What we haven’t discussed is how to find a model that gets the trade-off between bias and variance right. Two common techniques are subset selection and regularization.

The general idea for these strategies is that if you let the model use all the predictors you have (this is called a *full* or *saturated* model), it might be too big and overfit the data. You can get some hint of that if the cross-validated performance measure for a full model is not better than that for a null model or a single-predictor model. The latter ones are likely too simple and underfit (unless none or a single predictor describes your data), while the full model is likely too flexible and overfits. So you want to find a model that only uses the important/useful predictors.

You might have already seen subset selection in another statistics class. If that is the case, this is a refresher. Though you might not have used the approaches I describe here to perform subset selection. Subset selection is also called *variable selection* or *feature selection*. The overall idea is that you try models with different numbers of predictors included/excluded, you evaluate its performance (using cross-validation to get an honest estimate of model performance on new data), and pick the reduced/sub-model with the best performance. The performance is measured using whatever metric you decided to optimize on (SSR, RMSE, Accuracy, F1 score, etc.). There are different ways you can go about trying different sub-models, as described in the following sections.

In forward selection, you start with the null-model (no predictors), then evaluate all single-predictor models. The best-performing 1-predictor model is chosen, and you build all 2-predictor models on top of that chosen 1-predictor model. You choose the best 2-predictor model, then go on to add a third, and do that until adding any further predictors does not lead to a model that performs better (**again, measured by CV on the test set**) than the smaller model. E.g., you might find that adding any of the remaining 3rd predictors does not lead to a better-performing model than the 2-predictor model you already have. Then you stop here and pick the 2-predictor model as your best.

Let’s consider an example. Say you want to predict `BMI`

based on `age`

, `sex`

, `calorie intake`

, and `exercise level`

. Since BMI is continuous, we can consider a linear regression model, and we might try to minimize RMSE between model-predicted BMI and actual measurements. We start by computing RMSE for the null model. Then we compute RMSE for all 1-predictor models. Assume that all 1-predictor models have RMSE lower than the null model, and the lowest RMSE (all computed using cross-validation) of the 1-predictor models is the one that includes `calorie intake`

. Next, you try all 2-predictor models that include `calorie intake`

. Among those, a model with `calorie intake + exercise level`

has the lowest RMSE, and it’s also lower than the model with `calorie intake`

only. Next, you try all 3-predictor models that include the 2 chosen predictors. In this example, only 2 of those 3-predictor models are possible, namely `calorie intake + exercise level + age`

and `calorie intake + exercise level + sex`

. Let’s assume both of those models have RMSE that is larger than the 2 predictor model. Note, it is important to compute RMSE through cross-validation since the RMSE evaluated on the data used for fitting will always be lower for the bigger model. We thus found our best model, namely one that includes the 2 predictors `calorie intake + exercise level`

.

Backward selection is essentially the reverse to forward selection: You start with the model with all predictors, then evaluate all models with one predictor dropped. The smaller model with the best performance is your choice. You then drop each predictor at a time from that model. You continue until dropping predictors does not improve performance anymore. Then you keep the model with the best performance. Note that this method does not work on datasets that have more predictors than observations (e.g., most -omics type data).

Let’s consider the above example again. We now start with a model that includes all 4 predictors and compute cross-validated RMSE. Next, we drop each predictor and try all 3-predictor models. Let’s say the one that dropped `exercise level`

is the 3-predictor model with the lowest RMSE, lower than all other 3-predictor models, and lower than the full model. Next, we drop each variable from the current best candidate 3-predictor model, which is `age + sex + calorie intake`

. All of the 2-predictor models have RMSE larger than the 3-predictor model. We thus found our final model, with the 3 predictors `age + sex + calorie intake`

.

Note that it is also possible to do “forward-backward” selection: you can choose to start with the full or reduced model, but the method is allowed to add/remove predictors at every step depending on what improves performance most.

You might have noticed that in my example, forward selection and backward selection do not necessarily produce the same model. The reason for that is that neither forward- or backward-selection try all possible models (predictor combinations). Both approaches use what is called a *greedy* strategy. Once they included or excluded a predictor, it can never be dropped or included. Once you picked your 1-predictor model, all subsequent models will consist of that 1 predictor model - even if there might be a better performing 2-predictor model that doesn’t include the 1st predictor you picked. The same holds for dropping predictors in backward selection.

Also, as in my example, often (but not always), the final model determined through backward selection tends to include more predictors than a model determined through forward selection. In such a case, you have to decide which model to choose. Alternatively, you can try a different way of choosing the model. One approach is to try models with all combinations of predictors. This is called *exhaustive subset selection* or *best subset selection* (although note that “best subset selection” can also refer to a different method which still hasn’t gained much traction). If you have only a few predictors, you can try every combination of predictors. However, once you have more than 5-10 predictors (depending on the number of observations), trying every combination of sub-models would be too time-consuming. Therefore, most of the time, such an exhaustive feature selection is not feasible.

There have recently been a number of approaches that try to be somewhat smarter in finding the best model, without having to try every combination. The `caret`

package has several methods, such as Genetic Algorithms and Simulated Annealing. Those methods often do a good job finding a very good sub-model in a reasonable amount of time. While the only guarantee that you found the single best model is to try them all using exhaustive subset selection, these modern approaches usually find a model that is either the best or at least close to the best while doing so rather efficiently. See sections 18-22 of the caret manual for more information if you are interested. The `mlr/mlr3`

packages also have similar such selection methods (in fact, I find the `mlr/mlr3`

ones better). `tidymodels`

currently does not have any subset selection procedure implemented (see more below).

All the approaches above use cross-validated performance of some metric to determine which predictors should be included in the model. Another, still commonly used strategy is to do subset selection manually. The analyst usually starts with a simple model containing just one predictor (the most important one), and then decides based on model results (e.g., systematic deviations in the residuals) and/or based on scientific knowledge, if adding additional predictors is justified. For a skilled analyst with a lot of subject matter expertise, this can lead to models that are all-around very good. They have good performance, make scientific sense, and are still simple. Unfortunately, this requires a lot of knowledge, time and careful analysis. Also, the analyst needs to honestly report any model they tried, not just the one they end up with. Otherwise they are using and not reporting researcher’s degrees of freedom, which leads to overfitting. Most analyses are not done that way.

I mentioned this before in various places and it’s a theme throughout the course, that unless data is collected and analyzed using pre-specified protocols (e.g., a clinical trial for FDA approval) it is almost always the case that scientists have so much flexibility during the data collection and analysis stages that p-values and similar statistical measures are meaningless. The **Researcher’s degrees of freedom** issue is everywhere and essentially impossible to get rid of. For more on this topic, see the “Pitfalls” section on the General Resources page.

All of this applies to manual subset selection (which is why I’m thinking about it here and adding this note), but it is more general and applies to every step in a scientific project.

Sometimes, you have predictors of particular interest, e.g., your main exposure. What if the subset selection algorithm removes those variables? That can happen. One reason is that the predictor/variable you thought was important for the outcome is in fact not. Alternatively, it could be that your main variable correlates strongly with one or a combination of other variables, and the model routine kept those other variables instead. Some packages in R that do subset selection allow you to specify certain predictors that are kept, no matter what. This way you can keep the ones you know are important, and let the computer/stats decide if the other ones are useful or not. It is also worth thinking about what the model is trying to tell you. If it throws out the variable you are interested in, it means it might not be as important in influencing the outcome as you thought. This is a situation that requires careful thinking and judgement and then based on your best judgment, you proceed.

Most commonly, subset selection is used by starting with a full model that includes each predictor/feature in linear form only. However, there is no reason why one couldn’t build a full model that includes every predictor in more complicated forms. For instance, one could start with a model that includes each term linear and quadratic, e.g.

\[Y=b_0 + b_1X_1 + b_2X_1^2 + b_3 X_2 + b_4 X_2^2 + ...\]

Interaction terms (e.g., \(X_1 X_2\)) are of course also possible. One could then perform subset selection by removing each term (though if one removes the simple term, \(b_iX_i\), it is often a good idea to remove the higher-order terms which include \(X_i\), especially interaction terms). It is rarely done, and I don’t think I’ve seen an example of this in the literature. If interactions or higher-order terms are suspected, one often uses a different type of machine learning model that more naturally allows relations between predictors and outcome that are not linear (e.g. quadratic or other complicated forms). We’ll visit some of those models soon.

You might have learned other approaches to compare model performance, e.g., using AIC or similar metrics or F-tests. As previously mentioned, I consider cross-validation the most robust and reliable method and, thus, do not discuss the others. AIC or similar might sometimes be necessary if the cross-validation approach takes too long to run, e.g., if you have a lot of data or a complicated model. I don’t see a use for any of the other methods. (But I’m happy to change my mind if someone gives me some convincing reasons 😃.)

The subset selection methods described above can be applied to any generalized linear model (GLM), such as linear, logistic and many others. For subset selection to work, variables need to show up in such a form that they can be added or removed easily. We will soon look at machine learning models where the predictor variables are part of the model in a complicated manner and can’t be easily added/removed in a stepwise manner. Thus, subset selection as described here might not work for such models. However, the similar and related term of *feature selection* does apply. The term comes from ML, the idea is similar: Try to find those features/variables that are relevant. Many of those models have some internal mechanisms that allow the model to decide that a certain predictor is not useful and thus remove it. We’ll discuss that more soon.

So I just explained in detail how subset selection works and what it tries to accomplish (make simpler models and prevent overfitting), and now I’m going to tell you to maybe avoid it 🤦?! Here is the discussion. If you go to the statistical literature or online, you will find papers and posts telling you why subset selection is bad. Here is a nice paper that lists subset selection as one: *Five myths about variable selection.*. Many similar papers and discussions exist. If you want to see more, check Wikipedia or do a Google search on `problems with stepwise selection`

. You’ll find plenty. The general criticism is that subset selection doesn’t solve the problem of overfitting. Instead, by trying lots of different models, it might make it worse. Most of the criticism is directed at doing subset selection through AIC, F-test or other criteria, not using CV. I haven’t yet seen a convincing article that explains that subset selection if using CV is equally bad (if you find one, please let me know!) Even when using CV, subset selection might still lead to overfitting, especially if the amount of data you have is not large and thus CV estimates not too reliable. Probably because of that, and also because there are other methods, such as regularization (which we’ll cover next), that accomplish similar goals often more efficiently, the `tidymodels`

folks seem to not be interested in implementing subset selection anytime soon.

Because of this, and since I don’t want to jump around different R packages more than needed, we won’t actually practice subset selection. But it’s an approach worth knowing about. While it is an approach that - in my opinion - can be useful at times, you might want to consider other methods at reducing overfitting first, e.g., the regularization approach we discuss next.

For more details on subset selection, see section 6.1. of ISLR and sections 18-22 of the caret manual.

As just mentioned, since `tidymodels`

currently doesn’t support the subset selection approaches described here, we won’t try them. The next topic, regularization, does something similar to subset selection and is often more computationally efficient. We will try some of those methods.