Getting predicted values with confidence interval on Negative binomial regression with fixed effects I'm using the R fixest package to fit a negative binomial regression model with three fixed effect variables. The function I am using is fenebin().
I want to draw a plot with my predicted values and a confidence interval similarly to this post: https://fromthebottomoftheheap.net/2018/12/10/confidence-intervals-for-glms/
Since the fixest package provides the standard-error of the predicted value with GLM models, how can I get these standard errors?
Any help would be much appreciated!

Unfortunately, the current version (0.11.1) of fixest doesn't seem to implement the computation of prediction standard errors for GLMs. Although this has not been stated on the help page of predict (?predict.fixest), if you try to compute these standard errors using a GLM, you'll get an error message.
Here is an R demonstration of this using Poisson regression.
library(fixest)
library(dplyr)
wasp <- read.csv('darlingtonia.csv', sep = ",", skip = 1)
wasp <- mutate(wasp, lvisited = as.logical(visited))

res = feglm(lvisited ~ leafHeight, wasp, "poisson")
obs_fe = predict(res)

# generate fake new data
ndata <- with(wasp, data_frame(leafHeight = seq(min(leafHeight),
                                                max(leafHeight),
                                                length = 100)))
# try to extrat prediction se
ndata <- add_column(ndata,
                    wrong_se = predict(object = res,
                                       newdata = ndata,
                                       type = 'response',
                                       se.fit = TRUE)$se.fit)

Error in predict.fixest(object = res, newdata = ndata, type = "response",  : 
  The standard-error of the prediction is currently only available for OLS models, sorry.

