---
title: "Multiple regression principles"
subtitle: "Notes and in-class exercises"
format: 
  html:
    embed-resources: true
    toc: true
---


# Notes

## Learning goals

Working with multiple predictors in our plots and models can get complicated!
There are no *recipes* for this process.

BUT there are some guiding *principles* that assist in long-term retention, deeper understanding, and the ability to generalize our tools in new settings.

By the end of this lesson, you should be familiar with some general *principles* for... 

- incorporating additional quantitative or categorical predictors in a visualization
- how additional quantitative or categorical predictors impact the physical representation of a model
- interpreting quantitative or categorical coefficients in a multiple regression model

## Readings and videos

Please watch the following video **before** class.

- [Interpreting multivariate models](https://www.youtube.com/watch?v=DnXwu1OWgMM) ([slides](https://drive.google.com/file/d/1fiJL1IbReg6RJmjxRideYzYTyQtWR3Or/view?usp=sharing))

# Exercises

Let's revisit the bikeshare data:

```{r warning = FALSE, message = FALSE}
# Load packages & import data
library(readr)
library(ggplot2)
library(dplyr)

bikes <- read_csv("https://mac-stat.github.io/data/bikeshare.csv") %>% 
  rename(rides = riders_registered)
```

Our goal is to understand how / why registered ridership from day to day.

To this end, we'll build various **multiple linear regression** models of `rides` by different combinations of the possible predictors.

```{r}
# Check out the data
head(bikes)
```



## Exercise 1: Review visualization

Let's build a model of `rides` by `windspeed` (quantitative) and `weekend` status (categorical).

a. Write a model statement for this regression model.

> E[Y | X] = ...

b. Plot & describe, in words, the relationship between these 3 variables.


```{r}
# Plot of rides vs windspeed & weekend
# HINT: Start with a plot of rides vs windspeed, then add an aesthetic for weekend!

```

> Your description here


## Exercise 2: Review model

Let's build the model. Run the following code:

```{r}
bike_model_1 <- lm(rides ~ windspeed + weekend, data = bikes)
coef(summary(bike_model_1))
```

The model formula with our coefficient estimates filled in is therefore:

E[rides | windspeed, weekendTRUE] = 4738.38 - 63.97 * windspeed - 925.16 * weekendTRUE

This model formula is represented by 2 lines, one corresponding to weekends and the other to weekdays. Simplify the model formula above for weekdays and weekends:       
    
weekdays: rides = ___ - ___ windspeed   

weekends: rides = ___ - ___ windspeed

> Pause. Did you fill in the model formulas above? Go back and look at activity 08 if you need an example!

## Exercise 3: Review coefficient interpretation


a. The intercept coefficient, 4738.38, represents the *intercept* of the sub-model for weekdays, the reference category. What's its *contextual* interpretation?
    
> Your response here

b. The `windspeed` coefficient, -63.97, represents the *shared slope* of the weekend and weekday sub-models. What's its *contextual* interpretation?
    
> Your response here

c. The `weekendTRUE` coefficient, -925.16, represents the *change in intercept* for the weekend vs weekday sub-model. What's its *contextual* interpretation?

> Your response here


## Exercise 4: 2 categorical predictors -- visualization

Thus far, we've explored a couple examples of multiple regression models that have 2 predictors, 1 quantitative and 1 categorical.

So what happens when *both* predictors are categorical?!

To this end, let's model `rides` by `weekend` status and `season`.

The below code plots `rides` vs `season`.

Modify this code to *also* include information about `weekend`.

HINT: Remember the visualization *principle* that additional categorical predictors require some sort of grouping mechanism / mechanism that distinguishes between the 2 groups.

```{r}
# rides vs season
bikes %>% 
  ggplot(aes(y = rides, x = season)) + 
  geom_boxplot()

# rides vs season AND weekend
bikes %>%
  ggplot(aes(y = rides, x = season, ___ = ___)) +
  geom_boxplot()
```

> Pause. Did your modified code run? Can you see all the variables represented in the second plot above?

## Exercise 5: follow-up

a. Describe (in words) the relationship of ridership with season & weekend status.

> Your response here

b. A model of `rides` by `season` alone would be represented by only 4 expected outcomes, 1 for each season. Considering this and the plot above, how do you *anticipate* a model of `rides` by `season` and `weekend` status will be represented?  


    - 2 lines, 1 for each weekend status
    - 8 lines, 1 for each possible combination of season & weekend
    - 2 expected outcomes, 1 for each weekend status
    - 8 expected outcomes, 1 for each possible combination of season & weekend
    
> Copy/paste which one you think is correct here!



## Exercise 6: 2 categorical predictors -- build the model

Let's build the multiple regression model of `rides` vs `season` and `weekend`:
    
```{r}
bike_model_2 <- lm(rides ~ weekend + season, bikes)
coef(summary(bike_model_2))
```

Thus the model formula with coefficient estimates filled in is given by:

E[rides | weekend, season] = 4260.45 - 912.33 weekendTRUE - 116.38 seasonspring + 438.44 seasonsummer - 1719.06 seasonwinter

a. Use this model to predict the ridership on the following days:    
    
```{r}
# a fall weekday
4260.45 - 912.33*___ - 116.38*___  + 438.44*___ - 1719.06*___

# a winter weekday    
4260.45 - 912.33*___ - 116.38*___  + 438.44*___ - 1719.06*___

# a fall weekend day        
4260.45 - 912.33*___ - 116.38*___  + 438.44*___ - 1719.06*___

# a winter weekend day
4260.45 - 912.33*___ - 116.38*___  + 438.44*___ - 1719.06*___
```

> Pause. Did you fill in the code above? If you want to check that you did it correctly, look at the solutions on the course website, or call over Taylor or a preceptor to check your answer!

b. We only made 4 predictions here. How many *possible* predictions does this model produce? Is this consistent with your intuition in the previous exercise?

> Your response here




## Exercise 7: 2 categorical predictors -- interpret the model

Use your above predictions and visualization to fill in the below interpretations of the model coefficients.

Hint: What is the consequence of plugging in 0 or 1 for the different `weekend` and `season` categories?    

a. Interpreting 4260: On average, we expect there to be 4260 riders on (weekdays/weekends) during the (fall/spring/summer/winter).    

b. Interpreting -912: On average, *in any season*, we expect there to be 912 (more/fewer) riders on weekends than on ___.

An alternative interpretation: On average, we expect there to be 912 (more/fewer) riders on weekends than on ___, adjusting for season.

c. Interpreting -1719: On average, *on both weekdays and weekends*, we expect there to be 1719 (more/fewer) riders in winter than in ___.

An alternative interpretation: On average, we expect there to be 1719 (more/fewer) riders in winter than in ___, controlling for weekday status.


> Pause. For parts a through c, did you select the correct word, and fill in everything appropriately?





## Exercise 8: 2 quantitative predictors -- visualization   

Next, consider the relationship between `rides` and 2 *quantitative* predictors: `windspeed` and `temp_feel`.

Check out the plot of this relationship below.

This reflects the visualization *principle* that quantitative variables require some sort of numerical scaling mechanism -- `rides` and `windspeed` get numerical axes, and `temp_feel` gets a color scale.

![](https://mac-stat.github.io/images/155/bikes_multivar_3.png)

Modify the code below to recreate this plot.

```{r}
bikes %>%
  ggplot(aes(y = rides, x = windspeed, ___ = ___)) +
  geom_point()
```

> Pause. Does your recreated plot look the same as the image above / on the course website?


## Exercise 9: follow-up

Describe (in words) the relationship of ridership with windspeed & temperature.


> Your response here


## Exercise 10: 2 quantitative predictors -- modeling

Let's build the multiple regression model of `rides` vs `windspeed` and `temp_feel`:
    
```{r}
bike_model_3 <- lm(rides ~ windspeed + temp_feel, data = bikes)
coef(summary(bike_model_3))
```

Thus the model formula with coefficient estimates filled in is given by,

E[rides | windspeed, temp_feel] = -24.06 - 36.54 windspeed + 55.52 temp_feel


a. Interpret the intercept coefficient, -24.06, in context.

> Your response here

b. Interpret the `windspeed` coefficient, -36.54, in context. 

> Your response here

c. Interpret the `temp_feel` coefficient, 55.52, in context. 

> Your response here

    
    


## Exercise 11: Which is "best"?

We've now observed 3 different models of ridership, each having 2 predictors.

The R-squared values of these models, along with those of the simple linear regression models with each predictor alone, are summarized below.
    
model                     predictors                R-squared
------------------------- ------------------------- ------------
`bike_model_1`            `windspeed` & `weekend`   0.119
`bike_model_2`            `weekend` & `season`      0.349
`bike_model_3`            `windspeed` & `temp_feel` 0.310
`bike_model_4`            `windspeed`               0.047
`bike_model_5`            `temp_feel`               0.296
`bike_model_6`            `weekend`                 0.074
`bike_model_7`            `season`                  0.279


a. Which model does the best job of explaining the variability in ridership from day to day?

> Your response here

b. If you could only pick one predictor, which would it be?

> Your response here

c. What happens to R-squared when we add a second predictor to our model, and why does this make sense? For example, how does the R-squared for model 1 (with both windspeed and weekend) compare to those of model 4 (only windspeed) and model 6 (only weekend)?

> Your response here

d. Are 2 predictors always better than 1? Provide evidence and explain why this makes sense.

> Your response here



## Exercise 12: Principles of interpretation

These exercises have revealed some **principles** behind interpreting model coefficients, summarized below.

Review and confirm that these make sense.

---

**Principles of interpretation**

Consider a multiple linear regression model:

$$E[y | x_1, x_2, ..., x_p] = \beta_0 + \beta_1 x_1 + \beta_2 x_2 + ... + \beta_p x_p$$

We can interpret the coefficients as follows:    

- $\beta_0$ ("beta 0") is the y-intercept. It describes the average value of $y$ when $x_1, x_2,..., x_k$ are all 0, ie. when all quantitative predictors are set to 0 and the categorical predictors are set to their reference levels.    

- $\beta_i$ ("beta i") is the coefficient of $x_i$.    

    - If $x_i$ is quantitative, $\beta_i$ describes the average change in $y$ associated with a 1-unit increase in $x_i$ while at a fixed set of the other $x$.  
    
    - If $x_i$ represents a category of a categorical variable, $\beta_i$ describes the average difference in $y$ between this category and the reference category, while at a fixed set of the other $x$.    

---


> Pause. Does the above all make sense? Write down any remaining questions you may have here, that you could ask in office hours.



# Extra practice

The following exercises provide extra practice.
If you don't get to these during class, you're encouraged to try them outside class.


## Exercise 13: Practice 1

Consider the relationship of `rides` vs `weekend` and `weather_cat`.    

a. Construct a visualization of this relationship.   
b. Construct a model of this relationship.    
c. Interpret the first 3 model coefficients. 



## Exercise 14: Practice 2

Consider the relationship of `rides` vs `temp_feel` and `humidity`.    

a. Construct a visualization of this relationship.    
b. Construct a model of this relationship.    
c. Interpret the first 3 model coefficients.    
    



## Exercise 15: Practice 3

Consider the relationship of `rides` vs `temp_feel` and `weather_cat`.    

a. Construct a visualization of this relationship.    
b. Construct a model of this relationship.    
c. Interpret the first 3 model coefficients.    
    


## Exercise 16: CHALLENGE

We've explored models with 2 predictors.
What about 3 predictors?!
Consider the relationship of `rides` vs `temp_feel`, `humidity`, AND `weekend`.

a. Construct a visualization of this relationship.    
b. Construct a model of this relationship.    
c. Interpret each model coefficient.    
 
          


