How to display the number of groups in fixest regression tables using modelsummary?

I'm using the fixest package for fixed effects regressions and modelsummary to create tables. I'd like the tables to show the number of groups (e.g., number of unique values for the fixed effect variables). Here's a minimal example:

library(fixest)
library(modelsummary)

# Example fixed effects regression
reg1 <- feols(Sepal.Length ~
                Sepal.Width + Petal.Length + Petal.Width |
                Species,
              data = iris
)

# Create summary table
msummary(reg1)

This produces a nice table, but I'd like to add the number of groups in my fixed effect (in this example, it would be 3 species). I've tried using gof_map to include this information, but it doesn't seem to extract the number of groups by default:

msummary(reg1, 
         gof_map = c("nobs", "r.squared", "adj.r.squared", "n"))

Is there a way to make modelsummary display the number of groups for each fixed effect in my model? I'd prefer a solution that works with the standard modelsummary functionality, but I'm open to custom functions if needed.

Note: In my actual research, I'm working with panel data where the number of groups is important for interpreting the results.

You could use the gof_function argument to pull the `fixef_sizes` attribute from the model. For example, using the "mtcars" dataset:

mod = feols(mpg ~ disp | am + gear + carb, mtcars)
msummary(mod, gof_function=\(model) data.frame(t(model$fixef_sizes)))

      --- snip ---