Part A

Create a new variable for Age2. Then run a quadratic regression.

# make Age_sq variable
speed<-speed %>%
  mutate(Age_sq=Age^2)

# run quadratic regression
reg_quad<-lm(Amount ~ Age + Age_sq, data = speed) 
summary(reg_quad) 

## 
## Call:
## lm(formula = Amount ~ Age + Age_sq, data = speed)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -125.54  -44.96   -5.25   33.25  599.88 
## 
## Coefficients:
##               Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 146.746357   2.269451  64.662  < 2e-16 ***
## Age          -1.173833   0.125360  -9.364  < 2e-16 ***
## Age_sq        0.011357   0.001578   7.198 6.25e-13 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 56.08 on 31671 degrees of freedom
##   (36683 observations deleted due to missingness)
## Multiple R-squared:  0.005912,   Adjusted R-squared:  0.005849 
## F-statistic: 94.18 on 2 and 31671 DF,  p-value: < 2.2e-16

Part B

Try running the same regression using the alternate notation: lm(Y~X+I(X^2)). This method allows you to not have to create a new variable first. Do you get the same results?

reg_quad_alt<-lm(Amount ~ Age + I(Age^2), data = speed) 
summary(reg_quad_alt) 

## 
## Call:
## lm(formula = Amount ~ Age + I(Age^2), data = speed)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -125.54  -44.96   -5.25   33.25  599.88 
## 
## Coefficients:
##               Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 146.746357   2.269451  64.662  < 2e-16 ***
## Age          -1.173833   0.125360  -9.364  < 2e-16 ***
## I(Age^2)      0.011357   0.001578   7.198 6.25e-13 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 56.08 on 31671 degrees of freedom
##   (36683 observations deleted due to missingness)
## Multiple R-squared:  0.005912,   Adjusted R-squared:  0.005849 
## F-statistic: 94.18 on 2 and 31671 DF,  p-value: < 2.2e-16

Yes, this is the same output.

Part C

Write out the estimated regression equation.

^Amounti=146.75−1.17Agei+0.01Age2i

Part D

Is this model an improvement from the linear model?¹

Yes, a slight improvement. R2 went from 0.004 on the linear model to 0.006 on the quadratic model.

Part E

Write an equation for the marginal effect of Age on Amount.

The marginal effect is measured by the first derivative of the regression equation with respect to Age. But you can just remember the resulting formula and plug in the parameters:

dYdX=β1+2β2XdAmountdAge=−1.17+2(0.01)Age=−1.17+0.02Age

# PUT CODE HERE

Part F

Predict the marginal effect on Amount of being one year older when you are 18. How about when you are 40?

For 18 year olds:

dAmountdAge=−1.17+0.02(18)=−1.17+0.36=−$0.81

For 40 year olds:

dAmountdAge=−1.17+0.02(40)=−1.17+0.80=−$0.37

# Let's do this in R:

# we need broom
library(broom)
tidy_reg_quad<-tidy(reg_quad)

tidy_reg_quad

# save beta 1 
quad_beta_1<-tidy_reg_quad %>%
  filter(term=="Age") %>%
  pull(estimate)

# save beta 2
quad_beta_2<-tidy_reg_quad %>%
  filter(term=="Age_sq") %>%
  pull(estimate)

# create function to estimate marginal effects
marginal_effect<-function(x){
  return(quad_beta_1+2*quad_beta_2*x)
}

# run the function on the 18-year-old and the 40-year-old
marginal_effect(c(18,40))

## [1] -0.7649738 -0.2652574

# close enough, we had some rounding error

Part G

Our quadratic function is a U-shape. According to the model, at what age is the amount of the fine minimized?

We can set the derivative equal to 0, or you can just remember the formula and plug in the parameters:

dYdX=β1+2β2X0=β1+2β2X−β1=2β2X−12β1β2=Age∗−12−1.170.01=Age∗−12117≈Age∗58.5≈Age∗

# Let's do this in R:

-0.5*(quad_beta_1/quad_beta_2)

## [1] 51.67795

# again, some rounding error

Part H

