All installed packages were loaded to be used in the R Program.
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=stargazerlibrary(lmtest)## Loading required package: zoo##
## Attaching package: 'zoo'## The following objects are masked from 'package:base':
##
## as.Date, as.Date.numericlibrary(tidyverse)## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr 1.1.1 ✔ readr 2.1.4
## ✔ forcats 1.0.0 ✔ stringr 1.5.0
## ✔ ggplot2 3.4.2 ✔ tibble 3.2.1
## ✔ lubridate 1.9.2 ✔ tidyr 1.3.0
## ✔ purrr 1.0.1## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag() masks stats::lag()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errorslibrary(plm)##
## Attaching package: 'plm'
##
## The following objects are masked from 'package:dplyr':
##
## between, lag, leadtest.dat <- load("Ozone_Drought_Final.RData")
regions <- read.csv("region_code.csv")
df <- combinedAir.final
df$State.Code <- as.numeric(df$State.Code)
df$month = as.numeric(df$month)
df1 <- df %>%
mutate(no_drought = ifelse(USDM.categorical == "NoDrought", 1, 0),
moderate = ifelse(USDM.categorical == "ModerateDrought", 1, 0),
severe = ifelse(USDM.categorical == "SevereDrought", 1, 0))md2_random<- plm(Max.Ozone ~ moderate + severe + elevation + Longitude + Latitude,
data = df1,
index = c("GEOID", "Year","month"),
model = ("random"))## Warning in pdata.frame(data, index): duplicate couples (id-time) in resulting pdata.frame
## to find out which, use, e.g., table(index(your_pdataframe), useNA = "ifany")md2_fix<- plm(Max.Ozone ~ moderate + severe + elevation + Longitude + Latitude,
data = df1,
index = c("GEOID", "Year","month"),
effect = "twoway")## Warning in pdata.frame(data, index): duplicate couples (id-time) in resulting pdata.frame
## to find out which, use, e.g., table(index(your_pdataframe), useNA = "ifany")stargazer(md2_random, md2_fix, type = "text", title = "All Regions", align = TRUE,
add.lines = list(c("Model Name", "Random Effect Model", "Fixed Effect Model")),
column.labels = c("Random Effect Model", "Fixed Effect Model"),
dep.var.caption = "Dependent Variable: Ozone Max",
dep.var.labels.include = FALSE,
omit.table.layout = "#")##
## All Regions
## ================================================================
## Dependent Variable: Ozone Max
## ---------------------------------------------------
## Random Effect Model Fixed Effect Model
## ----------------------------------------------------------------
## moderate 1.309*** 1.416***
## (0.014) (0.014)
##
## severe 2.604*** 3.065***
## (0.020) (0.021)
##
## elevation 0.007*** 0.007***
## (0.00004) (0.00004)
##
## Longitude 0.142*** 0.323***
## (0.009) (0.032)
##
## Latitude 0.018 0.537***
## (0.024) (0.045)
##
## Constant 51.806***
## (1.232)
##
## ----------------------------------------------------------------
## Model Name Random Effect Model Fixed Effect Model
## Observations 6,608,332 6,608,332
## R2 0.011 0.009
## Adjusted R2 0.011 0.009
## F Statistic 51,505.070*** 11,699.010*** (df = 5; 6607385)
## ================================================================
## Note: *p<0.1; **p<0.05; ***p<0.01# Adding regions column
df2 <- df1 %>%
merge(regions, by = "State.Code")ne <- df2 %>%
filter(noaa_region == "northeast")
ne_random <- plm(Max.Ozone ~ moderate + severe + elevation + Longitude + Latitude,
data = ne,
index = c("GEOID", "Year", "month"),
model = ("random"))## Warning in pdata.frame(data, index): duplicate couples (id-time) in resulting pdata.frame
## to find out which, use, e.g., table(index(your_pdataframe), useNA = "ifany")ne_fix <- plm(Max.Ozone ~ moderate + severe + elevation + Longitude + Latitude,
data = ne,
index = c("GEOID", "Year", "month"),
effect = "twoway")## Warning in pdata.frame(data, index): duplicate couples (id-time) in resulting pdata.frame
## to find out which, use, e.g., table(index(your_pdataframe), useNA = "ifany")stargazer(ne_random, ne_fix, type = "text", title = "Northeast", align = TRUE,
add.lines = list(c("Model Name", "Random Effect Model", "Fixed Effect Model")),
column.labels = c("Random Effect Model", "Fixed Effect Model"),
dep.var.caption = "Dependent Variable: Ozone Max",
dep.var.labels.include = FALSE,
omit.table.layout = "#")##
## Northeast
## ============================================================
## Dependent Variable: Ozone Max
## -----------------------------------------------
## Random Effect Model Fixed Effect Model
## ------------------------------------------------------------
## moderate 0.497*** -0.495***
## (0.050) (0.053)
##
## severe -1.826*** -4.222***
## (0.162) (0.167)
##
## elevation 0.011*** 0.013***
## (0.0002) (0.0002)
##
## Longitude 1.855*** 4.292***
## (0.087) (0.147)
##
## Latitude -3.442*** -4.644***
## (0.141) (0.203)
##
## Constant 319.751***
## (10.608)
##
## ------------------------------------------------------------
## Model Name Random Effect Model Fixed Effect Model
## Observations 960,225 960,225
## R2 0.005 0.005
## Adjusted R2 0.005 0.005
## F Statistic 3,234.418*** 903.199*** (df = 5; 960049)
## ============================================================
## Note: *p<0.1; **p<0.05; ***p<0.01nr <- df2 %>%
filter(noaa_region == "northern_rockies")
nr_random <- plm(Max.Ozone ~ moderate + severe + elevation + Longitude + Latitude,
data = nr,
index = c("GEOID", "Year", "month"),
model = ("random"))## Warning in pdata.frame(data, index): duplicate couples (id-time) in resulting pdata.frame
## to find out which, use, e.g., table(index(your_pdataframe), useNA = "ifany")nr_fix <- plm(Max.Ozone ~ moderate + severe + elevation + Longitude + Latitude,
data = nr,
index = c("GEOID", "Year", "month"),
effect = "twoway")## Warning in pdata.frame(data, index): duplicate couples (id-time) in resulting pdata.frame
## to find out which, use, e.g., table(index(your_pdataframe), useNA = "ifany")stargazer(nr_random, nr_fix, type = "text", title = "Northern Rockies", align = TRUE,
add.lines = list(c("Model Name", "Random Effect Model", "Fixed Effect Model")),
column.labels = c("Random Effect Model", "Fixed Effect Model"),
dep.var.caption = "Dependent Variable: Ozone Max",
dep.var.labels.include = FALSE,
omit.table.layout = "#")##
## Northern Rockies
## ============================================================
## Dependent Variable: Ozone Max
## -----------------------------------------------
## Random Effect Model Fixed Effect Model
## ------------------------------------------------------------
## moderate 0.512*** 0.180***
## (0.048) (0.052)
##
## severe 0.803*** -0.018
## (0.076) (0.084)
##
## elevation 0.005*** 0.005***
## (0.0002) (0.0002)
##
## Longitude -0.111** -1.073***
## (0.045) (0.094)
##
## Latitude -0.270*** 0.705***
## (0.080) (0.142)
##
## Constant 35.327***
## (5.288)
##
## ------------------------------------------------------------
## Model Name Random Effect Model Fixed Effect Model
## Observations 226,906 226,906
## R2 0.032 0.004
## Adjusted R2 0.032 0.003
## F Statistic 1,181.555*** 171.206*** (df = 5; 226832)
## ============================================================
## Note: *p<0.1; **p<0.05; ***p<0.01nw <- df2 %>%
filter(noaa_region == "northwest")
nw_random <- plm(Max.Ozone ~ moderate + severe + elevation + Longitude + Latitude,
data = nw,
index = c("GEOID", "Year", "month"),
model = ("random"))## Warning in pdata.frame(data, index): duplicate couples (id-time) in resulting pdata.frame
## to find out which, use, e.g., table(index(your_pdataframe), useNA = "ifany")nw_fix <- plm(Max.Ozone ~ moderate + severe + elevation + Longitude + Latitude,
data = nw,
index = c("GEOID", "Year", "month"),
effect = "twoway")## Warning in pdata.frame(data, index): duplicate couples (id-time) in resulting pdata.frame
## to find out which, use, e.g., table(index(your_pdataframe), useNA = "ifany")stargazer(nw_random, nw_fix, type = "text", title = "Northwest", align = TRUE,
add.lines = list(c("Model Name", "Random Effect Model", "Fixed Effect Model")),
column.labels = c("Random Effect Model", "Fixed Effect Model"),
dep.var.caption = "Dependent Variable: Ozone Max",
dep.var.labels.include = FALSE,
omit.table.layout = "#")##
## Northwest
## ============================================================
## Dependent Variable: Ozone Max
## -----------------------------------------------
## Random Effect Model Fixed Effect Model
## ------------------------------------------------------------
## moderate 0.471*** 0.228**
## (0.082) (0.104)
##
## severe -0.993*** -2.230***
## (0.190) (0.209)
##
## elevation 0.011*** 0.012***
## (0.0002) (0.0002)
##
## Longitude -2.543*** -4.407***
## (0.137) (0.171)
##
## Latitude -2.535*** -3.378***
## (0.277) (0.385)
##
## Constant -157.097***
## (18.929)
##
## ------------------------------------------------------------
## Model Name Random Effect Model Fixed Effect Model
## Observations 105,338 105,338
## R2 0.046 0.041
## Adjusted R2 0.046 0.041
## F Statistic 4,064.467*** 906.615*** (df = 5; 105285)
## ============================================================
## Note: *p<0.1; **p<0.05; ***p<0.01ohv <- df2 %>%
filter(noaa_region == "ohio_valley")
ohv_random <- plm(Max.Ozone ~ moderate + severe + elevation + Longitude + Latitude,
data = ohv,
index = c("GEOID", "Year", "month"),
model = ("random"))## Warning in pdata.frame(data, index): duplicate couples (id-time) in resulting pdata.frame
## to find out which, use, e.g., table(index(your_pdataframe), useNA = "ifany")ohv_fix <- plm(Max.Ozone ~ moderate + severe + elevation + Longitude + Latitude,
data = ohv,
index = c("GEOID", "Year", "month"),
effect = "twoway")## Warning in pdata.frame(data, index): duplicate couples (id-time) in resulting pdata.frame
## to find out which, use, e.g., table(index(your_pdataframe), useNA = "ifany")stargazer(ohv_random, ohv_fix, type = "text", title = "Ohio Valley", align = TRUE,
add.lines = list(c("Model Name", "Random Effect Model", "Fixed Effect Model")),
column.labels = c("Random Effect Model", "Fixed Effect Model"),
dep.var.caption = "Dependent Variable: Ozone Max",
dep.var.labels.include = FALSE,
omit.table.layout = "#")##
## Ohio Valley
## =============================================================
## Dependent Variable: Ozone Max
## ------------------------------------------------
## Random Effect Model Fixed Effect Model
## -------------------------------------------------------------
## moderate 1.872*** 0.153***
## (0.045) (0.046)
##
## severe 3.628*** 0.627***
## (0.105) (0.108)
##
## elevation 0.005*** 0.004***
## (0.0003) (0.0003)
##
## Longitude 0.627*** 3.304***
## (0.069) (0.199)
##
## Latitude -0.917*** -1.994***
## (0.117) (0.230)
##
## Constant 132.649***
## (8.026)
##
## -------------------------------------------------------------
## Model Name Random Effect Model Fixed Effect Model
## Observations 1,006,480 1,006,480
## R2 0.008 0.001
## Adjusted R2 0.008 0.0005
## F Statistic 3,183.263*** 134.398*** (df = 5; 1006272)
## =============================================================
## Note: *p<0.1; **p<0.05; ***p<0.01sth <- df2 %>%
filter(noaa_region == "south")
sth_random <- plm(Max.Ozone ~ moderate + severe + elevation + Longitude + Latitude,
data = sth,
index = c("GEOID", "Year", "month"),
model = ("random"))## Warning in pdata.frame(data, index): duplicate couples (id-time) in resulting pdata.frame
## to find out which, use, e.g., table(index(your_pdataframe), useNA = "ifany")sth_fix <- plm(Max.Ozone ~ moderate + severe + elevation + Longitude + Latitude,
data = sth,
index = c("GEOID", "Year", "month"),
effect = "twoway")## Warning in pdata.frame(data, index): duplicate couples (id-time) in resulting pdata.frame
## to find out which, use, e.g., table(index(your_pdataframe), useNA = "ifany")stargazer(sth_random, sth_fix, type = "text", title = "South", align = TRUE,
add.lines = list(c("Model Name", "Random Effect Model", "Fixed Effect Model")),
column.labels = c("Random Effect Model", "Fixed Effect Model"),
dep.var.caption = "Dependent Variable: Ozone Max",
dep.var.labels.include = FALSE,
omit.table.layout = "#")##
## South
## ============================================================
## Dependent Variable: Ozone Max
## -----------------------------------------------
## Random Effect Model Fixed Effect Model
## ------------------------------------------------------------
## moderate 1.363*** 0.961***
## (0.035) (0.038)
##
## severe 2.938*** 2.183***
## (0.052) (0.059)
##
## elevation 0.014*** 0.014***
## (0.001) (0.001)
##
## Longitude 0.450*** -0.014
## (0.069) (0.169)
##
## Latitude 0.500*** 3.482***
## (0.070) (0.197)
##
## Constant 63.687***
## (6.890)
##
## ------------------------------------------------------------
## Model Name Random Effect Model Fixed Effect Model
## Observations 951,486 951,486
## R2 0.015 0.003
## Adjusted R2 0.015 0.003
## F Statistic 4,518.735*** 507.127*** (df = 5; 951335)
## ============================================================
## Note: *p<0.1; **p<0.05; ***p<0.01se <- df2 %>%
filter(noaa_region == "southeast")
se_random <- plm(Max.Ozone ~ moderate + severe + elevation + Longitude + Latitude,
data = se,
index = c("GEOID", "Year", "month"),
model = ("random"))## Warning in pdata.frame(data, index): duplicate couples (id-time) in resulting pdata.frame
## to find out which, use, e.g., table(index(your_pdataframe), useNA = "ifany")se_fix <- plm(Max.Ozone ~ moderate + severe + elevation + Longitude + Latitude,
data = se,
index = c("GEOID", "Year", "month"),
effect = "twoway")## Warning in pdata.frame(data, index): duplicate couples (id-time) in resulting pdata.frame
## to find out which, use, e.g., table(index(your_pdataframe), useNA = "ifany")stargazer(se_random, se_fix, type = "text", title = "Southeast", align = TRUE,
add.lines = list(c("Model Name", "Random Effect Model", "Fixed Effect Model")),
column.labels = c("Random Effect Model", "Fixed Effect Model"),
dep.var.caption = "Dependent Variable: Ozone Max",
dep.var.labels.include = FALSE,
omit.table.layout = "#")##
## Southeast
## ==============================================================
## Dependent Variable: Ozone Max
## -------------------------------------------------
## Random Effect Model Fixed Effect Model
## --------------------------------------------------------------
## moderate 3.868*** 2.831***
## (0.034) (0.037)
##
## severe 5.240*** 3.592***
## (0.051) (0.055)
##
## elevation 0.006*** 0.009***
## (0.0004) (0.0005)
##
## Longitude -0.043 -0.026
## (0.060) (0.169)
##
## Latitude 0.651*** -0.398**
## (0.060) (0.197)
##
## Constant 15.664***
## (6.029)
##
## --------------------------------------------------------------
## Model Name Random Effect Model Fixed Effect Model
## Observations 982,391 982,391
## R2 0.030 0.008
## Adjusted R2 0.030 0.008
## F Statistic 20,731.880*** 1,660.708*** (df = 5; 982203)
## ==============================================================
## Note: *p<0.1; **p<0.05; ***p<0.01sw <- df2 %>%
filter(noaa_region == "southwest")
sw_random <- plm(Max.Ozone ~ moderate + severe + elevation + Longitude + Latitude,
data = sw,
index = c("GEOID", "Year", "month"),
model = ("random"))## Warning in pdata.frame(data, index): duplicate couples (id-time) in resulting pdata.frame
## to find out which, use, e.g., table(index(your_pdataframe), useNA = "ifany")sw_fix <- plm(Max.Ozone ~ moderate + severe + elevation + Longitude + Latitude,
data = sw,
index = c("GEOID", "Year", "month"),
effect = "twoway")## Warning in pdata.frame(data, index): duplicate couples (id-time) in resulting pdata.frame
## to find out which, use, e.g., table(index(your_pdataframe), useNA = "ifany")stargazer(sw_random, sw_fix, type = "text", title = "Southwest", align = TRUE,
add.lines = list(c("Model Name", "Random Effect Model", "Fixed Effect Model")),
column.labels = c("Random Effect Model", "Fixed Effect Model"),
dep.var.caption = "Dependent Variable: Ozone Max",
dep.var.labels.include = FALSE,
omit.table.layout = "#")##
## Southwest
## ==============================================================
## Dependent Variable: Ozone Max
## -------------------------------------------------
## Random Effect Model Fixed Effect Model
## --------------------------------------------------------------
## moderate -0.444*** 0.089**
## (0.034) (0.040)
##
## severe 3.715*** 4.295***
## (0.043) (0.055)
##
## elevation 0.008*** 0.008***
## (0.0001) (0.0001)
##
## Longitude 0.255*** 1.148***
## (0.063) (0.082)
##
## Latitude 0.510*** 1.668***
## (0.069) (0.100)
##
## Constant 42.653***
## (7.484)
##
## --------------------------------------------------------------
## Model Name Random Effect Model Fixed Effect Model
## Observations 683,530 683,530
## R2 0.045 0.032
## Adjusted R2 0.045 0.032
## F Statistic 22,580.700*** 4,493.770*** (df = 5; 683441)
## ==============================================================
## Note: *p<0.1; **p<0.05; ***p<0.01um <- df2 %>%
filter(noaa_region == "upper_midwest")
um_random <- plm(Max.Ozone ~ moderate + severe + elevation + Longitude + Latitude,
data = um,
index = c("GEOID", "Year", "month"),
model = ("random"))## Warning in pdata.frame(data, index): duplicate couples (id-time) in resulting pdata.frame
## to find out which, use, e.g., table(index(your_pdataframe), useNA = "ifany")um_fix <- plm(Max.Ozone ~ moderate + severe + elevation + Longitude + Latitude,
data = um,
index = c("GEOID", "Year", "month"),
effect = "twoway")## Warning in pdata.frame(data, index): duplicate couples (id-time) in resulting pdata.frame
## to find out which, use, e.g., table(index(your_pdataframe), useNA = "ifany")stargazer(um_random, um_fix, type = "text", title = "Upper Midwest", align = TRUE,
add.lines = list(c("Model Name", "Random Effect Model", "Fixed Effect Model")),
column.labels = c("Random Effect Model", "Fixed Effect Model"),
dep.var.caption = "Dependent Variable: Ozone Max",
dep.var.labels.include = FALSE,
omit.table.layout = "#")##
## Upper Midwest
## ============================================================
## Dependent Variable: Ozone Max
## -----------------------------------------------
## Random Effect Model Fixed Effect Model
## ------------------------------------------------------------
## moderate 0.149** -1.140***
## (0.062) (0.065)
##
## severe 0.283 -2.467***
## (0.207) (0.210)
##
## elevation 0.029*** 0.001
## (0.002) (0.002)
##
## Longitude 1.245*** 11.215***
## (0.080) (0.356)
##
## Latitude 1.309*** 6.775***
## (0.130) (0.212)
##
## Constant 87.542***
## (8.051)
##
## ------------------------------------------------------------
## Model Name Random Effect Model Fixed Effect Model
## Observations 382,355 382,355
## R2 0.006 0.005
## Adjusted R2 0.006 0.005
## F Statistic 646.255*** 393.032*** (df = 5; 382246)
## ============================================================
## Note: *p<0.1; **p<0.05; ***p<0.01wst <- df2 %>%
filter(noaa_region == "west")
wst_random <- plm(Max.Ozone ~ moderate + severe + elevation + Longitude + Latitude,
data = wst,
index = c("GEOID", "Year", "month"),
model = ("random"))## Warning in pdata.frame(data, index): duplicate couples (id-time) in resulting pdata.frame
## to find out which, use, e.g., table(index(your_pdataframe), useNA = "ifany")wst_fix <- plm(Max.Ozone ~ moderate + severe + elevation + Longitude + Latitude,
data = wst,
index = c("GEOID", "Year", "month"),
effect = "twoway")## Warning in pdata.frame(data, index): duplicate couples (id-time) in resulting pdata.frame
## to find out which, use, e.g., table(index(your_pdataframe), useNA = "ifany")stargazer(wst_random, wst_fix, type = "text", title = "West", align = TRUE,
add.lines = list(c("Model Name", "Random Effect Model", "Fixed Effect Model")),
column.labels = c("Random Effect Model", "Fixed Effect Model"),
dep.var.caption = "Dependent Variable: Ozone Max",
dep.var.labels.include = FALSE,
omit.table.layout = "#")##
## West
## ===============================================================
## Dependent Variable: Ozone Max
## --------------------------------------------------
## Random Effect Model Fixed Effect Model
## ---------------------------------------------------------------
## moderate 1.056*** 1.171***
## (0.031) (0.042)
##
## severe 0.993*** 1.449***
## (0.036) (0.061)
##
## elevation 0.006*** 0.006***
## (0.0001) (0.0001)
##
## Longitude 0.305*** 0.283***
## (0.044) (0.045)
##
## Latitude 0.036 0.109
## (0.069) (0.070)
##
## Constant 74.249***
## (5.543)
##
## ---------------------------------------------------------------
## Model Name Random Effect Model Fixed Effect Model
## Observations 1,290,716 1,290,716
## R2 0.013 0.011
## Adjusted R2 0.013 0.011
## F Statistic 15,613.030*** 2,977.455*** (df = 5; 1290634)
## ===============================================================
## Note: *p<0.1; **p<0.05; ***p<0.01