Discussion_2

Author

Angelo Marinaro

Packages List

library(AER)
Loading required package: car
Loading required package: carData
Loading required package: lmtest
Loading required package: zoo

Attaching package: 'zoo'
The following objects are masked from 'package:base':

    as.Date, as.Date.numeric
Loading required package: sandwich
Loading required package: survival
library(dplyr)

Attaching package: 'dplyr'
The following object is masked from 'package:car':

    recode
The following objects are masked from 'package:stats':

    filter, lag
The following objects are masked from 'package:base':

    intersect, setdiff, setequal, union
library(tidyr)
library(lubridate)

Attaching package: 'lubridate'
The following objects are masked from 'package:base':

    date, intersect, setdiff, union
library(haven)
library(ggplot2)
library(zoo)
library(forecast)
Registered S3 method overwritten by 'quantmod':
  method            from
  as.zoo.data.frame zoo 
library(stargazer)

Please cite as: 
 Hlavac, Marek (2022). stargazer: Well-Formatted Regression and Summary Statistics Tables.
 R package version 5.2.3. https://CRAN.R-project.org/package=stargazer 
library(glmnet)
Loading required package: Matrix

Attaching package: 'Matrix'
The following objects are masked from 'package:tidyr':

    expand, pack, unpack
Loaded glmnet 4.1-8
library(vars)
Loading required package: MASS

Attaching package: 'MASS'
The following object is masked from 'package:dplyr':

    select
Loading required package: strucchange
Loading required package: urca
library(corrplot)
corrplot 0.92 loaded
library(car)

I. Choosing Data set

#library(help = "datasets")
#data(package = "AER")
#?mtcars
dt <- mtcars
head(dt)
                   mpg cyl disp  hp drat    wt  qsec vs am gear carb
Mazda RX4         21.0   6  160 110 3.90 2.620 16.46  0  1    4    4
Mazda RX4 Wag     21.0   6  160 110 3.90 2.875 17.02  0  1    4    4
Datsun 710        22.8   4  108  93 3.85 2.320 18.61  1  1    4    1
Hornet 4 Drive    21.4   6  258 110 3.08 3.215 19.44  1  0    3    1
Hornet Sportabout 18.7   8  360 175 3.15 3.440 17.02  0  0    3    2
Valiant           18.1   6  225 105 2.76 3.460 20.22  1  0    3    1

This data set is cross sectional since it does not contain a time variable. The data contains fuel consumption and other aspects for 32 automobiles for the year 1974.

The variables are:

  • mpg : miles per gallon

  • cyl: number of cylinders of the engine

  • disp: displacement in cubic inches

  • hp: gross horsepower

  • drat: rear axle ratio

  • wt: weight (1000 lbs)

  • qsec: 1/4 mile time

  • vs: engine shape

  • am: transmission

  • gear: number of forward gears

  • carb: number of carburetors

II. Choosing variables

cor1 <- cor(dt)
corrplot(cor1, type = "upper", order = "original", tl.col = "black" ,
         tl.srt = 90, tl.cex = .8)

The variables that seem to have a stronger relationship with miles per gallon are cyl, disp, hp, drat, and wt. Trying to account for multicollinearity, I am choosing cyl and wt as my variables for my regression. They both show a negative correlation to miles per gallon and show a positive correlation between themselves.

hist(dt$mpg, 
     xlab="Miles per Gallon",
     main = "Miles per Gallon")

hist(dt$cyl,
     xlab="Cylinders",
     main = "Cylinders")

hist(dt$wt,
     main = "Weight", 
     xlab="Weight (1000 lbs)")

mtcars_sub<- mtcars[, c("cyl", "wt", "mpg")]
pairs(mtcars_sub, main = "Scatterplot Matrix of mtcars")

\[ \begin{align*} Miles\ per \ gallon_i = \beta_0 + \beta_1 Cylinder + \beta_2Weight + \epsilon_i \end{align*} \]

III. Running OLS using lm

reg1 <- lm(data = dt, 
           formula = mpg ~ cyl + wt)

summary(reg1)

Call:
lm(formula = mpg ~ cyl + wt, data = dt)

Residuals:
    Min      1Q  Median      3Q     Max 
-4.2893 -1.5512 -0.4684  1.5743  6.1004 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)  39.6863     1.7150  23.141  < 2e-16 ***
cyl          -1.5078     0.4147  -3.636 0.001064 ** 
wt           -3.1910     0.7569  -4.216 0.000222 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 2.568 on 29 degrees of freedom
Multiple R-squared:  0.8302,    Adjusted R-squared:  0.8185 
F-statistic: 70.91 on 2 and 29 DF,  p-value: 6.809e-12
coefficients <- coef(reg1)
print(coefficients)
(Intercept)         cyl          wt 
  39.686261   -1.507795   -3.190972

Interpretation:

Since I did not perform any transformations to the data, the model is a level-level, meaning that the interpretation for the variables are:

  • For a unit increase in the number of cylinders, miles per gallon reduce by 1.51 units

  • For a unit increase in weight of the car, miles per gallon decrease by 3.19 units

IV. Running OLS matrix algebra

y <- as.vector(mtcars_sub$mpg)

X <- as.matrix(mtcars_sub[-ncol(mtcars_sub)])

int <- rep(x = 1, 
           times =  length(y))

X <- cbind(int, X)
coeff <- round(solve(t(X) %*% X) %*% t(X) %*% y,
               4)
 
coeff
       [,1]
int 39.6863
cyl -1.5078
wt  -3.1910

Interpretation:

Since I did not perform any transformations to the data, the model is a level-level, meaning that the interpretation for the variables are:

  • For a unit increase in the number of cylinders, miles per gallon reduce by 1.51 units

  • For a unit increase in weight of the car, miles per gallon decrease by 3.19 units

Both outputs from the lm() function and by doing matrix algebra yield the same results.

V. Bonus

N   <- nrow(X)
k   <- ncol(X) - 1

beta_hat    <-    solve( t(X) %*% X )  %*%  t(X) %*% y

sigma_sq  <- 
  residual_variance <- (N-k-1)^-1 * sum((y - X %*% beta_hat)^2)
sigma_sq
[1] 6.592137
residual_std_error  <- sqrt(residual_variance)

var_betaHat       <- sigma_sq * solve(t(X) %*% X)
round(var_betaHat, 
      digits = 1)
     int  cyl   wt
int  2.9 -0.3 -0.3
cyl -0.3  0.2 -0.2
wt  -0.3 -0.2  0.6
coeff_std_errors  <- sqrt(diag(var_betaHat))
coeff_std_errors
      int       cyl        wt 
1.7149840 0.4146883 0.7569065