# 8 Panel Data

## 8.1 Seminar

Create a new script with the following lines at the top, and save it as `seminar8.R`

```
rm(list = ls())
setwd("~/PUBL0055")
```

### 8.1.1 Required Packages

Let’s start off by installing the `plm`

package that allows us to estimate linear models for panel data. You only need to do this once, so don’t put it in your script, instead type it in the R console.

`install.packages("plm")`

Now let’s load these packages:

```
library(plm)
library(lmtest)
library(texreg)
```

### 8.1.2 More Guns, Less Crime

Gun rights advocate John Lott argues in his book More Guns, Less Crime that crime rates in the United States decrease when gun ownership restrictions are relaxed. The data used in Lott’s research compares violent crimes, robberies, and murders across 50 states to determine whether the so called “shall” laws that remove discretion from license granting authorities actually decrease crime rates. So far 41 states have passed these “shall” laws where a person applying for a licence to carry a concealed weapon doesn’t have to provide justification or “good cause” for requiring a concealed weapon permit.

Let’s load the dataset used by Lott and see if we can test the arguments made by gun rights advocates.

`guns <- read.csv("guns.csv")`

The variables we’re interested in are described below.

Variable | Description |
---|---|

mur | Murder rate (incidents per 100,000) |

shall | = 1 if the state has a shall-carry law in effect in that year = 0 otherwise |

incarc_rate | Incarceration rate in the state in the previous year (sentenced prisoners per 100,000 residents; value for the previous year) |

pm1029 | Percent of state population that is male, ages 10 to 29 |

stateid | ID number of states (Alabama = 1, Alaska = 2, etc.) |

year | Year (1977-1999) |

We will focus on murder rates in this example but you could try the same with variables measuring violent crimes or robberies as well.

Let’s create a factor variable representing whether a state has passed “shall” law or not. The variable already exists as `0`

or `1`

but we want to convert it to a factor for our analysis.

`guns$shall <- factor(guns$shall, levels = c(0, 1), labels =c("NO", "YES"))`

### 8.1.3 Fixed Effects

Let’s estimate a fixed effect model on panel data using the `plm()`

function with `shall`

, `incarc_rate`

, and `pm1029`

as the independent variables.

```
state_effects <- plm(
mur ~ shall + incarc_rate + pm1029,
data = guns,
index = c("stateid", "year"),
effect = "individual"
)
summary(state_effects)
```

```
Oneway (individual) effect Within Model
Call:
plm(formula = mur ~ shall + incarc_rate + pm1029, data = guns,
effect = "individual", index = c("stateid", "year"))
Balanced Panel: n = 51, T = 23, N = 1173
Residuals:
Min. 1st Qu. Median 3rd Qu. Max.
-21.102428 -0.958945 0.016047 1.082008 29.031961
Coefficients:
Estimate Std. Error t-value Pr(>|t|)
shallYES -1.4513886 0.3154300 -4.6013 4.678e-06 ***
incarc_rate 0.0174551 0.0011261 15.4998 < 2.2e-16 ***
pm1029 0.9582993 0.0859610 11.1481 < 2.2e-16 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Total Sum of Squares: 12016
Residual Sum of Squares: 9800
R-Squared: 0.18444
Adj. R-Squared: 0.14581
F-statistic: 84.3526 on 3 and 1119 DF, p-value: < 2.22e-16
```

The `state_effects`

model shows that all three of our independent variables are statistically significant, with `shall`

decreasing murder rates by `1.45`

incidents per `100000`

members of the population. The effects of incarceration rate and percentage of male population between `10`

and `29`

years old are also statistically significant.

Before drawing any conclusions let’s make sure whether there are any state effects in our model using `plmtest()`

.

`plmtest(state_effects, effect="individual")`

```
Lagrange Multiplier Test - (Honda) for balanced panels
data: mur ~ shall + incarc_rate + pm1029
normal = 47.242, p-value < 2.2e-16
alternative hypothesis: significant effects
```

The p-value suggests the presence of state effects. In addition to state fixed effects, a number of factors could affect the murder rate that are not specific to an individual state. We can model these time fixed effects using the `effect = "time"`

argument in `plm()`

.

```
time_effects <- plm(
mur ~ shall + incarc_rate + pm1029,
data = guns,
index = c("stateid", "year"),
effect = "time"
)
summary(time_effects)
```

```
Oneway (time) effect Within Model
Call:
plm(formula = mur ~ shall + incarc_rate + pm1029, data = guns,
effect = "time", index = c("stateid", "year"))
Balanced Panel: n = 51, T = 23, N = 1173
Residuals:
Min. 1st Qu. Median 3rd Qu. Max.
-21.68350 -2.04596 -0.31955 1.76758 35.20084
Coefficients:
Estimate Std. Error t-value Pr(>|t|)
shallYES -0.2521605 0.3032163 -0.8316 0.4058
incarc_rate 0.0412157 0.0007704 53.4993 <2e-16 ***
pm1029 0.2148597 0.1407613 1.5264 0.1272
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Total Sum of Squares: 65571
Residual Sum of Squares: 18141
R-Squared: 0.72334
Adj. R-Squared: 0.71731
F-statistic: 999.623 on 3 and 1147 DF, p-value: < 2.22e-16
```

