# Clear the workspace
rm(list = ls()) # Clear environment
gc() # Clear unused memory## used (Mb) gc trigger (Mb) limit (Mb) max used (Mb)
## Ncells 520270 27.8 1151834 61.6 NA 669302 35.8
## Vcells 958697 7.4 8388608 64.0 16384 1839795 14.1 cat("\f") # Clear the console dev.off # Clear the charts## function (which = dev.cur())
## {
## if (which == 1)
## stop("cannot shut down device 1 (the null device)")
## .External(C_devoff, as.integer(which))
## dev.cur()
## }
## <bytecode: 0x7f95e9109110>
## <environment: namespace:grDevices># Create a list of all packages to be used
packages <- c("knitr",
"AER",
"stargazer",
"plm")
#Install Packages
for (i in 1:length(packages)) {# Iterate through each package in the 'packages' vector.
if (!packages[i] %in% rownames(installed.packages())) {#check if package is installed
install.packages(packages[i], #install package if it isn't already
repos = "http://cran.rstudio.com/",
dependencies = TRUE)}
library(packages[i], character.only = TRUE)} # Load the package into the R session.## 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##
## 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# Remove packages and i
rm(packages, i)
# Import "Guns" dataset from AER
data("Guns")
# Create a duplicate table that keeps data in original format
guns <- Guns| Variable | Description |
|---|---|
| state | factor indicating state |
| year | factor indicating year |
| violent | violent crime rate (incidents per 100,000 members of the population) |
| murder | murder rate (incidents per 100,000) |
| robbery | robbery rate (incidents per 100,000) |
| prisoners | incarceration rate in the state in the previous year (sentenced prisoners per 100,000 residents; value for the previous year) |
| afam | percent of state population that is African-American, ages 10 to 64 |
| cauc | percent of state population that is Caucasian, ages 10 to 64 |
| male | percent of state population that is male, ages 10 to 29 |
| population | state population, in millions of people |
| income | real per capita personal income in the state (US dollars) |
| density | population per square mile of land area, divided by 1,000 |
| law | factor. Does the state have a shall carry law in effect in that year? |
# Observe the structure of the dataset
str(guns)## 'data.frame': 1173 obs. of 13 variables:
## $ year : Factor w/ 23 levels "1977","1978",..: 1 2 3 4 5 6 7 8 9 10 ...
## $ violent : num 414 419 413 448 470 ...
## $ murder : num 14.2 13.3 13.2 13.2 11.9 10.6 9.2 9.4 9.8 10.1 ...
## $ robbery : num 96.8 99.1 109.5 132.1 126.5 ...
## $ prisoners : int 83 94 144 141 149 183 215 243 256 267 ...
## $ afam : num 8.38 8.35 8.33 8.41 8.48 ...
## $ cauc : num 55.1 55.1 55.1 54.9 54.9 ...
## $ male : num 18.2 18 17.8 17.7 17.7 ...
## $ population: num 3.78 3.83 3.87 3.9 3.92 ...
## $ income : num 9563 9932 9877 9541 9548 ...
## $ density : num 0.0746 0.0756 0.0762 0.0768 0.0772 ...
## $ state : Factor w/ 51 levels "Alabama","Alaska",..: 1 1 1 1 1 1 1 1 1 1 ...
## $ law : Factor w/ 2 levels "no","yes": 1 1 1 1 1 1 1 1 1 1 ...# Make all data numeric
Guns[] <- lapply(Guns, as.numeric)
#show how data was transformed
kable(head(guns), caption = "ORIGINAL FORMAT")| year | violent | murder | robbery | prisoners | afam | cauc | male | population | income | density | state | law |
|---|---|---|---|---|---|---|---|---|---|---|---|---|
| 1977 | 414.4 | 14.2 | 96.8 | 83 | 8.384873 | 55.12291 | 18.17441 | 3.780403 | 9563.148 | 0.0745524 | Alabama | no |
| 1978 | 419.1 | 13.3 | 99.1 | 94 | 8.352101 | 55.14367 | 17.99408 | 3.831838 | 9932.000 | 0.0755667 | Alabama | no |
| 1979 | 413.3 | 13.2 | 109.5 | 144 | 8.329575 | 55.13586 | 17.83934 | 3.866248 | 9877.028 | 0.0762453 | Alabama | no |
| 1980 | 448.5 | 13.2 | 132.1 | 141 | 8.408386 | 54.91259 | 17.73420 | 3.900368 | 9541.428 | 0.0768288 | Alabama | no |
| 1981 | 470.5 | 11.9 | 126.5 | 149 | 8.483435 | 54.92513 | 17.67372 | 3.918531 | 9548.351 | 0.0771866 | Alabama | no |
| 1982 | 447.7 | 10.6 | 112.0 | 183 | 8.514000 | 54.89621 | 17.51052 | 3.925229 | 9478.919 | 0.0773185 | Alabama | no |
kable(head(Guns),caption = "FORMAT FOR SUMMARY")| year | violent | murder | robbery | prisoners | afam | cauc | male | population | income | density | state | law |
|---|---|---|---|---|---|---|---|---|---|---|---|---|
| 1 | 414.4 | 14.2 | 96.8 | 83 | 8.384873 | 55.12291 | 18.17441 | 3.780403 | 9563.148 | 0.0745524 | 1 | 1 |
| 2 | 419.1 | 13.3 | 99.1 | 94 | 8.352101 | 55.14367 | 17.99408 | 3.831838 | 9932.000 | 0.0755667 | 1 | 1 |
| 3 | 413.3 | 13.2 | 109.5 | 144 | 8.329575 | 55.13586 | 17.83934 | 3.866248 | 9877.028 | 0.0762453 | 1 | 1 |
| 4 | 448.5 | 13.2 | 132.1 | 141 | 8.408386 | 54.91259 | 17.73420 | 3.900368 | 9541.428 | 0.0768288 | 1 | 1 |
| 5 | 470.5 | 11.9 | 126.5 | 149 | 8.483435 | 54.92513 | 17.67372 | 3.918531 | 9548.351 | 0.0771866 | 1 | 1 |
| 6 | 447.7 | 10.6 | 112.0 | 183 | 8.514000 | 54.89621 | 17.51052 | 3.925229 | 9478.919 | 0.0773185 | 1 | 1 |
# Create a Summary table to check if data is balanced
summary_data <- data.frame(
Mean = round(colMeans(Guns, na.rm = TRUE),2),
Median = round(sapply(Guns, median, na.rm = TRUE),2),
SD = round(sapply(Guns, sd, na.rm = TRUE),2),
Min = round(sapply(Guns, min, na.rm = TRUE),2),
Max = round(sapply(Guns, max, na.rm = TRUE),2),
NAs = round((colSums(is.na(Guns)) / nrow(Guns)) * 100, 2))
#create the table
num_observations <- nrow(Guns)
caption_summary <- paste("SUMMARY DATA","(", num_observations, "OBSERVATIONS)")
kable(summary_data, caption = caption_summary)| Mean | Median | SD | Min | Max | NAs | |
|---|---|---|---|---|---|---|
| year | 12.00 | 12.00 | 6.64 | 1.00 | 23.00 | 0 |
| violent | 503.07 | 443.00 | 334.28 | 47.00 | 2921.80 | 0 |
| murder | 7.67 | 6.40 | 7.52 | 0.20 | 80.60 | 0 |
| robbery | 161.82 | 124.10 | 170.51 | 6.40 | 1635.10 | 0 |
| prisoners | 226.58 | 187.00 | 178.89 | 19.00 | 1913.00 | 0 |
| afam | 5.34 | 4.03 | 4.89 | 0.25 | 26.98 | 0 |
| cauc | 62.95 | 65.06 | 9.76 | 21.78 | 76.53 | 0 |
| male | 16.08 | 15.90 | 1.73 | 12.21 | 22.35 | 0 |
| population | 4.82 | 3.27 | 5.25 | 0.40 | 33.15 | 0 |
| income | 13724.80 | 13401.55 | 2554.54 | 8554.88 | 23646.71 | 0 |
| density | 0.35 | 0.08 | 1.36 | 0.00 | 11.10 | 0 |
| state | 26.00 | 26.00 | 14.73 | 1.00 | 51.00 | 0 |
| law | 1.24 | 1.00 | 0.43 | 1.00 | 2.00 | 0 |
Time Component: Year
Entity: State
Based on the summary table from above, the dataset appears to be balanced and does not have any apparent issues. One thing to note is that all factor variables were converted to numeric for the table. Thus, each distinct factor was converted to a numerical value in order of occurrence (Reason why 1977 was converted to 1). Additionally, the reason why there are 51 distinct observations for “state” is because it includes Washington DC.
\[ violent= \beta_0+ \beta_1(prisoners) + \beta_2(afam) + \beta_3(cauc) + \beta_4(male) +\beta_5(population) + \beta_6(income) + \beta_7(density) + \beta_8(law) + u_i \]
The regression above regresses one potential dependent variable from the dataset (violent crime rate) across all potential explanatory variables in the set (incarceration rate, demographics, etc). The idea is to see which variables, especially gun laws, have the most influence on violent crime.
# Run regression without fixed effects
reg_plain <- lm(data = guns, violent ~
prisoners + afam + cauc + male + population + income + density + law)
# Present regression
stargazer(reg_plain, title = "Without Fixed Effects",
align = TRUE, type = "text", header = FALSE)##
## Without Fixed Effects
## ===============================================
## Dependent variable:
## ---------------------------
## violent
## -----------------------------------------------
## prisoners 0.816***
## (0.044)
##
## afam 11.292
## (6.864)
##
## cauc 2.702
## (3.454)
##
## male 9.317**
## (4.441)
##
## population 18.506***
## (1.055)
##
## income 0.001
## (0.003)
##
## density 94.669***
## (5.428)
##
## lawyes -92.721***
## (13.426)
##
## Constant -178.778
## (224.010)
##
## -----------------------------------------------
## Observations 1,173
## R2 0.724
## Adjusted R2 0.722
## Residual Std. Error 176.313 (df = 1164)
## F Statistic 381.103*** (df = 8; 1164)
## ===============================================
## Note: *p<0.1; **p<0.05; ***p<0.01There are a few variables from this regression that have coefficients that contradict expectations for this data. The variables with contradicting direction are Caucasian and income. There is a stigma and a general assumption that predominantly white and wealthier areas have a higher rate of white collar crime (tax evasion and embezzlement) and a lower violent crime rate. However, since these variables and afam are not statistically significant, this model provides no conclusive evidence of any of the conclusions for said variables.
This model also yielded some expected results. First, is that gun laws reduce violent crime rates. By limiting the means of procuring a weapon, states also limit the potential to commit a crime. Also, population density makes sense to have a higher crime rate. More urban areas tend to have higher crime rates and more wealth disparity. The incarceration rate also makes sense but could potentially be a reverse causality as well. As more crimes are committed, the number of people being incarcerated is bound to increase.
This model may have omitted variable bias given that crime rates may vary at different economic periods and with different state legislation beyond gun laws. These could be controlled for by using fixed effects
\[ \text{violent}_{ist} - \bar{\text{violent}}_{i} = \beta_1 \times (\text{prisoners}_{ist} - \bar{\text{prisoners}}_{i}) + \beta_2 \times (\text{afam}_{ist} - \bar{\text{afam}}_{i}) + \beta_3 \times (\text{cauc}_{ist} - \bar{\text{cauc}}_{i}) + \beta_4 \times (\text{male}_{ist} - \bar{\text{male}}_{i}) + \beta_5 \times (\text{income}_{ist} - \bar{\text{income}}_{i}) + \beta_6 \times (\text{density}_{ist} - \bar{\text{density}}_{i}) + \beta_7 \times (\text{law}_{ist} - \bar{\text{law}}_{i}) +\beta_8(population_{ist}-\bar{population_i}) \epsilon_{ist} \]
Speaking to equation:
The equation above displays how fixed effects works. The main differences are that it doesn’t include a constant and each variable as well as each variable being measured in “\((x_{ist}-\bar{x_i})\)”.
\(x_{ist}\) represents the variable in “s” state at “t” time. A second regression was run that dropped the “s” fixed effect, making it \(x_{it}\).
\(\bar{x_i}\) represents the mean for each variable in the regression. Thus, the equation for fixed effects is the value of a variable at a certain point minus the mean for that variable.
Speaking to reasoning:
The reason for two separate fixed effects models is as follows:
State and Year:
State: As mentioned earlier, there could be legislation or other background factors in individual states that influence crime factors. A hurricane could displace residents or potential tax hikes could limit purchasing power and push them towards desperation.
Year: This variable is the focal point because there may be some periods that result in an increase in crime trends. These may include recessions and rises in political factions.
Just year:
# Run fixed effects regressions
#Fixed both state and year
reg_fe1 <- plm(violent ~
prisoners + afam + cauc + male + population + income + density + law,
data = guns,
model = "within",
index = c("state", "year"))
#Fixed for only year
reg_fe2 <- plm(violent ~
prisoners + afam + cauc + male + population + income + density + law,
data = guns,
model = "within",
index = "year")
# Present regression
stargazer(reg_plain, reg_fe1, reg_fe2, title = "OLS to Fixed Effects Comparison",
align = TRUE, type = "text", header = FALSE)##
## OLS to Fixed Effects Comparison
## ================================================================================================
## Dependent variable:
## ----------------------------------------------------------------------------
## violent
## OLS panel
## linear
## (1) (2) (3)
## ------------------------------------------------------------------------------------------------
## prisoners 0.816*** 0.097* 0.972***
## (0.044) (0.057) (0.046)
##
## afam 11.292 12.808 14.867**
## (6.864) (10.883) (7.268)
##
## cauc 2.702 10.323*** 4.279
## (3.454) (3.110) (3.638)
##
## male 9.317** -23.864*** -9.150
## (4.441) (3.925) (6.983)
##
## population 18.506*** 12.244** 17.487***
## (1.055) (5.347) (1.006)
##
## income 0.001 -0.004 0.009***
## (0.003) (0.004) (0.003)
##
## density 94.669*** -155.954*** 76.775***
## (5.428) (52.118) (5.561)
##
## lawyes -92.721*** -18.575 -49.962***
## (13.426) (11.563) (13.426)
##
## Constant -178.778
## (224.010)
##
## ------------------------------------------------------------------------------------------------
## Observations 1,173 1,173 1,173
## R2 0.724 0.204 0.749
## Adjusted R2 0.722 0.162 0.743
## Residual Std. Error 176.313 (df = 1164)
## F Statistic 381.103*** (df = 8; 1164) 35.645*** (df = 8; 1114) 426.247*** (df = 8; 1142)
## ================================================================================================
## Note: *p<0.1; **p<0.05; ***p<0.01I believe that the biggest problem is that I still included all explanatory variables for fixed effects by state. Since I only really wanted to observe the influence of gun laws, adding the other variables may over-complicate the model. Additionally, they may already be accounted for using fixed effects on the state level. When accounting for year only, however, there should be more variation given that they are the sole explanatory values for each variable.
state and year:
Magnitude for incarceration rate drops significantly
Increased positive magnitude of caucasian and now significant
Income flips signs but still insignificant
Male population flips to negative and becomes more significant
Density flips and is significant
Magnitude for gun laws drops significantly and becomes insignificant
This was a surprising outcome but you could make some conclusions off of it. One explanation could be that it may just be cancelled out by accounting for fixed effects by state. Since fixed effects control for all factors in a state, gun laws may already be accounted for.
Assuming that it wasn’t cancelled out, you may conclude that violent crimes can exist without the presence of guns. So if a person is going to commit a violent crime they may find an alternative weapon, or just slip through gun laws. This does not mean that the laws are ineffective however. This only suggests that the violent crime rate is not effected but gun laws may still reduce the outcome of each case. For instance, someone could pull a gun and kill someone or they may simply punch them. Both will increase the violent crime rate by the same amount, but one is objectively more severe than the other.
Just year:
This model seems a little more practical if you are seeking to measure all variables since, as mentioned earlier, fixed effects for state may push out or kind of cancel out variables. Here are the notable changes in regards to the original OLS:
Male population flips to negative and becomes insignificant
Positive value for income becomes significant
Gun laws magnitude is still reduced