feglm(y ~ x + z,
df, family = 'probit')
This does not seem to work in R. This is odd. The fixest package is really great. But for two step regressions such as the Heckman in a panel context, would involve a first stage probit. It doesn't seem to work with or without fixed effects. I would hate to have to move back to the plm package just to perform a probit. Any insight on whether people have to switch back and forth between these two packages would also be helpful.
This website claims that this works (https://stata2r.github.io/fixest/#models), although I cannot see any documentation explaining this phenomenon.