I'm trying to do as much as I can with R instead of Stata, but Stata has positive network externalities in the economics profession and it's hard to move away if you have coauthors using it. Anyway, one of the most common regressions I have to run is a fixed effects regression with clustered standard errors. Suppose that Y is your dependent variable, X is an explanatory variable and F is a categorical variable that defines your fixed effects. You also want to cluster your standard errors according to F. If F has a lot of categories, demeaning instead of including the fixed effects explicitly has to be used so as not to have an excessive number of regressors. This would be the code in Stata:

xtset F

xtreg Y X, fe cluster(F) robust dfadj

The "dfadj" option is something I discovered recently that is needed to get the correct standard errors (it stands for "degrees of freedom adjustment" and at least for me, it was quite hidden in the manuals).

To produce the same results in R, you have to use the felm function of the lfe package. This is the code:

model <- felm(Y ~ X | F | 0 | F)

summary(model)