Confidence interval for the difference of two negative binomial rates I have a model with a negative binomial distribution using the glm.nb function from R:
Call:
MASS::glm.nb(formula = Counts ~ Gender + offset(log(Offset_Days)), 
    data = r_dataset, init.theta = 1.023811633, link = log)

Deviance Residuals: 
    Min       1Q   Median       3Q      Max  
-0.5973  -0.5973  -0.5594  -0.5594   3.7897  

Coefficients:
            Estimate Std. Error  z value Pr(>|z|)    
(Intercept) -7.67759    0.02613 -293.787  < 2e-16 ***
GenderMale   0.14226    0.03450    4.124 3.73e-05 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for Negative Binomial(1.0238) family taken to be 1)

    Null deviance: 12494  on 22136  degrees of freedom
Residual deviance: 12477  on 22135  degrees of freedom
AIC: 22558

Number of Fisher Scoring iterations: 1


              Theta:  1.0238 
          Std. Err.:  0.0863 

 2 x log-likelihood:  -22551.7100

From here I get that the ratio of the two gender groups shows the males having rates about 14% greater than females.
My question is if there is a way to derive a confidence interval for the difference of the rates? From the model I get that the rates for the male and female groups are 0.1948688 0.1690280 respectively. From this data I can get a confidence interval of the ratio, but I'm wondering if there is a way to get a confidence interval for the difference, 0.1948688 - 0.1690280 = 0.0258408?
Also, I only have one covariate here, but in the final model there will be several, both continuous and discrete.
 A: Two ways I see.
Marginal Effects
You're basically asking about a contrast.  The {marginaleffects} library handles this elegantly.
library(tidyverse)
library(rsample)
library(marginaleffects)

y <- c(7, 5, 4, 7, 5, 2, 11, 5, 5, 4, 2, 3, 4, 3, 5, 9, 6, 7, 10, 6, 12,
       6, 3, 5, 3, 9, 13, 0, 6, 1, 2, 0, 1, 0, 0, 4, 5, 1, 5, 3, 3, 4)

sex <- sample(c("M","F"), size = length(y), replace = T)

d <- tibble(y, sex)


# Marginal Effects

fit <- MASS::glm.nb(y ~ ., data = d)

marginaleffects::marginaleffects(fit, variables = 'sex') %>% 
  summary


  Term Contrast  Effect Std. Error z value Pr(>|z|)  2.5 % 97.5 %
1  sex    M - F -0.6364      1.024 -0.6216  0.53419 -2.643   1.37

Model type:  negbin 
Prediction type:  response 

Here, the marginaleffects function call computes the Male/Female contrast.
Bootstrapping
# bootstrap

compute_differences <- function(data, ...){
  d <- analysis(data)
  fit <- MASS::glm.nb(y ~ ., data = d)
  
  # Predict for males
  newdata <-  tibble(sex = c('F','M'))
  preds <- predict(fit, newdata = newdata, type='response')
  
  tibble(term = 'difference', estimate = diff(preds))
}


bootstraps(d, times = 250, apparent=TRUE) %>% 
  mutate(
    difference = map(splits, compute_differences)
  ) %>% 
  int_bca(difference, .fn = compute_differences)

 term       .lower .estimate .upper .alpha .method
  <chr>       <dbl>     <dbl>  <dbl>  <dbl> <chr>  
1 difference  -2.85     -1.03   1.26   0.05 BCa   


Note that both of these methods only really work conditional on all the other variables being the same.  If Sex is associated with a 14% change in the rate, then 14% of 100 is a lot more than 14% of 10.
A: *

*Draw, with replacement, a random sample from your observations. The size of the sample should be the same as that of your observations.


*For this random sample, fit the model as you have above.


*For the fitted model in step 2, compute the difference in rates. Call this a bootstrap replication of the difference.


*Repeat the thre steps above many times, yielding many bootstrap replications.


*The middle 95 % of the bootstrap replications (when they are placed in order of size) is the 95 % confidence interval of the difference.
I mention this method first because it is very general and will work regardless of what number you're interested in. This method is also more accurate in the sense that it preserves the range and distribution of your number of interest.
However, specifically in your case, you have another approximative option: the variance of a difference between two numbers is the sum of the variances of each number. With the variance and a normality assumption that might not hold, you can get a confidence interval.