The `incarc_rate`

variable is the only statistically significant variable in the time fixed effects model.

Now let’s run `plmtest`

on the `time_effects`

model to verify if time fixed effects are indeed present in the model.

`plmtest(time_effects, effect="time")`

```
Lagrange Multiplier Test - time effects (Honda) for balanced
panels
data: mur ~ shall + incarc_rate + pm1029
normal = 16.104, p-value < 2.2e-16
alternative hypothesis: significant effects
```

The *p-value* tells us that we can reject the null hypothesis so we know that there are time fixed effects present in our model.

We already confirmed the presence of state fixed effects in the first model we estimated. Now, in order to control for both state AND time fixed effects, we need to estimate a model using the `effect = "twoways"`

argument.

```
twoway_effects <- plm(
mur ~ shall + incarc_rate + pm1029,
data = guns,
index = c("stateid", "year"),
effect = "twoways"
)
summary(twoway_effects)
```

```
Twoways effects Within Model
Call:
plm(formula = mur ~ shall + incarc_rate + pm1029, data = guns,
effect = "twoways", index = c("stateid", "year"))
Balanced Panel: n = 51, T = 23, N = 1173
Residuals:
Min. 1st Qu. Median 3rd Qu. Max.
-19.2097691 -0.9748749 -0.0069663 1.0119176 27.1354552
Coefficients:
Estimate Std. Error t-value Pr(>|t|)
shallYES -0.5640474 0.3325054 -1.6964 0.0901023 .
incarc_rate 0.0209756 0.0011252 18.6411 < 2.2e-16 ***
pm1029 0.7326357 0.2189770 3.3457 0.0008485 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Total Sum of Squares: 11263
Residual Sum of Squares: 8519.4
R-Squared: 0.24357
Adj. R-Squared: 0.19186
F-statistic: 117.746 on 3 and 1097 DF, p-value: < 2.22e-16
```

In a twoway fixed effects model `shall`

is no longer significant and the effect of male population between `10`

and `29`

years old has decreased from `0.95`

to `0.73`

incidents per `100,000`

population.

The results of all three models are shown below.

```
screenreg(
list(state_effects, time_effects, twoway_effects),
custom.model.names = c("State Fixed Effects", "Time Fixed Effects", "Twoway Fixed Effects")
)
```

```
==========================================================================
State Fixed Effects Time Fixed Effects Twoway Fixed Effects
--------------------------------------------------------------------------
shallYES -1.45 *** -0.25 -0.56
(0.32) (0.30) (0.33)
incarc_rate 0.02 *** 0.04 *** 0.02 ***
(0.00) (0.00) (0.00)
pm1029 0.96 *** 0.21 0.73 ***
(0.09) (0.14) (0.22)
--------------------------------------------------------------------------
R^2 0.18 0.72 0.24
Adj. R^2 0.15 0.72 0.19
Num. obs. 1173 1173 1173
==========================================================================
*** p < 0.001, ** p < 0.01, * p < 0.05
```

### 8.1.4 Serial Correlation

For time series data we need to address the potential for serial correlation in the error term. We will test for serial correlation with Breusch-Godfrey test using `pbgtest()`

and provide solutions for correcting it if necessary.

`pbgtest(twoway_effects)`

```
Breusch-Godfrey/Wooldridge test for serial correlation in panel
models
data: mur ~ shall + incarc_rate + pm1029
chisq = 765.16, df = 23, p-value < 2.2e-16
alternative hypothesis: serial correlation in idiosyncratic errors
```

The null hypothesis for the Breusch-Godfrey test is that there is no serial correlation. The `p-value`

from the test tells us that we can reject the null hypothesis and confirms the presence of serial correlation in our error term.

We can correct for serial correlation using `coeftest()`

similar to how we corrected for heteroskedastic errors. We’ll use the `vcovHC()`

function for obtaining a heteroskedasticity-consistent covariance matrix, but since we’re interested in correcting for autocorrelation as well, we will specify `method = "arellano"`

which corrects for both heteroskedasticity and autocorrelation.

```
twoway_effects_hac <- coeftest(
twoway_effects,
vcov = vcovHC(twoway_effects, method = "arellano", type = "HC3")
)
screenreg(
list(twoway_effects, twoway_effects_hac),
custom.model.names = c("Twoway Fixed Effects", "Twoway Fixed Effects (HAC)")
)
```

