CS109A Introduction to Data Science

Lab 3: Scikit-learn for Regression¶

Harvard University
Fall 2018
Instructors: Pavlos Protopapas and Kevin Rader
Authors: David Sondak, Will Claybaugh, Eleni Kaxiras

Table of Contents¶

1. Array creation and reshape
2. Some plotting
3. Simple linear regression
4. $k$-nearest neighbors

Learning Goals¶

Overall description and goal for the lab.

By the end of this lab, you should be able to:

• Understand array reshaping
• Review how to make plots
• Feel comfortable with simple linear regression
• Feel comfortable with $k$ nearest neighbors

This lab corresponds to lecture 4 and maps on to homework 2 (and beyond).

Simple Linear Regression¶

Linear regression and its many extensions are a workhorse of the statistics and data science community, both in application and as a reference point for other models. Most of the major concepts in machine learning can be and often are discussed in terms of various linear regression models. Thus, this section will introduce you to building and fitting linear regression models and some of the process behind it, so that you can 1) fit models to data you encounter 2) experiment with different kinds of linear regression and observe their effects 3) see some of the technology that makes regression models work.

Linear regression with a toy dataset¶

We first examine a toy problem, focusing our efforts on fitting a linear model to a small dataset with three observations. Each observation consists of one predictor $x_i$ and one response $y_i$ for $i = 1, 2, 3$,

To be very concrete, let's set the values of the predictors and responses.

There is no line of the form $\beta_0 + \beta_1 x = y$ that passes through all three observations, since the data are not collinear. Thus our aim is to find the line that best fits these observations in the least-squares sense, as discussed in lecture.

• Make two numpy arrays out of this data, x_train and y_train
• Check the dimentions of these arrays
• Try to reshape them into a different shape
• Make points into a very simple scatterplot
• Make a better scatterplot

Formulae¶

Linear regression is special among the models we study beuase it can be solved explicitly. While most other models (and even some advanced versions of linear regression) must be solved itteratively, linear regression has a formula where you can simply plug in the data.

For the single predictor case it is: \begin{align} \beta_1 &= \frac{\sum_{i=1}^n{(x_i-\bar{x})(y_i-\bar{y})}}{\sum_{i=1}^n{(x_i-\bar{x})^2}}\\ \beta_0 &= \bar{y} - \beta_1\bar{x}\ \end{align}

Where $\bar{y}$ and $\bar{x}$ are the mean of the y values and the mean of the x values, respectively.

From the re-aranged second equation we can see that the best-fit line passes through $(\bar{x},\bar{y})$, the center of mass of the data

From any of the first equations, we can see that the slope of the line has to do with whether or not an x value that is above/below the center of mass is typically paired with a y value that is likewise above/below, or typically paired with one that is opposite.

Building a model from scratch¶

In this part, we will solve the equations for simple linear regression and find the best fit solution to our toy problem.

The snippets of code below implement the linear regression equations on the observed predictors and responses, which we'll call the training data set. Let's walk through the code.

We have to reshape our arrrays to 2D. We will see later why.

• make an array with shape (2,3)
• reshape it to a size that you want

• Why the empty brackets? (The numerator and denominator are scalars, as expected.)

Turn the code from the above cells into a function called simple_linear_regression_fit, that inputs the training data and returns beta0 and beta1.

To do this, copy and paste the code from the above cells below and adjust the code as needed, so that the training data becomes the input and the betas become the output.

def simple_linear_regression_fit(x_train: np.ndarray, y_train: np.ndarray) -> np.ndarray:

    return

Check your function by calling it with the training data from above and printing out the beta values.

• Let's run this function and see the coefficients

• Do the values of beta0 and beta1 seem reasonable?
• Plot the training data using a scatter plot.
• Plot the best fit line with beta0 and beta1 together with the training data.

The values of beta0 and beta1 seem roughly reasonable. They capture the positive correlation. The line does appear to be trying to get as close as possible to all the points.

Building a model with statsmodels and sklearn¶

Now that we can concretely fit the training data from scratch, let's learn two python packages to do it all for us:

• statsmodels and
• scikit-learn (sklearn).

Our goal is to show how to implement simple linear regression with these packages. For an important sanity check, we compare the $\beta$ values from statsmodels and sklearn to the $\beta$ values that we found from above with our own implementation.

