R: panel data and fixest package

Question

I am trying to incorporate panel data analysis to a negative binomial with fixest package. I am unsure whether I am implementing the 'panel.id' argument correctly. According to the fixest's CRAN there are to identical ways to incorporate panel data analysis. First, I could set up a panel data set with function panel.

Thus.

dfpanel <- panel(df1, ~id + year)

Alternatively, I could work with my original dataframe but use panel identifier argument panel.id=~id+year.

In theory, both the first two models should calculate the same standard errors.

m1 <- fenegbin(incident~ treatdummy + poverty + rural|state,
            data=dfpanel)

m2 <- fenegbin(incident~ treatdummy + poverty + rural|state,
panel.id=~id+year,
data=df1)

m3 <- fenegbin(incident~ treatdummy + poverty + rural|state,
data=df1)

Just to compare, I am running a third model with neither panel.id not panel data. If I am implementing panel analysis correct, m1 and m2 should have the same SEs

Nevertheless, I am getting that m1 and m2 have different SEs. m1 and m3 have the same SE, suggesting that I am missing something about how to implement panel data with fixest correctly.

Here is the output

> etable(m1, m2, m3)
                               m1                  m2                m3
Dependent Var.:          homicide            homicide          homicide
                                                                       
minedummy        0.3134* (0.1406)   0.3134** (0.1077)  0.3134* (0.1406)
imarg           -0.3467* (0.1504) -0.3467*** (0.0600) -0.3467* (0.1504)
ruralurban      1.116*** (0.2388)   1.116*** (0.0918) 1.116*** (0.2388)
Fixed-Effects:  ----------------- ------------------- -----------------
state                         Yes                 Yes               Yes
_______________ _________________ ___________________ _________________
S.E.: Clustered         by: state              by: id         by: state
Observations               30,888              30,888            30,888

What am I doing wrong?

Answer 1

Going off of @Parfait's comment, the SE cluster groupings are the culprit. You can make the SEs agree by explicitly specifying the cluster argument in m2. See the example below:

data("InstInnovation")
library(fixest)
II_panel = panel(InstInnovation, ~company+year)

form = as.formula(patents ~ cites | industry)
m1 = fenegbin(form, 
              data=II_panel)

m2 = fenegbin(form, data=InstInnovation, panel.id=~company+year)

m3 = fenegbin(form, data=InstInnovation)

# same as m2, but specify cluster explicitly
m4 = fenegbin(form, data=InstInnovation, panel.id=~company+year, cluster='industry')

etable(m1, m2, m3, m4) yields:

Dependent Var.:            patents            patents            patents            patents
                                                                                           
cites           0.0023*** (0.0006) 0.0023*** (0.0005) 0.0023*** (0.0006) 0.0023*** (0.0006)
Fixed-Effects:  ------------------ ------------------ ------------------ ------------------
industry                       Yes                Yes                Yes                Yes
_______________ __________________ __________________ __________________ __________________
S.E.: Clustered       by: industry        by: company       by: industry       by: industry
Observations                 6,208              6,208              6,208              6,208
Squared Cor.               0.05756            0.05756            0.05756            0.05756
Pseudo R2                  0.11051            0.11051            0.11051            0.11051
BIC                       37,037.2           37,037.2           37,037.2           37,037.2
Over-dispersion            0.43466            0.43466            0.43466            0.43466
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Where you can see that m1, m3, and m4 all have the same standard errors. That said, I agree that this is somewhat confusing behavior and may be worth raising in the Git issues as suggested.