Inconsistent results when moving FEs to main formula.

[Garcese]

First I'd like to say I'm a big fan of the fixest package, it is my go-to R package for all things regression. I ran into an issue, however, when trying to recover the constant in my regression with fixed effects. I have the following regression:

reg <- feols(dependent ~ x + y + z | FE, data = mydat, cluster = ~FE)

I wanted to see the constant for this model. I saw in a stack exchange post that this is possible by doing:

reg$sumFE %>% mean()

Intuitively, this makes sense, and it did seem to work. However, to verify this, I tried placing the FE inside the main regression formula instead, since feols() outputs the constant (and it's standard error, which the other method does not give) when no FEs are specified. So, I ran this:

reg <- feols(dependent ~ x + y + z + FE, data = mydat, cluster = ~FE)

This not only gave me a different constant, but I also realized the SEs for all the variables were different from the original regression, although the coefficients remained the same. My assumption would be that it wouldn't matter whether I proprely specified the FE as a FE, or just put it as a regular categorical variable in the main regression, but it seems that something different is going on under the hood depending on where in the function the FE is included. For most of the variables, the SE has changed extremely little, but for others it changed a lot.

There is probably a good reason why the outputs are different, and it's probably an obvious reason that is explained clearly in the docs. That being said, if someone could address this here, that would be greatly appreciated. I am most concerned about the changing SEs, but if you also have a solution for getting the correct constant (and it's SE), that would be excellent.

Also, I tried this without clustered SEs, and still got different SEs, so I don't think that has anything to do with it. If you have any questions about the data I'm using, or would like to actually see the differing outputs, please feel free to ask.

Thanks!

[Garcese]

I'd like to write a partial answer to my own question. After messing around some more, I've noticed two things:

1. When you have not specified the cluster argument, and you include fixed effects after a | in the feols() call, the function automatically clusters SEs by the first specified FE. This was the cause of the large discrepancies I was seeing between regression results, as I did not notice the automatic clustered SEs. I suppose I should have included some sample results, I myself would have probably noticed the difference just from pasting them here! However, there are still some marginal differences, which I explain in (2).

2. When you are clustering SEs by a FE, the SEs will change very slightly depending on whether you place the clustered FE in the main formula or after a | in the feols() call. Only moving the clustered FE will change the SEs, moving the other FEs out of or into the main regression changes nothing. They change very, very minimally, from my experimentation, enough to change p-values by maybe 0.001 +/-. I assume, here, there is something different going on under the hood that is causing these minute differences, but because the differences are so small, they really do not seem to be a concern.

Additionally, the reason reg$sumFE %>% mean() was not matching up to the same regression when the FEs are placed in the main formula was because of the obvious reason that the FEs are not being demeaned out, so you need to average all the fixed effects by the number of times they appear in the data and add that to constant to achieve the same result as reg$sumFE %>% mean(). Another econometric blunder by yours truly.

So that solves 95% of my questions. I'm going to leave this open, just to see if someone knows what is causing the minute differences in SEs I outlined in 2. I have a separate question about retrieving the constant and it's SE, but I am going to open another issue to address that.

[s3alfisc]

Hi @Garcese, do your errors vanish if you impose the same small sample corrections on the fixed effects and dummy regressions? I.e. what happens if you set 'ssc = ssc(fixef.K = "none")'? It might be that fixest computes a different number of effective regressions. I would assume that for the dummy regression, fixed.K is always none, whereas for the fixed effects regression, the default is different?

[Garcese]

Thanks for the reply @s3alfisc. So, I just re-ran the regression I was working with to test this, adding the ssc = ssc(fixef.K = "none") argument. I was still noticing small differences in the standard errors of my covariates. For example, when I placed the clustered FE into the main formula, I got a SE of 0.042898, and when specified as a FE (not in main formula), I got 0.042389 for the same covariate.

These SEs did not match with either of the regressions when I did not use ssc = ssc(fixef.K = "none"). I got SEs of 0.043445 and 0.042916. It seems like when you include the clustered FE after a | (not in main formula), the SEs get slightly smaller. The coefficient is much larger than the SEs in this case - these small differences equate to very small p-value changes.

Hopefully this helps. I can provide sample data if you or anyone would like to investigate this further. You should be able to replicate it yourself with any data though, following the procedure listed above, but I can also provide data and code if wanted.

[lrberge]

Hi everyone!

The difference in SEs comes from the degrees of freedom correction and is normal. When SEs are clustered the parameters from the fixed-effects nested in the clusters are removed.

That is, the K in the adjustment below is smaller, leading to smaller SEs.

The nestedness is not checked on regular variables, hence you get all the parameters in the K when you include the FEs as regular dummies. That's the difference.

Here is an illustration (btw Alex's intuition was right but the option was "full" and not "none" ;-)):

base = setNames(iris, c("y", "x1", "x2", "x3", "species"))
library(fixest)

est_noFE = feols(y ~ x1 + x2 + species, base, vcov = ~species)
est_FE = feols(y ~ x1 + x2 | species, base, vcov = ~species)
est_FE_full = feols(y ~ x1 + x2 | species, base, vcov = ~species + ssc(fixef.K = "full"))

etable(est_noFE, est_FE, est_FE_full, drop = "^species")
#>                         est_noFE           est_FE      est_FE_full
#> Dependent Var.:                y                y                y
#> 
#> Constant         2.390* (0.4510)
#> x1               0.4322 (0.1624)  0.4322 (0.1613)  0.4322 (0.1624)
#> x2              0.7756* (0.1274) 0.7756* (0.1265) 0.7756* (0.1274)
#> Fixed-Effects:  ---------------- ---------------- ----------------
#> species                       No              Yes              Yes
#> _______________ ________________ ________________ ________________
#> S.E.: Clustered      by: species      by: species      by: species
#> Observations                 150              150              150
#> R2                       0.86331          0.86331          0.86331
#> Within R2                     --          0.64151          0.64151
#> ---
#> Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

[Garcese]

Thanks for the answer! I knew there had to be something I was missing. It was such a small difference, I almost didn't catch it at first.