For the purposes of this lab, statsmodels and sklearn do the same thing. More generally though, statsmodels tends to be easier for inference [finding the values of the slope and intercept and dicussing uncertainty in those values], whereas sklearn has machine-learning algorithms and is better for prediction [guessing y values for a given x value]. (Note that both packages make the same guesses, it's just a question of which activity they provide more support for.

Note: statsmodels and sklearn are different packages! Unless we specify otherwise, you can use either one.

Below is the code for statsmodels. Statsmodels does not by default include the column of ones in the $X$ matrix, so we include it manually with sm.add_constant.

Besides the beta parameters, results_sm contains a ton of other potentially useful information.

Now let's turn our attention to the sklearn library.

We should feel pretty good about ourselves now, and we're ready to move on to a real problem!

The shape of things in scikit-learn¶

Before diving right in to a "real" problem, we really ought to discuss more of the details of sklearn. We do this now. Along the way, we'll import the real-world dataset.

Scikit-learn is the main python machine learning library. It consists of many learners which can learn models from data, as well as a lot of utility functions such as train_test_split. It can be used in python by the incantation import sklearn.

In scikit-learn, an estimator is a Python object that implements the methods fit(X, y) and predict(T)

Let's see the structure of scikit-learn needed to make these fits. .fit always takes two arguments:

estimator.fit(Xtrain, ytrain)

We will consider two estimators in this lab: LinearRegression and KNeighborsRegressor.

Critically, Xtrain must be in the form of an array of arrays (or a 2x2 array) with the inner arrays each corresponding to one sample, and whose elements correspond to the feature values for that sample (visuals coming in a moment).

ytrain on the other hand is a simple array of responses. These are continuous for regression problems.

Practice with sklearn¶

We begin by loading up the mtcars dataset and cleaning it up a little bit.

Next, let's split the dataset into a training set and test set.

Now we have training and test data. We still need to select a predictor and a response from this dataset. Keep in mind that we need to choose the predictor and response from both the training and test set. You will do this in the exercises below. However, we provide some starter code for you to get things going.

Notice the shape of y_train.

Another way to see the shape is to use the shape method.

This is not an "array of arrays". That's okay! Remember, sklearn requires an array of arrays only for the predictor array! You will have to pay close attention to this in the exercises later.

For now, let's discuss two ways out of this debacle. All we'll do is get y_train to be an array of arrays. This doesn't hurt anything because sklearn doesn't care too much about the shape of y_train.

First, let's reshape y_train to be an array of arrays using the reshape method. We want the first dimension of y_train to be size $25$ and the second dimension to be size $1$.

Notice that y_train.shape[0] gives the size of the first dimension.

There's an even easier way to get the correct shape right from the beginning.

Finally, there is a nice shortcut to reshaping an array. numpy can infer a dimension based on the other dimensions specified.

In this case, we said the second dimension should be size $1$. Since the requirement of the reshape() method is that the requested dimensions be compatible, numpy decides the the first dimension must be size $25$.

What would the .shape return if we did y_train.values.reshape(-1,5)?

Okay, enough of that. The whole reason we went through that whole process was to show you how to reshape your data into the correct format.

IMPORTANT: Remember that your response variable ytrain can be a vector but your predictor variable xtrain must be an array!

Simple linear regression with automobile data¶

We will now use sklearn to predict automobile mileage per gallon (mpg) and evaluate these predictions. We already loaded the data and split them into a training set and a test set.

We need to choose the variables that we think will be good predictors for the dependent variable mpg. 

• Pick one variable to use as a predictor for simple linear regression. Create a markdown cell below and discuss your reasons.
• Justify your choice with some visualizations.
• Is there a second variable you'd like to use? For example, we're not doing multiple linear regression here, but if we were, is there another variable you'd like to include if we were using two predictors?

• Use sklearn to fit the training data using simple linear regression.
• Use the model to make mpg predictions on the test set.
• Plot the data and the prediction.
• Print out the mean squared error for the training set and the test set and compare.

Hints:

• Use the following to perform the analysis:

  from sklearn.linear_model import LinearRegression
  from sklearn.model_selection import train_test_split
  from sklearn.metrics import mean_squared_error
  

$k$-nearest neighbors¶

Great, so we did a simple linear regression on the car data.

Now that you're familiar with sklearn, you're ready to do a KNN regression. Let's use $5$ nearest neighbors.

Lets vary the number of neighbors and see what we get.

Notice how the $1$-NN goes through every point on the training set but utterly fails elsewhere. Lets look at the scores on the training set.

Why do we get a perfect $R^2$ at k=1?

• Make the same plot as above on the test set.
• What is the best $k$?