The process of predictive modeling

By looking at some of the different characterizations of models, we've already hinted at various steps of the predictive modeling process. In this section, we will present these steps in a sequence and make sure we understand how each of these contributes to the overall success of the endeavor.

Defining the model's objective

In a nutshell, the first step of every project is to figure out precisely what the desired outcome is, as this will help steer us to make good decisions throughout the course of the project. In a predictive analytics project, this question involves drilling into the type of prediction that we want to make and understanding the task in detail. For example, suppose we are trying to build a model that predicts employee churn for a company. We first need to define this task precisely, while trying to avoid making the problem overly broad or overly specific. We could measure churn as the percentage of new full time hires that defect from the company within their first six months. Notice that once we properly define the problem, we have already made some progress in thinking about what data we will have to work with. For example, we won't have to collect data from part-time contractors or interns. This task also means that we should collect data from our own company only, but at the same time recognize that our model might not necessarily be applicable to make predictions for the workforce of a different company. If we are only interested in churn, it also means that we won't need to make predictions about employee performance or sick days (although it wouldn't hurt to ask the person for whom we are building the model, to avoid surprises in the future).

Once we have a precise enough idea of the model we want to build, the next logical question to ask is what sort of performance we are interested in achieving, and how we will measure this. That is to say, we need to define a performance metric for our model and then a minimum threshold of acceptable performance. We will go into substantial detail on how to assess the performance of models in this book. For now, we want to emphasize that, although it is not unusual to talk about assessing the performance of a model after we have trained it on some data, in practice it is important to remember that defining the expectations and performance target for our model is something that a predictive modeler should discuss with the stakeholders of a project at the very beginning. Models are never perfect and it is easy to spiral into a mode of forever trying to improve performance. Clear performance goals are not only useful in guiding us to decide which methods to use, but also in knowing when our model is good enough.

Finally, we also need to think about the data that will be available to us when the time comes to collect it, and the context in which the model will be used. For example, suppose we know that our employee churn model will be used as one of the factors that determine whether a new applicant in our company will be hired. In this context, we should only collect data from our existing employees who were available before they were hired. We cannot use the result of their first performance review, as this data won't be available for a prospective applicant.

Collecting the data

Training a model to make predictions is often a data-intensive venture, and if there is one thing that you can never have too much of in this business, it is data. Collecting the data can often be the most time- and resource-consuming part of the entire process, which is why it is so critical that the first step of defining the task and identifying the right data to be collected is done properly. When we learn about how a model such as logistic regression works we often do this by way of an example dataset and this is largely the approach we'll follow in this book. Unfortunately, we don't have a way to simulate the process of collecting the data, and it may seem that most of the effort is spent on training and refining a model. When learning about models using existing datasets, we should bear in mind that a lot of effort has usually gone into collecting, curating, and preprocessing the data. We will look at data preprocessing more closely in a subsequent section.

While we are collecting data, we should always keep in mind whether we are collecting the right kind of data. Many of the sanity checks that we perform on data during preprocessing also apply during collection, in order for us to spot whether we have made a mistake early on in the process. For example, we should always check that we measure features correctly and in the right units. We should also make sure that we collect data from sources that are sufficiently recent, reliable, and relevant to the task at hand. In the employee churn model we described in the previous section, as we collect information about past employees we should ensure that we are consistent in measuring our features. For example, when measuring how many days a person has been working in our company, we should consistently use either calendar days or business days. We must also check that when collecting dates, such as when a person joined or left the company, we invariably either use the US format (month followed by day) or the European format (day followed by month) and do not mix the two, otherwise a date like 03/05/2014 will be ambiguous. We should also try to get information from as broad a sample as possible and not introduce a hidden bias in our data collection. For example, if we wanted a general model for employee churn, we would not want to collect data from only female employees or employees from a single department.

How do we know when we have collected enough data? Early on when we are collecting the data and have not built and tested any model, it is impossible to tell how much data we will eventually need, and there aren't any simple rules of thumb that we can follow. We can, however, anticipate that certain characteristics of our problem will require more data. For example, when building a classifier that will learn to predict from one of three classes, we may want to check whether we have enough observations representative of each class.

The greater the number of output classes we have, the more data we will need to collect. Similarly, for regression models, it is also useful to check that the range of the output variable in the training data corresponds to the range that we would like to predict. If we are building a regression model that covers a large output range, we will also need to collect more data compared to a regression model that covers a smaller output range under the same accuracy requirements.

Another important factor to help us estimate how much data we will need, is the desired model performance. Intuitively, the higher the accuracy that we need for our model, the more data we should collect. We should also be aware that improving model performance is not a linear process. Getting from 90 to 95% accuracy can often require more effort and a lot more data, compared to making the leap from 70 to 90%. Models that have fewer parameters or are simpler in their design, such as linear regression models, often tend to need less data than more complex models such as neural networks. Finally, the greater the number of features that we want to incorporate into our model, the greater the amount of data we should collect. In addition, we should be aware of the fact that this requirement for additional data is also not going to be linear. That is to say, building a model with twice the number of features often requires much more than twice the amount of original data. This should be readily apparent, if we think of the number of different combinations of inputs our model will be required to handle. Adding twice the number of dimensions results in far more than twice the number of possible input combinations. To understand this, suppose we have a model with three input features, each of which takes 10 possible values. We have 10³ = 1000 possible input combinations. Adding a single extra feature that also takes 10 values raises this to 10,000 possible combinations, which is much more than twice the number of our initial input combinations.

There have been attempts to obtain a more quantifiable view of whether we have enough data for a particular dataset, but we will not have time to cover them in this book. A good place to start learning more about this area of predictive modeling is to study learning curves. In a nutshell, with this approach we build consecutive models on the same dataset by starting off with a small portion of the data and successively adding more. The idea is that if throughout this process the predictive accuracy on testing data always improves without tapering off, we probably could benefit from obtaining more data. As a final note for the data collection phase, even if we think we have enough data, we should always consider how much it would cost us (in terms of time and resources) in order to get more data, before making the choice to stop collecting and begin modeling.

Picking a model

Once we are clear on the prediction task, and we have the right kind of data, the next step is to pick our first model. To begin with, there is no best model overall, not even a best model using a few rules of thumb. In most cases, it makes sense to start off with a simple model, such as a Naïve Bayes model or a logistic regression in the case of a classification task, or a linear model in the case of regression. A simple model will give us a starting baseline performance, which we can then strive to improve. A simple model to start off with might also help in answering useful questions, such as how each feature contributes to the result, that is, how important is each feature and is the relationship with the output positively or negatively correlated. Sometimes, this kind of analysis itself warrants the production of a simple model first, followed by a more complex one, which will be used for the final prediction.

Sometimes a simple model might give us enough accuracy for the task at hand so that we won't need to invest more effort in order to give us a little bit extra. On the other hand, a simple model will often end up being inadequate for the task, requiring us to pick something more complicated. Choosing a more complex model over a simpler one is not always a straightforward decision, even if we can see that the accuracy of the complex model will be much better. Certain constraints, such as the number of features we have or the availability of data, may prevent us from moving to a more complex model. Knowing how to choose a model involves understanding the various strengths and limitations of the models in our toolkit. For every model we encounter in this book, we will pay particular attention to learning these points. In a real-world project, to help guide our decision, we often go back to the task requirements and ask a few questions, such as:

• What type of task do we have? Some models are only suited for particular tasks such as regression, classification, or clustering.
• Does the model need to explain its predictions? Some models, such as decision trees, are better at giving insights that are easily interpretable to explain why they made a particular prediction.
• Do we have any constraints on prediction time?
• Do we need to update the model frequently and is training time, therefore, important?
• Does the model work well if we have highly correlated features?
• Does the model scale well for the number of features and amount of data that we have available? If we have massive amounts of data, we may need a model whose training procedure can be parallelized to take advantage of parallel computer architectures, for example.

In practice, even if our first analysis points towards a particular model, we will most likely want to try out a number of options before making our final decision.

Pre-processing the data

Before we can use our data to train our model, we typically need to preprocess it. In this section, we will discuss a number of common preprocessing steps that we usually perform. Some of these are necessary in order to detect and resolve problems in our data, while others are useful in order to transform our data and make them applicable to the model we have chosen.

Exploratory data analysis

Once we have some data and have decided to start working on a particular model, the very first thing we'll want to do is to look at the data itself. This is not necessarily a very structured part of the process; it mostly involves understanding what each feature measures and getting a sense of the data we have collected. It is really important to understand what each feature represents and the units in which it is measured. It is also a really good idea to check the consistent use of units. We sometimes call this investigative process of exploring and visualizing our data exploratory data analysis.

An excellent practice is to use the summary() function of R on our data frame to obtain some basic metrics for each feature, such as the mean and variance, as well as the largest and smallest values. Sometimes, it is easy to spot that a mistake has been made in data collection through inconsistencies in the data. For example, for a regression problem, multiple observations with identical feature values but wildly different outputs may (depending on the application) be a signal that there are erroneous measurements. Similarly, it is a good idea to know whether there are any features that have been measured in the presence of significant noise. This may sometimes lead to a different choice of model or it may mean that the feature should be ignored.

Tip

Another useful function used to summarize features in a data frame is the describe() function in the psych package. This returns information about how skewed each feature is, as well as the usual measures of a location (such as the mean and median) and dispersion (such as the standard deviation).

An essential part of exploratory data analysis is to use plots to visualize our data. There is a diverse array of plots that we can use depending on the context. For example, we might want to create box plots of our numerical features to visualize ranges and quartiles. Bar plots and mosaic plots are useful to visualize the proportions of our data under different combinations of values for categorical input features. We won't go into further detail on information visualization, as this is a field in its own right.

Tip

R is an excellent platform to create visualizations. The base R package provides a number of different functions to plot data. Two excellent packages to create more advanced plots are lattice and ggplot2. Good references for these two, which also cover principles used to make effective visualizations, are Lattice: Multivariate Data Visualization with R and ggplot2: Elegant Graphics for Data Analysis, both of which are published by Springer under the Use R! series.

Feature transformations

Often, we'll find that our numerical features are measured on scales that are completely different from each other. For example, we might measure a person's body temperature in degrees Celsius, so the numerical values will typically be in the range of 36-38. At the same time, we might also measure a person's white blood cell count per microliter of blood. This feature generally takes values in the thousands. If we are to use these features as an input to an algorithm, such as kNN, we'd find that the large values of the white blood cell count feature dominate the Euclidean distance calculation. We could have several features in our input that are important and useful for classification, but if they were measured on scales that produce numerical values much smaller than one thousand, we'd essentially be picking our nearest neighbors mostly on the basis of a single feature, namely the white blood cell count. This problem comes up often and applies to many models, not just kNN. We handle this by transforming (also referred to as scaling) our input features before using them in our model.

We'll discuss three popular options for feature scaling. When we know that our input features are close to being normally distributed, one possible transformation to use is Z-score normalization, which works by subtracting the mean and dividing it by the standard deviation:

E(x) is the expectation or mean of x, and the standard deviation is the square root of the variance of x, written as Var(x). Notice that as a result of this transformation, the new feature will be centered on a mean of zero and will have unit variance. Another possible transformation, which is better when the input is uniformly distributed, is to scale all the features and outputs so that they lie within a single interval, typically the unit interval [0,1]:

A third option is known as the Box-Cox transformation. This is often applied when our input features are highly skewed (asymmetric) and our model requires the input features to be normally distributed or symmetrical at the very least:

As λ is in the denominator, it must take a value other than zero. The transformation is actually defined for a zero-valued λ: in this case, it is given by the natural logarithm of the input feature, ln(x). Notice that this is a parameterized transform and so there is a need to specify a concrete value of λ. There are various ways to estimate an appropriate value for λ from the data itself. Indicatively, we'll mention a technique to do this, known as cross-validation, which we will encounter later on in this book in Chapter 5, Support Vector Machines.

Note

The original reference for the Box-Cox transformation is a paper published in 1964 by the Journal of the Royal Statistical Society, titled An analysis of Transformations and authored by G. E. P. Box and D. R. Cox.

To get a feel for how these transformations work in practice, we'll try them out on the Sepal.Length feature from our iris dataset. Before we do this, however, we'll introduce the first R package that we will be working with, caret.

The caret package is a very useful package that has a number of goals. It provides a number of helpful functions that are commonly used in the process of predictive modeling, from data preprocessing and visualization, to feature selection and resampling techniques. It also features a unified interface for many predictive modeling functions and provides functionalities for parallel processing.

Note

The definitive reference for predictive modeling using the caret package is a book called Applied Predictive Modeling, written by Max Kuhn and Kjell Johnson and published by Springer. Max Kuhn is the principal author of the caret package itself. The book also comes with a companion website at http://appliedpredictivemodeling.com.

When we transform our input features on the data we use to train our model, we must remember that we will need to apply the same transformation to the features of later inputs that we will use at prediction time. For this reason, transforming data using the caret package is done in two steps. In the first step, we use the preProcess() function that stores the parameters of the transformations to be applied to the data, and in the second step, we use the predict() function to actually compute the transformation. We tend to use the preProcess() function only once, and then the predict() function every time we need to apply the same transformation to some data. The preProcess() function takes a data frame with some numerical values as its first input, and we will also specify a vector containing the names of the transformations to be applied to the method parameter. The predict() function then takes the output of the previous function along with the data we want to transform, which in the case of the training data itself may well be the same data frame. Let's see all this in action:

> library("caret")
> iris_numeric <- iris[1:4]
> pp_unit <- preProcess(iris_numeric, method = c("range"))
> iris_numeric_unit <- predict(pp_unit, iris_numeric)
> pp_zscore <- preProcess(iris_numeric, method = c("center", "scale"))
> iris_numeric_zscore <- predict(pp_zscore, iris_numeric)
> pp_boxcox <- preProcess(iris_numeric, method = c("BoxCox"))
> iris_numeric_boxcox <- predict(pp_boxcox, iris_numeric)

Tip

Downloading the example code:

You can download the example code files from your account at http://www.packtpub.com for all the Packt Publishing books you have purchased. If you purchased this book elsewhere, you can visit http://www.packtpub.com/support and register to have the files e-mailed directly to you.

We've created three new versions of the numerical features of the iris data, with the difference being that in each case we used a different transformation. We can visualize the effects of our transformations by plotting the density of the Sepal.Length feature for each scaled data frame using the density() function and plotting the results, as shown here:

Notice that the Z-score and unit interval transformations preserve the overall shape of the density while shifting and scaling the values, whereas the Box-Cox transformation also changes the overall shape, resulting in a density that is less skewed than the original.

Encoding categorical features

Many models, from linear regression to neural networks, require all the inputs to be numerical, and so we often need a way to encode categorical fields on a numerical scale. For example, if we have a size feature that takes values in the set {small, medium, large}, we may want to represent this with the numerical values 1, 2, and 3, respectively. In the case of ordered categories, such as the size feature just described, this mapping probably makes sense.

The number 3 is the largest on this scale and this corresponds to the large category, which is further away from the small category, represented by the number 1 than it is from the medium category, represented by the value 2. Using this scale is only one possible mapping, and in particular, it forces the medium category to be equidistant from the large and small categories, which may or may not be appropriate based on our knowledge about the specific feature. In the case of unordered categories, such as brands or colors, we generally avoid mapping them onto a single numerical scale. For example, if we mapped the set {blue, green, white, red, orange} to the numbers one through five, respectively, then this scale is arbitrary and there is no reason why red is closer to white and far from blue. To overcome this, we create a series of indicator features, I_i, which take the following form:

We need as many indicator features as we have categories, so for our color example, we would create five indicator features. In this case, I₁, might be:

In this way, our original color feature will be mapped to five indicator features and for every observation, only one of these indicator features takes the value 1 and the rest will be 0 as each observation will involve one color value in our original feature. Indicator features are binary features as they only take on two values: 0 and 1.

Note

We may often encounter an alternative approach that uses only n-1 binary features to encode n levels of a factor. This is done by choosing one level to be the reference level and it is indicated where each one of the n-1 binary features takes the value 0. This can be more economical on the number of features and avoids introducing a linear dependence between them, but it violates the property that all features are equidistant from each other.

Missing data

Sometimes, data contains missing values, where, for certain observations, some features were unavailable or could not properly be measured. For example, suppose that in our iris dataset, we lost the measurement for a particular observation's petal length. We would then have a missing value for this flower sample in the Petal.Length feature. Most models do not have an innate ability to handle missing data. Typically, a missing value appears in our data as a blank entry or the symbol NA. We should check whether missing values are actually present in our data, but have been erroneously assigned a value, such as 0, which is often a very legitimate feature value.

Before deciding how to handle missing data, especially when our approach will be to simply throw away observations with missing values, we should recognize that the particular values that are missing might follow a pattern. Concretely, we often distinguish between different so-called mechanisms for missing values. In the ideal Missing Completely At Random (MCAR) scenario, missing values occur independently from the true values of the features in which they occur, as well as from all other features. In this scenario, if we are missing a value for the length of a particular iris flower petal, then this occurs independently from how long the flower petal actually was and the value of any other feature, such as whether the observation was from the versicolor species or the setosa species. The Missing At Random (MAR) scenario is a less ideal situation. Here, a missing value is independent of the true value of the feature in question, but may be correlated with another feature. An example of this scenario is when missing petal length values mostly occur in the setosa samples in our iris dataset, as long as they still occur independently of the true petal length values.

In the Missing Not At Random (MNAR) scenario, which is the most problematic case, there is some sort of a pattern that explains when values might be missing based on the true values of the feature itself. For example, if we had difficulty in measuring very small petal lengths and ended up with missing values as a result, simply removing the incomplete samples would result in a sample of observations with above average petal lengths, and so our data would be biased.

There are a number of ways to handle missing values, but we will not delve deeply into this problem in this book. In the rare cases where we have missing values, we will exclude them from our datasets, but be aware that in a real project, we would investigate the source of the missing values in order to be sure that we can do this safely. Another approach is to attempt to guess or impute the missing values. The kNN algorithm itself is one way to do this by finding the nearest neighbors of a sample with a missing value in one feature. This is done by using a distance computation that excludes the dimension that contains the missing value. The missing value is then computed as the mean of the values of the nearest neighbors in this dimension.

Note

The interested reader can find a detailed treatment of how to handle missing values in Statistical Analysis with Missing Data, Second Edition, by Roderick J. A. Little and Donald B. Rubin, published by Wiley.

Outliers

Outliers are also a problem that often needs to be addressed. An outlier is a particular observation that is very far from the rest of the data in one or more of its features. In some cases, this may represent an actual rare circumstance that is a legitimate behavior for the system we are trying to model. In other cases, it may be that there has been an error in measurement. For example, when reporting the ages of people, a value of 110 might be an outlier, which could happen because of a reporting error on an actual value of 11. It could also be the result of a valid, albeit extremely rare measurement. Often, the domain of our problem will give us a good indication of whether outliers are likely to be measurement errors or not, and if so, as part of preprocessing the data, we will often want to exclude outliers from our data completely. In Chapter 2, Linear Regression, we will look at outlier exclusion in more detail.

Removing problematic features

Preprocessing a dataset can also involve the decision to drop some of the features if we know that they will cause us problems with our model. A common example is when two or more features are highly correlated with each other. In R, we can easily compute pairwise correlations on a data frame using the cor() function:

> cor(iris_numeric)
             Sepal.Length Sepal.Width Petal.Length Petal.Width
Sepal.Length    1.0000000  -0.1175698    0.8717538   0.8179411
Sepal.Width    -0.1175698   1.0000000   -0.4284401  -0.3661259
Petal.Length    0.8717538  -0.4284401    1.0000000   0.9628654
Petal.Width     0.8179411  -0.3661259    0.9628654   1.0000000

Here, we can see that the Petal.Length feature is very highly correlated with the Petal.Width feature, with the correlation exceeding 0.96. The caret package offers the findCorrelation() function, which takes a correlation matrix as an input, and the optional cutoff parameter, which specifies a threshold for the absolute value of a pairwise correlation. This then returns a (possibly zero length) vector that shows the columns to be removed from our data frame due to correlation. The default setting of cutoff is 0.9:

> iris_cor <- cor(iris_numeric)
> findCorrelation(iris_cor)
[1] 3
> findCorrelation(iris_cor, cutoff = 0.99)
integer(0)
> findCorrelation(iris_cor, cutoff = 0.80)
[1] 3 4

An alternative approach to removing correlation is a complete transformation of the entire feature space as is done in many methods for dimensionality reduction, such as Principal Component Analysis (PCA) and Singular Value Decomposition (SVD). We'll see the former shortly, and the latter we'll visit in Chapter 11, Recommendation Systems.

In a similar vein, we might want to remove features that are linear combinations of each other. By linear combination of features, we mean a sum of features where each feature is multiplied by a scalar constant. To see how caret deals with these, we will create a new iris data frame with two additional columns, which we will call Cmb and Cmb.N, as follows:

> new_iris <- iris_numeric
> new_iris$Cmb <- 6.7 * new_iris$Sepal.Length –
                  0.9 * new_iris$Petal.Width
> set.seed(68)
> new_iris$Cmb.N <- new_iris$Cmb +   
                    rnorm(nrow(new_iris), sd = 0.1)
> options(digits = 4)
> head(new_iris,n = 3)
  Sepal.Length Sepal.Width Petal.Length Petal.Width   Cmb Cmb.N
1          5.1         3.5          1.4         0.2 33.99 34.13
2          4.9         3.0          1.4         0.2 32.65 32.63
3          4.7         3.2          1.3         0.2 31.31 31.27

As we can see, Cmb is a perfect linear combination of the Sepal.Length and Petal.Width features. Cmb.N is a feature that is the same as Cmb, but with some added Gaussian noise with a mean of zero and a very small standard deviation (0.1), so that the values are very close to those of Cmb. The caret package can detect exact linear combinations of features, though not if the features are noisy, using the findLinearCombos() function:

> findLinearCombos(new_iris)
$linearCombos
$linearCombos[[1]]
[1] 5 1 4

$remove
[1] 5

As we can see, the function only suggests that we should remove the fifth feature (Cmb) from our data frame, because it is an exact linear combination of the first and fourth features. Exact linear combinations are rare, but can sometimes arise when we have a very large number of features and redundancy occurs between them. Both correlated features as well as linear combinations are an issue with linear regression models, as we shall soon see in Chapter 2, Linear Regression. In this chapter, we'll also see a method of detecting features that are very nearly linear combinations of each other.

A final issue that we'll look at for problematic features, is the issue of having features that do not vary at all in our dataset, or that have near zero variance. For some models, having these types of feature does not cause us problems. For others, it may create problems and we'll demonstrate why this is the case. As in the previous example, we'll create a new iris data frame, as follows:

> newer_iris <- iris_numeric
> newer_iris$ZV <- 6.5
> newer_iris$Yellow <- ifelse(rownames(newer_iris) == 1, T, F
> head(newer_iris, n = 3)
  Sepal.Length Sepal.Width Petal.Length Petal.Width  ZV Yellow
1          5.1         3.5          1.4         0.2 6.5   TRUE
2          4.9         3.0          1.4         0.2 6.5  FALSE
3          4.7         3.2          1.3         0.2 6.5  FALSE

The ZV column has the constant number of 6.5 for all observations. The Yellow column is a fictional column that records whether an observation had some yellow color on the petal. All the observations, except the first, are made to have this feature set to FALSE and so this is a near zero variance column. The caret package uses a definition of near zero variance that checks whether the number of unique values that a feature takes as compared to the overall number of observations is very small, or whether the ratio of the most common value to the second most common value (referred to as the frequency ratio) is very high. The nearZeroVar() function applied to a data frame returns a vector containing the features that have zero or near zero variance. By setting the saveMetrics parameter to TRUE, we can see more information about the features in our data frame:

> nearZeroVar(newer_iris)
[1] 5 6
> nearZeroVar(newer_iris, saveMetrics = T)
             freqRatio percentUnique zeroVar   nzv
Sepal.Length     1.111       23.3333   FALSE FALSE
Sepal.Width      1.857       15.3333   FALSE FALSE
Petal.Length     1.000       28.6667   FALSE FALSE
Petal.Width      2.231       14.6667   FALSE FALSE
ZV               0.000        0.6667    TRUE  TRUE
Yellow         149.000        1.3333   FALSE  TRUE

Here, we can see that the ZV column has been identified as a zero variance column (which is also by definition a near zero variance column). The Yellow column does have a nonzero variance, but its high frequency ratio and low unique value percentage make it a near zero variance column. In practice, we tend to remove zero variance columns, as they don't have any information to give to our model. Removing near zero variance columns, however, is tricky and should be done with care. To understand this, consider the fact that a model for species prediction, using our newer iris dataset, might learn that if a sample has yellow in its petals, then regardless of all other predictors, we would predict the setosa species, as this is the species that corresponds to the only observation in our entire dataset that had the color yellow in its petals. This might indeed be true in reality, in which case, the yellow feature is informative and we should keep it. On the other hand, the presence of the color yellow on iris petals may be completely random and non-indicative of species, but also an extremely rare event. This would explain why only one observation in our dataset had the yellow color in its petals. In this case, keeping the feature is dangerous because of the aforementioned conclusion. Another potential problem with keeping this feature will become apparent when we look at splitting our data into training and test sets, as well as other cases of data splitting, such as cross-validation, described in Chapter 5, Support Vector Machines. Here, the issue is that one split in our data may lead to unique values for a near zero variance column, for example, only FALSE values for our Yellow iris column.

Feature engineering and dimensionality reduction

The number and type of features that we use with a model is one of the most important decisions that we will make in the predictive modeling process. Having the right features for a model will ensure that we have sufficient evidence on which to base a prediction. On the flip side, the number of features that we work with is precisely the number of dimensions that the model has. A large number of dimensions can be the source of several complications. High dimensional problems often suffer from data sparsity, which means that because of the number of dimensions available, the range of possible combinations of values across all the features grows so large that it is unlikely that we will ever collect enough data in order to have enough representative examples for training. In a similar vein, we often talk about the curse of dimensionality. This describes the fact that, because of the overwhelmingly large space of possible inputs, data points that we have collected are likely to be far away from each other in the feature space. As a result, local methods, such as k-nearest neighbors, that make predictions using observations in the training data that are close to the point for which we are trying to make a prediction, will not work as well in high dimensions. A large feature set is also problematic in that it may significantly increase the time needed to train (and predict, in some cases) our model.

Consequently, there are two types of process that feature engineering involves. The first of these, which grows the feature space, is the design of new features based on features within our data. Sometimes, a new feature that is a product or ratio of two original features might work better. There are many ways to combine existing features into new ones, and often it is expert knowledge from the problem's particular application domain that might help guide us. In general though, this process takes experience and a lot of trial and error. Note that there is no guarantee that adding a new feature will not degrade performance. Sometimes, adding a feature that is very noisy or highly correlated with an existing feature may actually cause us to lose accuracy.

The second process in feature engineering is feature reduction or shrinkage, which reduces the size of the feature space. In the previous section on data preprocessing, we looked at how we can detect individual features that may be problematic for our model in some way. Feature selection refers to the process in which the subset of features that are the most informative for our target output is selected from the original pool of features. Some methods, such as tree-based models, have built-in feature selection, as we shall see in Chapter 6, Tree-based Methods. In Chapter 2, Linear Regression, we'll also explore methods to perform feature selection for linear models.

Another way to reduce the overall number of features, a concept known as dimensionality reduction, is to transform the entire set of features into a completely new set of features that are fewer in number. A classic example of this is Principal Component Analysis (PCA).

In a nutshell, PCA creates a new set of input features, known as principal components, all of which are linear combinations of the original input features. For the first principal component, the linear combination weights are chosen in order to capture the maximum amount of variation in the data. If we could visualize the first principal component as a line in the original feature space, this would be the line in which the data varies the most. It also happens to be the line that is closest to all the data points in the original feature space. Every subsequent principal component attempts to capture a line of maximum variation, but in such a way that the new principal component is uncorrelated with the previous ones already computed. Thus, the second principal component selects the linear combination of original input features that have the highest degree of variation in the data, while being uncorrelated with the first principal component.

The principal components are ordered naturally in descending order according to the amount of variation that they capture. This allows us to perform dimensionality reduction in a simple manner by keeping the first N components, where we choose N so that the components chosen incorporate a minimum amount of the variance from the original dataset. We won't go into the details of the underlying linear algebra necessary to compute the principal components.

Instead, we'll direct our attention to the fact that this process is sensitive to the variance and scale of the original features. For this reason, we often scale our features before carrying out this process. To visualize how useful PCA can be, we'll once again turn to our faithful iris dataset. We can use the caret package to carry out PCA. To do this, we specify pca in the method parameter of the preProcess() function. We can also use the thresh parameter, which specifies the minimum variance we must retain. We'll explicitly use the value 0.95 so that we retain 95% of the variance of the original data, but note that this is also the default value of this parameter:

> pp_pca <- preProcess(iris_numeric,  method = c("BoxCox", "center", "scale", "pca"), thresh = 0.95)
> iris_numeric_pca <- predict(pp_pca, iris_numeric)
> head(iris_numeric_pca, n = 3)
     PC1     PC2
1 -2.304 -0.4748
2 -2.151  0.6483
3 -2.461  0.3464

As a result of this transformation, we are now left with only two features, so we can conclude that the first two principal components of the numerical iris features incorporate over 95% of the variation in the data.

If we are interested in learning the weights that were used to compute the principal components, we can inspect the rotation attribute of the pp_pca object:

> options(digits = 2)
> pp_pca$rotation
               PC1    PC2
Sepal.Length  0.52 -0.386
Sepal.Width  -0.27 -0.920
Petal.Length  0.58 -0.049
Petal.Width   0.57 -0.037

This means that the first principal component, PC1, was computed as follows:

Sometimes, instead of directly specifying a threshold for the total variance captured by the principal components, we might want to examine a plot of each principal component and its variance. This is known as a scree plot, and we can build this by first performing PCA and indicating that we want to keep all the components. To do this, instead of specifying a variance threshold, we set the pcaComp parameter, which is the number of principal components we want to keep. We will set this to 4, which includes all of them, remembering that the total number of principal components is the same as the total number of original features or dimensions we started out with. We will then compute the variance and cumulative variance of these components and store it in a data frame. Finally, we will plot this in the figure that follows, noting that the numbers in brackets are cumulative percentages of variance captured:

> pp_pca_full <- preProcess(iris_numeric,  method = c("BoxCox", "center", "scale", "pca"), pcaComp = 4)
> iris_pca_full <- predict(pp_pca_full, iris_numeric)
> pp_pca_var <- apply(iris_pca_full, 2, var)
> iris_pca_var <- data.frame(Variance = round(100 * pp_pca_var / sum(pp_pca_var), 2), CumulativeVariance = round(100 * cumsum(pp_pca_var) / sum(pp_pca_var), 2))
> iris_pca_var
    Variance CumulativeVariance
PC1    73.45              73.45
PC2    22.82              96.27
PC3     3.20              99.47
PC4     0.53             100.00

As we can see, the first principal component accounts for 73.45% of the total variance in the iris dataset, while together with the second component, the total variance captured is 96.27%. PCA is an unsupervised method for dimensionality reduction that does not make use of the output variable even when it is available. Instead, it looks at the data geometrically in the feature space. This means that we cannot ensure that PCA will give us a new feature space that will perform well in our prediction problem, beyond the computational advantages of having fewer features. These advantages might make PCA a viable choice even when there is reduction in model accuracy as long as this reduction is small and acceptable for the specific task. As a final note, we should point out that weights of the principal components, often referred to as loadings, are unique within a sign flip as long as they have been normalized. In cases where we have perfectly correlated features or perfect linear combinations we will obtain a few principal components that are exactly zero.

Training and assessing the model

In our earlier discussion of parametric models, we saw that they come with a procedure to train the model using a set of training data. Nonparametric models will typically either perform lazy learning, in which case there really isn't an actual training procedure at all beyond memorizing the training data, or, as in the case of splines, will perform local computations on the training data.

Either way, if we are to assess the performance of our model, we need to split our data into a training set and a test set. The key idea is that we want to assess our model based on how we expect it to perform on unseen future data. We do this by using the test set, which is a portion (typically 15-30%) of the data we collected and set aside for this purpose and haven't used during training. For example, one possible divide is to have a training set with 80% of the observations in our original data, and a test set with the remaining 20%. The reason why we need a test set is that we cannot use the training set to fairly assess our model performance, since we fit our model to the training data and it does not represent data that we haven't seen before. From a prediction standpoint, if our goal was to maximize performance on our training data alone, then the best thing to do would be to simply memorize the input data along with the desired output values and our model would thus be a simple look-up table!

A good question to ask would be. How do we decide how much data to use for training and testing? There is a trade-off that is involved here that makes the answer to this question nontrivial. On the one hand, we would like to use as much data as possible in our training set, so that the model has more examples from which to learn. On the other hand, we would like to have a large test set so that we can test our trained model using many examples in order to minimize the variance of our estimate of the model's predictive performance. If we only have a handful of observations in our test set, then we cannot really generalize about how our model performs on unseen data overall.

Another factor that comes into play is how much starting data we have collected. If we have very little data, we may have to use a larger amount in order to train our model, such as an 85-15 split. If we have enough data, then we might consider a 70-30 split so that we can get a more accurate prediction on our test set.

To split a dataset using the caret package, we can use the createDataPartition() function to create a sampling vector containing the indices of the rows we will use in our training set. These are selected by randomly sampling the rows until a specified proportion of the rows have been sampled, using the p parameter:

> set.seed(2412)
> iris_sampling_vector <- createDataPartition(iris$Species, p = 0.8, list = FALSE)

Tip

It is good practice, when reporting the results of a statistical analysis involving a random number generation, to apply the set.seed() function on a randomly chosen but fixed number. This function ensures that the random numbers that are generated from the next function call involving random number generation will be the same every time the code is run. This is done so that others who read the analysis are able to reproduce the results exactly. Note that if we have several functions in our code that perform random number generation, or the same function is called multiple times, we should ideally apply set.seed() before each one of them.

Using our sampling vector, which we created for the iris dataset, we can construct our training and test sets. We'll do this for a few versions of the iris dataset that we built earlier when we experimented with different feature transformations:

> iris_train     <- iris_numeric[iris_sampling_vector,]
> iris_train_z   <- iris_numeric_zscore[iris_sampling_vector,]
> iris_train_pca <- iris_numeric_pca[iris_sampling_vector,]
> iris_train_labels <- iris$Species[iris_sampling_vector]
>
> iris_test        <- iris_numeric[-iris_sampling_vector,]
> iris_test_z      <- iris_numeric_zscore[-iris_sampling_vector,]
> iris_test_pca    <- iris_numeric_pca[-iris_sampling_vector,]
> iris_test_labels <- iris$Species[-iris_sampling_vector]

We are now in a position to build and test three different models for the iris dataset. These are in turn, the unnormalized model, a model where the input features have been centered and scaled with a Z-score transformation, and the PCA model with two principal components. We could use our test set in order to measure the predictive performance of each of these models after we build them; however, this would mean that, in our final estimate of unseen accuracy, we will have used the test set in the model selection, thus producing a biased estimate. For this reason, we often maintain a separate split of the data, usually as large as the test set, known as the validation set. This is used to tune model parameters, such as k in kNN, and among different encodings and transformations of the input features before using the test set to predict unseen performance. In Chapter 5, Support Vector Machines, we'll discuss an alternative to this approach known as cross-validation.

Once we have split our data, trained our model by following the relevant training procedure that it requires, and tuned our model parameters, we then have to assess its performance on the test set. Typically, we won't find the same performance on our test set as on our training set. Sometimes, we may even find that the performance we see when we deploy our model does not match what we expected to see, based on the performance on our training or test sets. There are a number of possible reasons for this disparity in performance. The first of these is that the data we may have collected may either not be representative of the process that we are modeling, or that there are certain combinations of feature inputs that we simply did not encounter in our training data. This could produce results that are inconsistent with our expectations. This situation can happen both in the real world, but also with our test set if it contains outliers, for example. Another common situation is the problem of model overfitting.

Overfitting is a problem in which some models, especially more flexible models, perform well on their training dataset, but perform significantly worse on an unseen test set. This occurs when a model matches the observations in the training data too closely and fails to generalize on unseen data. Put differently, the model is picking up on spurious details and variations in a training dataset, which are not representative of the underlying population as a whole. Overfitting is one of the key reasons why we do not choose our model based on its performance on the training data. Other sources of discrepancy between training and test data performance are model bias and variance. Together, these actually form a well-known trade-off in statistical modeling known as the bias-variance tradeoff.

The variance of a statistical model refers to how much the model's predicted function would change, should a differently chosen training set (but generated from the exact same process or system that we are trying to predict as the original) be used to train the model. A low variance is desired because essentially, we don't want to predict a very different function with a different training set that is generated from the same process. Model bias refers to the errors inherently introduced in the predicted function, as a result of the limitation as to what functional forms the specific model can learn. For example, linear models introduce bias when trying to approximate nonlinear functions because they can only learn linear functions. The ideal scenario for a good predictive model is to have both a low variance and a low bias. It is important for a predictive modeler to be aware of the fact that there is a bias-variance trade-off that arises from the choice of models. Models that are typically more complex because of the fact that they make fewer assumptions on the target function are prone to less bias, but higher variance than simpler but more restrictive models, such as linear models. This is because more complex models are able to approximate the training data more closely due to their flexibility, but as a result, they are more sensitive to changes in training data. This, of course, is also related to the problem of overfitting that complex models often exhibit.

We can actually see the effects of overfitting by first training some kNN models on our iris datasets. There are a number of packages that offer an implementation of the kNN algorithm, but we will use the knn3() function provided by the caret package with which we are familiar. To train a model using this function, all we have to do is provide it with a data frame that contains the numerical input features, a vector of output labels, and k, the number of nearest neighbors we want to use for the prediction:

> knn_model     <- knn3(iris_train, iris_train_labels, k = 5)
> knn_model_z   <- knn3(iris_train_z, iris_train_labels, k = 5)
> knn_model_pca <- knn3(iris_train_pca, iris_train_labels, k = 5)

To see the effect of different values of k, we will use the iris PCA model that is conveniently available in two dimensions for us to visualize and repeatedly train:

In the preceding plots, we have used different symbols to denote data points corresponding to different species. The lines shown in the plots correspond to the decision boundaries between the different species, which are the class labels of our output variable. Notice that using a low value of k, such as 1, captures local variations in the data very closely and, as a result, the decision boundaries are very irregular. A higher value of k uses many neighbors to create a prediction, resulting in a smoothing effect and smoother decision boundaries. Tuning k in kNN is an example of tuning a model parameter to balance the effect of overfitting.

We haven't mentioned any specific performance metrics in this section. There are different measures of model quality relevant to regression and classification, and we will address these after we wrap up our discussion on the predictive modeling process.

Repeating with different models and final model selection

During the first iteration of this process (and this is very much an iterative process!), we usually arrive at this stage, having trained and assessed a simple model. Simple models usually allow us to get to a quick and dirty solution with minimum effort, thus giving us an early sense of how far away we are from a model that will make predictions with reasonable accuracy. Simple models are also great at giving us a baseline level of performance against which we can benchmark the performance of future models. As modelers, we often acquire a preference toward one method over others, but it is important to remember that it is generally well worth the effort to try out different approaches to a problem and use the data to help us decide which one we should end up using.

Before picking the final model, it is worth considering whether it might be a good idea to use more than one model to solve our problem. In Chapter 7, Ensemble Methods, we spend an entire chapter studying techniques that involve many models working together to boost the predictive accuracy of the overall system.

Deploying the model

Once we have chosen the final model to use, we want to finalize its implementation so that the end users can use it reliably. Programmers refer to this process as deploying to production. This is where sound software engineering principles become extremely important. The following guidelines offer some useful advice:

• The model should be optimized to improve the speed at which it computes predictions. For example, this means ensuring that any features that are computed at runtime are done so efficiently.
• The model should be well documented. The final input features should be clearly defined, and the method and data used for training should be stored so that it can easily be retrained if changes need to be made. The original performance on the training and test set should also be stored as a reference for subsequent improvements.
• The model's performance should be monitored over time. This is important, not only as a means of verifying that the model works as intended, but also in order to catch any potential data shifts. If the process that is being modeled changes over time, it is likely that our model performance will degrade and this will signal the need for a new model to be trained.
• The software used to implement the model should be properly tested using standard unit and integration tests. Often, we will use a lot of existing R packages whose functions have already undergone testing, but the final deployment of a model may require us to write some additional code ourselves, such as for feature computation.
• The deployed model should be able to handle errors in the input. For example, if some of the input features are missing, it should be made very clear to the user why the model is unable to make a prediction through appropriate error messages. Errors and warnings should also be logged, especially if the model is deployed for continuous predictions in real-time settings.