Standard error on the mean (rare events/low rates)

Question

How do I calculate the standard error on, $r$, an aggregated rate? (for a sample of $n$)

I will go through two simple cases first (1 and 2) to put in perspective what I mean by "aggregated rate" which will be explained in case 3.

1. Binary outcome - expected number of successes

Starting simple.

i     x
      (1/0)
-----------
1     0
2     0
3     1
...   ...
N     0

If the probability of an attack happening to a subject is $p$, and there are $n$ subjects, then the expected number of attacks is $np$ with a standard error on the expected value $\sqrt{p(1-p)/n}$. If $p=1$%, then the vast majority of subjects will not have the event happen to them ($x=0$).

2. Poisson outcome - expected count of events

Slightly more complex - instead of looking at a binary outcome, the counts can be examined.

i     x
      Pois(0.01)
----------------
1     0
2     1
3     2
...   ...
N     0

The mean number of times a subject is attacked is $0.01$. A fairly valid model would be to treat the number of attacks as a Poisson distribution with $\lambda = 0.01$. If a sample of $n$ is taken from this underlying population, the mean number of attacks is $\lambda=0.01$ attacks subject$^{-1}$ with a standard error on the mean of $1/\sqrt{\lambda n}$. When $\lambda = 0.01$, the vast majority of subjects will have experienced $x=0$ attacks.

3. Discrete-continuous hybrid - expected rate of events

More complex - instead of looking at counts, looking at rates.

i     x           w      r
      Pois(0.01)  miles  rate
--------------------------------
1     0           10     0
2     1            5     0.2
3     2           20     0.1
...   ...         ...    ...
N     0           50     0

Instead of looking at attack numbers, I am now looking at attack rates (e.g., attacks per mile, or attacks per visit to London, etc). Again, the vast majority of subjects have a rate of $0$ attacks mile$^{-1}$ owing to the rarity of the event in the first place.

The aggregated rate would be:

$$ r = \frac{\sum^N_{i=1} x_i}{\sum^N_{i=1} w_i} $$

But how do I calculate the standard error on, $r$, an aggregated rate? (for a sample of $n$)

I found a similar question here, but could not quite find what I was looking for.

Answer 1

Here's an example of how you might do this with a Poisson regression and a logarithmic offset for area in Stata (though this should be doable in any statistics package). The model looks like:

$$E[ attacks] = \exp\left(\alpha+ 1 \cdot \ln{miles} \right)= \exp\left(\alpha \right) \cdot miles,$$

so the expected value of the rate is $$E \left[ \frac{attacks}{miles} \right] = \exp\left(\alpha \right).$$

The intuition is that the expected number of attacks is proportional to area, but the rate is constant.

Once you have the intercept and its standard error, it is pretty easy to calculate the rate and its standard error using the delta method.

. clear

. /* Fake Data */
. set seed 101011979

. set obs 1000
number of observations (_N) was 0, now 1,000

. generate miles    = mod(_n, 10)*10 + 10

. generate attacks  = rpoisson(.05)*miles

. generate rate     = attacks/mile

. generate ln_miles = ln(miles)

. 
. /* Poisson Model */
. constraint define 1 _b[ln_miles] = 1

. poisson attacks ln_miles, constraint(1) nolog

Poisson regression                              Number of obs     =      1,000
                                                Wald chi2(0)      =          .
Log likelihood = -7393.2626                     Prob > chi2       =          .

 ( 1)  [attacks]ln_miles = 1
------------------------------------------------------------------------------
     attacks |      Coef.   Std. Err.      z    P>|z|     [95% Conf. Interval]
-------------+----------------------------------------------------------------
    ln_miles |          1  (constrained)
       _cons |  -3.205332   .0211762  -151.36   0.000    -3.246836   -3.163827
------------------------------------------------------------------------------

. 
. generate pred_attacks    = exp(_b[_cons])*miles                   // predicted number of attacks given the miles                

. generate pred_attacks_se = exp(_b[_cons])*_se[_cons]*miles // SE of that prediction using delta method

. generate pred_rate       = exp(_b[_cons])                  // predicted number of attacks with miles set to 1 (aka the rate), constant  

. generate pred_rate_se    = exp(_b[_cons])*_se[_cons]       // std. error of the rate using delta method

. 
. summarize pred_attacks pred_rate attacks rate

    Variable |        Obs        Mean    Std. Dev.       Min        Max
-------------+---------------------------------------------------------
pred_attacks |      1,000        2.23    1.165162   .4054545   4.054545
   pred_rate |      1,000    .0405455           0   .0405455   .0405455
     attacks |      1,000        2.23    13.17196          0        200
        rate |      1,000        .038    .1964551          0          2