You are correct in noting that, adding factor(FIPS)
in side the regression formula should be equivalent to specifiying factor(FIPS) + YEARxQT
in the felm()
fixed effects block. I am less certain that your reg2
model that throws up the NA's is working correctly however, as that would typically throw up an error message before getting that result. Without actually seeing your data, I can only presume your issue is related to a bug, or has something very specific to do with your data.
Fortunately, there are lots of different ways to do this estimation, and it might be worth it to try this using different packages to corroborate the correct estimation and better understand what is going on with your specific use of felm()
.
For example, just to check, I downloaded some sample panel data and ran a FE estimation using a binary for the main regressor. I run a Least Square Dummy Variable model, a FE using the plm
package, and then two versions with felm
, both placing the entity fixed effect inside and outside of felm
's FE block as you did in your example. In all four cases the coefficient estimate was exactly the same (which it should be).
library(foreign)
library(plm)
library(lfe)
Panel <- read.dta("http://dss.princeton.edu/training/Panel101.dta")
modlsdv <-lm(y ~ op + factor(country) + factor(year), data=Panel)
modplm <- plm(y ~ op + factor(year), data=Panel, index=c("country", "year"), model="within")
modfelm1 <- felm(y ~ op | country + year | 0 | country, data=Panel)
#in this last model, I move the entity FE outside the FE block
modfelm2 <- felm(y ~ op + factor(country) | year | 0 | country, data=Panel)
I only show the coefficient of interest as all the factored-out variables take up a lot of space, but you can see below that they are all equivalent (apart from the SE's in the felm
models being different because of how they were clustered).
────────────────────────────────────────────────────────────────────────────
LSDV plm felm1 felm2
─────────────────────────────────────────────────────────────
(Intercept) 419655690.33
3
(1320579554.
769)
op -1434277502. -1434277502. -1434277502. -1434277502.
244 244 244 244
* *
(771947977.3 (771947977.3 (524652684.0 (594276313.4
01) 01) 46) 03)
─────────────────────────────────────────────────────────────
N 70 70 70 70
R2 0.386 0.264 0.386 0.386
logLik -1609.649 -1609.649 -1609.649
AIC 3255.298 3255.298 3255.298
──────────────────────────────────────────────────────────────────────────
reg1 <- felm(HPI~FLOODED+ lag(dat.reg[,'HPI'],1) + factor(FIPS) | YEARxQT | 0 | FIPS, data= dat.reg)
. The one that returned NA:reg2 <- felm(HPI~FLOODED+ lag(dat.reg[,'HPI'],1) | FIPS+YEARxQT | 0 | FIPS, data= dat.reg)
. $\endgroup$reg2
call give you a Warning message? $\endgroup$