```
=============================================================
Twoway Fixed Effects Twoway Fixed Effects (HAC)
-------------------------------------------------------------
shallYES -0.56 -0.56
(0.33) (0.48)
incarc_rate 0.02 *** 0.02 *
(0.00) (0.01)
pm1029 0.73 *** 0.73
(0.22) (0.54)
-------------------------------------------------------------
R^2 0.24
Adj. R^2 0.19
Num. obs. 1173
=============================================================
*** p < 0.001, ** p < 0.01, * p < 0.05
```

We can see that with heteroskedasticity and autocorrelation consistent (HAC) standard errors, the percent of male population (10 - 29 yr old) is no longer a significant predictor in our model.

### 8.1.5 Cross Sectional Dependence

If a federal law imposed restrictions on gun ownership or licensing requirements then the changes would likely affect all 50 states. This is an example of Cross Sectional Dependence and not accounted for in a fixed effect model. Other scenarios could also trigger cross sectional dependence that we should take into consideration. For example, security policies and law enforcement efforts might change after an extraordinary event (think of mass shootings or terrorist attacks) thus influencing law enforcement practices in all states. We can check for cross sectional dependence using the Pesaran cross sectional dependence test or `pcdtest()`

.

`pcdtest(twoway_effects)`

```
Pesaran CD test for cross-sectional dependence in panels
data: mur ~ shall + incarc_rate + pm1029
z = 3.9121, p-value = 9.148e-05
alternative hypothesis: cross-sectional dependence
```

As we’ve seen with other tests, the null hypothesis is that there is no cross sectional dependence. The p-value, however tells that there is indeed cross-sectional dependence and we need to correct it.

**Driscoll and Kraay (1998) (SCC)**: We use the cross-sectional and serial correlation (SCC) method by Driscoll and Kraay for obtaining heteroskedasticity and autocorrelation consistent errors that are also robust to cross-sectional dependence. We can get SCC corrected covariance matrix using the `vcovSCC()`

function.

```
twoway_effects_scc <- coeftest(
twoway_effects,
vcov = vcovSCC(twoway_effects, type="HC3", cluster = "group")
)
twoway_effects_scc
```

```
t test of coefficients:
Estimate Std. Error t value Pr(>|t|)
shallYES -0.564047 0.542698 -1.0393 0.29888
incarc_rate 0.020976 0.010321 2.0324 0.04236 *
pm1029 0.732636 0.551066 1.3295 0.18396
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
```

```
screenreg(
list(twoway_effects, twoway_effects_scc),
custom.model.names = c("Twoway Fixed Effects", "Twoway Fixed Effects (SCC)")
)
```

```
=============================================================
Twoway Fixed Effects Twoway Fixed Effects (SCC)
-------------------------------------------------------------
shallYES -0.56 -0.56
(0.33) (0.54)
incarc_rate 0.02 *** 0.02 *
(0.00) (0.01)
pm1029 0.73 *** 0.73
(0.22) (0.55)
-------------------------------------------------------------
R^2 0.24
Adj. R^2 0.19
Num. obs. 1173
=============================================================
*** p < 0.001, ** p < 0.01, * p < 0.05
```

### 8.1.6 Exercises

- Create a new file called
`assignment8.R`

in your`PUBL0055`

folder and write all the solutions in it. Load the Comparative Political Dataset

`CPDS_1960-2013_stata.dta`

from your`PUBL0055`

folder or from the datasets repository:`cpds <- read.dta("https://uclspp.github.io/datasets/data/CPDS_1960-2013_stata.dta")`

Variable Description vturn Voter turnout in election judrev Judicial review (existence of an independent body which decides whether laws conform to the constitution).

Coded 0 = no, 1 = yes.ud Net union membership as a proportion wage and salary earners in employment (union density) unemp Unemployment rate, percentage of civilian labour force. unemp_pmp Cash expenditure for unemployment benefits as a percentage of GDP (public and mandatory private). realgdpgr Growth of real GDP, percent change from previous year. debt Gross general government debt (financial liabilities) as a percentage of GDP. socexp_t_pmp Total public and mandatory private social expenditure as a percentage of GDP. - Estimate a model for the electoral fractionalization of the party system as coded by the
`rae_ele`

variable in the dataset using all the variables listed above.- Estimate a fixed effect model and test for country and time fixed effects.
- Run the necessary tests to check whether country and time fixed effects are present.

- Estimate a twoway model and compare to the previous country and time fixed effect models.
- Test for serial correlation and cross sectional dependence in the twoway model.
- If either serial correlation or cross sectional dependence is present, use the methods discussed in the seminar to obtain heteroskedastic and autocorrelation consistent standard errors.
- Compare the HAC and spatially robust standard errors with the twoway model estimated earlier.
- Display the results in publication-ready tables and discuss the substantively significant findings.
- Save your script, which should now include the answers to all the exercises.
Source your script, i.e. run the entire script all at once. Fix the script if you get any error messages.