Create a scatterplot between Amount and Age and add a a layer of a linear regression (as always), and an additional layer of your predicted quadratic regression curve. The regression curve, just like any regression line, is a geom_smooth() layer on top of the geom_point() layer. We will need to customize geom_smooth() to geom_smooth(method="lm", formula="y~x+I(x^2) (copy/paste this verbatim)! This is the same as a regression line (method="lm"), but we are modifying the formula to a polynomial of degree 2 (quadratic): y=a+bx+cx2.

ggplot(data=speed)+
  aes(x = Age,
      y = Amount)+
  geom_point()+
  geom_smooth(method="lm",
              formula="y~poly(x,2)",
              color="red")

## Warning: Removed 36683 rows containing non-finite values (stat_smooth).

## Warning: Removed 36683 rows containing missing values (geom_point).

Part I

It’s quite hard to see the quadratic curve with all those data points. Redo another plot and this time, only keep the quadratic stat_smooth() layer and leave out the geom_point() layer. This will only plot the regression curve.

ggplot(data=speed)+
  aes(x = Age,
      y = Amount)+
  geom_smooth(method="lm",
              formula="y~poly(x,2)",
              color="red")

## Warning: Removed 36683 rows containing non-finite values (stat_smooth).

Part A

Creating new variables as necessary, run a linear-log model of Amount on MPHover. Write down the estimated regression equation, and interpret the coefficient on MPHover (^β1). Make a scatterplot with the regression line.²

# create log of MPHover
speed<-speed %>%
  mutate(log_mph=log(MPHover))

# Run linear-log regression
linear_log_reg<-lm(Amount ~ log_mph, data = speed)
summary(linear_log_reg)

## 
## Call:
## lm(formula = Amount ~ log_mph, data = speed)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -190.61  -16.44    8.56   20.52  425.33 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -200.0975     1.9401  -103.1   <2e-16 ***
## log_mph      115.7544     0.6922   167.2   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 40.99 on 31672 degrees of freedom
##   (36683 observations deleted due to missingness)
## Multiple R-squared:  0.4689, Adjusted R-squared:  0.4689 
## F-statistic: 2.796e+04 on 1 and 31672 DF,  p-value: < 2.2e-16

# note we could have done this without creating the variable
# just take log() inside the regression:
linear_log_reg_alt<-lm(Amount ~ log(MPHover), data = speed)
summary(linear_log_reg)

## 
## Call:
## lm(formula = Amount ~ log_mph, data = speed)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -190.61  -16.44    8.56   20.52  425.33 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -200.0975     1.9401  -103.1   <2e-16 ***
## log_mph      115.7544     0.6922   167.2   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 40.99 on 31672 degrees of freedom
##   (36683 observations deleted due to missingness)
## Multiple R-squared:  0.4689, Adjusted R-squared:  0.4689 
## F-statistic: 2.796e+04 on 1 and 31672 DF,  p-value: < 2.2e-16

^Amounti=−200.10+115.75ln(MPHoveri)

A 1% increase in speed (over the speed limit) increases the fine by 115.75100≈$1.16

Part B

Creating new variables as necessary, run a log-linear model of Amount on MPHover. Write down the estimated regression equation, and interpret the coefficient on MPHover (^β1). Make a scatterplot with the regression line.³

# create log of Amount
speed<-speed %>%
  mutate(log_Amount=log(Amount))

# Run log-linear regression
log_linear_reg<-lm(log_Amount ~ MPHover, data = speed)
summary(log_linear_reg)

## 
## Call:
## lm(formula = log_Amount ~ MPHover, data = speed)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -4.2353 -0.1895  0.0929  0.2733  1.2970 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 3.8790452  0.0058716   660.6   <2e-16 ***
## MPHover     0.0484962  0.0003256   148.9   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.3355 on 31672 degrees of freedom
##   (36683 observations deleted due to missingness)
## Multiple R-squared:  0.4119, Adjusted R-squared:  0.4119 
## F-statistic: 2.219e+04 on 1 and 31672 DF,  p-value: < 2.2e-16

# again we could have done this without creating the variable
# just take log() inside the regression:
log_linear_reg_alt<-lm(log(Amount) ~ MPHover, data = speed)
summary(log_linear_reg_alt)

## 
## Call:
## lm(formula = log(Amount) ~ MPHover, data = speed)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -4.2353 -0.1895  0.0929  0.2733  1.2970 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 3.8790452  0.0058716   660.6   <2e-16 ***
## MPHover     0.0484962  0.0003256   148.9   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.3355 on 31672 degrees of freedom
##   (36683 observations deleted due to missingness)
## Multiple R-squared:  0.4119, Adjusted R-squared:  0.4119 
## F-statistic: 2.219e+04 on 1 and 31672 DF,  p-value: < 2.2e-16

^ln(Amounti)=3.87+0.05MPHoveri

For every 1 MPH in speed (over the speed limit), expected fine increases by 0.05×100%=5%

Part C

Creating new variables as necessary, run a log-log model of Amount on MPHover. Write down the estimated regression equation, and interpret the coefficient on MPHover (^β1). Make a scatterplot with the regression line.⁴

# Run log-log regression
log_log_reg<-lm(log_Amount ~ log_mph, data = speed)
summary(log_log_reg)

## 
## Call:
## lm(formula = log_Amount ~ log_mph, data = speed)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -4.1417 -0.1727  0.1034  0.2319  2.3669 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  2.30860    0.01564   147.6   <2e-16 ***
## log_mph      0.86196    0.00558   154.5   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.3304 on 31672 degrees of freedom
##   (36683 observations deleted due to missingness)
## Multiple R-squared:  0.4297, Adjusted R-squared:  0.4297 
## F-statistic: 2.386e+04 on 1 and 31672 DF,  p-value: < 2.2e-16

# again we could have done this just taking log()s inside the regression:
log_log_reg_alt<-lm(log(Amount) ~ log(MPHover), data = speed)
summary(log_log_reg_alt)

## 
## Call:
## lm(formula = log(Amount) ~ log(MPHover), data = speed)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -4.1417 -0.1727  0.1034  0.2319  2.3669 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept)   2.30860    0.01564   147.6   <2e-16 ***
## log(MPHover)  0.86196    0.00558   154.5   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.3304 on 31672 degrees of freedom
##   (36683 observations deleted due to missingness)
## Multiple R-squared:  0.4297, Adjusted R-squared:  0.4297 
## F-statistic: 2.386e+04 on 1 and 31672 DF,  p-value: < 2.2e-16

^ln(Amounti)=2.31+0.86ln(MPHoveri)

For every 1% increase in speed (over the speed limit), expected fine increases by 0.86%.

Part D

Which of the three log models has the best fit?⁵

We can compare the R2’s of the three models or compare scatterplots with the regression lines. I will make a table of the three regressions with huxreg for easy comparison of fit:

library(huxtable)

## 
## Attaching package: 'huxtable'

## The following object is masked from 'package:dplyr':
## 
##     add_rownames

## The following object is masked from 'package:purrr':
## 
##     every

## The following object is masked from 'package:ggplot2':
## 
##     theme_grey

huxreg("Linear-Log" = linear_log_reg,
       "Log-Linear" = log_linear_reg,
       "Log-Log" = log_log_reg,
       coefs = c("Constant" = "(Intercept)",
                 "MPH Over" = "MPHover",
                 "Log(MPH Over)" = "log_mph"),
       statistics = c("N" = "nobs",
                      "R-Squared" = "r.squared",
                      "SER" = "sigma"),
       number_format = 2)

It appears the linear-log model has the best fit with the highest R2 out of the three, but not by very much.

If we wanted to compare scatterplots:

ggplot(data = speed)+
  aes(x = log_mph,
      y = Amount)+
  geom_point()+
  geom_smooth(method="lm", color="red")+
  labs(title="Linear-Log Model")

## Warning: Removed 36683 rows containing non-finite values (stat_smooth).

## Warning: Removed 36683 rows containing missing values (geom_point).

ggplot(data = speed)+
  aes(x = MPHover,
      y = log_Amount)+
  geom_point()+
  geom_smooth(method="lm", color="red")+
  labs(title="Log-Linear Model")

## Warning: Removed 36683 rows containing non-finite values (stat_smooth).

## Warning: Removed 36683 rows containing missing values (geom_point).

ggplot(data = speed)+
  aes(x = log_mph,
      y = log_Amount)+
  geom_point()+
  geom_smooth(method="lm", color="red")+
  labs(title="Log-Log Model")

## Warning: Removed 36683 rows containing non-finite values (stat_smooth).

## Warning: Removed 36683 rows containing missing values (geom_point).

Part A

Create new standardized variables for Amount, Age, and MPHover.⁶

# make standardized variables
speed<-speed %>%
  mutate(std_Amount = scale(Amount),
         std_Age = scale(Age),
         std_MPH = scale(MPHover))

Part B

Run a regression of standardized Amount on standardized Age and MPHover. Interpret ^β1 and ^β2 Which variable has a bigger effect on Amount?

# make standardized variables
std_reg<-lm(std_Amount ~ std_Age + std_MPH, data = speed)
summary(std_reg)

## 
## Call:
## lm(formula = std_Amount ~ std_Age + std_MPH, data = speed)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -5.4896 -0.3517  0.1366  0.4579  4.0234 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -0.234516   0.004205 -55.768   <2e-16 ***
## std_Age      0.005986   0.004220   1.418    0.156    
## std_MPH      0.622823   0.003496 178.130   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.7053 on 31671 degrees of freedom
##   (36683 observations deleted due to missingness)
## Multiple R-squared:  0.5026, Adjusted R-squared:  0.5026 
## F-statistic: 1.6e+04 on 2 and 31671 DF,  p-value: < 2.2e-16