# Section 0: Getting set up

Load the following packages:

```{r, echo = FALSE, eval = TRUE}
# You probably already have these packages installed, so let's just load them
library(tidyverse)
library(modelr)
library(knitr)
library(broom)

options(scipen = 999) # disable scientific notation

# For labs, we want to see all our code
knitr::opts_chunk$set(echo = TRUE)
```

# Section 1: Interactions in OLS in `R`

In this section we will focus on estimating an *interaction* model in `R`. Recall from lecture that an interaction model is appropriate when you think the effect of some variable $x$ on your dependent variable $y$ depends on a third variable $z$. For example, Graff-Zivin and Neidell (2012) find that the productivity of female agricultural workers is less sensitive to air pollution than that of men's. This suggests that *gender* influences the relationship between *productivity* and *pollution*.

To learn how to estimate interaction models, we will use an example from the automobile industry (from last week) which has data on fuel efficiency and automobile characteristics for cars of two vintages: 1999 cars, and 2008 cars. Last week we found that vintage influenced fuel efficiency: more recent cars have higher miles per gallon. We also found that our estimates were biased if we did not also control for engine size, as higher engine sizes lead to lower fuel efficiency, and they are also correlated with vintage.

## Parallel slopes model:

The model we estimated last week was:

$$hwy_i =\beta_{0}+\beta_{1} \cdot displ_i +\beta_{2} \cdot \text vintage_i+\varepsilon_i$$

The results of this "parallel slopes" model from last week looked like this:

```{r, echo=TRUE}
mpg <- mpg %>% mutate(year = as.factor(year)) # Recall, we do this to ensure our year variable is treated as a categorical variable

mod <- lm(hwy ~ displ + year, data = mpg)

mpg %>% 
  ggplot(aes(x = displ, y = hwy, color = year)) +
  geom_point() +
  geom_line(data = augment(mod), aes(y = .fitted, color = year)) + 
  labs(x = "Engine Displacement, in litres",
       y = "Highway Miles per gallon") +
  scale_colour_discrete("Year")
```

These slopes being **parallel** comes from our regression equation having just one coefficient on engine size; regardless of your vintage, we estimate the same impact of engine size on fuel economy.

This may not be accurate. If technology is improving over time, the effect of engine size on fuel economy could decline over time -- having a larger engine lowers fuel economy, but maybe we are getting better at not sacrificing miles per gallon as we increase engine size. In this subsection, we want to see if this true or or not. We will use the `mpg` dataset and is pre-loaded in `R`.

## Interaction model

We hypothesize that the relationship between miles per gallon and engine size may not actually be the same across the two vintages. That is, we want a model in which those "parallel slopes" above are allowed to differ by vintage.

We use an interaction model to achieve this: $$hwy_i=\beta_{0}+\beta_{1} \cdot displ_i + \beta_{2} \cdot vintage_i + \beta_{3} \cdot displ_i \cdot vintage_i + \varepsilon_i$$ The `lm` package makes estimating interactions easy. Just use a colon (`:`) between the two variables you want to interact; the colon tells `lm` to multiply the two variables by one another.

Complete the following:

1.  Using the model specified above, use `lm()` to estimate all three coefficients using the `mpg` data. Make sure you are treating `year` as a categorical variable! Use `summary(lm())` to print the regression results.

```{r}

```

2.  What does the intercept tell you, in words?

3.  What does the coefficient on `displ` tell you, in words?

4.  What does the coefficient estimate on the 2008 vintage indicator variable (`year2008`) tell you, in words?

5.  What does the coefficient on `displ:year2008` tell you, in words?[^1]

[^1]: Hint: This is tricky. Start by writing down the relationship between miles per gallon and engine size when the car vintage is 2008. Then write the same thing for car vintage 1999. Can you see what role this coefficient plays?

# Section 2: Plotting an interaction model

Interaction models are difficult to interpret. Writing out the algebra of what each independent variable and its associated coefficient mean is really valuable, but visualizing the data and resulting model can also help. Here, let's remake our scatter plot above, but overlay our regression results from the interaction model.

**Intuition check before you dive in:** What do you expect to see in terms of intercepts and slopes of your vintage-specific regression lines?

We will use `ggplot2` with `augment()` to obtain predictions, in combination with `geom_line()`.

```{r}

```

# Section 3: Adjusted $R^2$

We saw in the previous exercise that an interaction model resulted in different slopes for each group. However, is this more complex model actually "better"? We have explored using $R^2$ as a measure of model fit. It is a common one and you should be comfortable computing and interpreting it. Maybe it can help us answer the question of which model to use in this case.

However, as we discussed in lecture, $R^2$ *mechanically* increases as you add more independent variables to your regression. Adjusted $R^2$ is designed to "penalize" your measure of model fit for those additional independent variables.

These two measures are: $$R^{2}=1 - \frac{\sum_i e_i^2}{\sum_i (y_i - \bar y)^2}$$ $$ \overline{R}^2 = 1 - \dfrac{\sum_i e_i^2/(n-k-1)}{\sum_i \left( y_i - \overline{y} \right)^2/(n-1)} $$

In this case, the interaction model added a new variable to our regression, so we should be comparing adjusted $R^2$ results from a model without the interaction to a model with it. You can see the adjusted $R^2$ in `summary(lm())`, but you can also access it by saving your `lm` object and calling `summary(mod)$adj.r.squared`.

```{r}
# no interaction

# with the interaction

```

Since the adjusted $R^2$ value declined slightly, we conclude that the interaction model did *not* meaningfully improve model fit.

Note that if you compare $R^2$ across models, without accounting for the additional independent variables, you get a slightly different picture:

```{r}
# no interaction

# with the interaction

```

For multiple regression, you should generally rely on adjusted $R^2$ and not $R^2$.
