1
$\begingroup$

I am fairly new to stats and the field of confidence intervals, so I apologize if the question seems trivial.

I am working on a physics experiment and I am measuring the time to occurrence of a particular event. I know occurrence of this event follows a Poisson process.

My question is, if I calculate the average time to occurrence based on 10 such events, how do I calculate the upper and lower error bound (in other words, the confidence interval) of the estimated average. For instance, if I want the estimated average time to occurrence to be within $\pm 20\%$ of the actual mean based on a confidence level, how many samples do I need to record?

Someone told me that I need to use inverse chi-squared distribution, but I would like to have a more in depth tutorial/info in this.

Thank you very much for your time and help.

$\endgroup$
2
$\begingroup$

The first thing is that with a constant intensity Poisson-process with rate $\lambda$ per unit time, the time between events has an exponential distribution with mean $\mu=1/\lambda$.

Second, you actually ask two separate questions, so let's deal with them one at a time.

if I calculate the average time to occurrence based on 10 such events, how do I calculate the upper and lower error bound (in other words, the confidence interval) of the estimated average.

If you average ten such events, the mean is distributed as $\text{gamma}(10,\mu/10)$ (this is the gamma parameterized in the shape-scale form, not the shape-rate form)

If you want a $1-\alpha$ interval for the population mean, then one approach is to construct a pivotal quantity for the mean, which is to say a quantity $Q$ which is a function of the data and the parameter whose distribution doesn't depend on the parameter.

Note that $\hat{\mu}/\mu$ will also have a gamma distribution, it will be distributed as $\text{gamma}(10,1/10)$, so that is a pivotal quantity (and works just fine for this), but it may be easier to work with its inverse $\mu/\hat{\mu}$, which will be distributed as an inverse-gamma distribution, specifically $\text{inverse-gamma}(10,1/10)$.

Similarly, for a sample of size $n$, the distribution of $\mu/\hat{\mu}$ will be $\text{inverse-gamma}(n,1/n)$.

(The inverse chi-square is a particular case of this distribution, and it, too can be used for this purpose, with a suitable rescaling. So the advice you were given was correct, if not immediately helpful to you.)

(I'll assume you want a symmetric (in probability) interval.)

Let $l$ be the $\alpha/2$ quantile of an $\text{inverse-gamma}(n,1/n)$ distribution, and let $u$ be the $1-\alpha/2$ quantile - that is, the values that cut off a lower and upper tail area of $\alpha/2$.

Then $P(l<Q<u) = 1-\alpha$, so a $1-\alpha$ interval for $Q$ is $(l,u)$. Hence we can work out the interval for $\mu$:

$$l<\mu/\hat{\mu}<u$$

$$l\,\hat{\mu}<\mu<u\,\hat{\mu}$$

where $\hat{\mu}=\bar x$, which is to say a $1-\alpha$ interval for $\mu$ is $(l\,\bar x,u\, \bar x)$

Imagine the sample size was 10, and that the mean time was 3.21 seconds. Further, imagine you wanted a 95% interval.

I'm going to work in R, but anything that will give you gamma, inverse gamma, chi-square or inverse-chi-square quantiles will work.

So here's the 0.025 and 0.975 quantiles of a gamma(10,1/10):

 qgamma(.025,10,scale=.1)
[1] 0.4795389
> qgamma(.975,10,scale=.1)
[1] 1.70848

The corresponding cutoffs for an inverse gamma will be the inverses of these, flipped around:

> 1/qgamma(.975,10,scale=.1)
[1] 0.5853155
> 1/qgamma(.025,10,scale=.1)
[1] 2.085337

So that's our $l$ and $u$ respectively.

Hence a 95% interval for $\mu$ in this example is $(3.21\times 0.5853, 3.21\times 2.0853)$, or $(1.879, 6.694)$

That interval is far bigger that $\pm 20\%$, of course.

If you want to do it in terms of an inverse-chi-square instead, note that a $χ²(\nu)$ is a $\text{gamma}(\nu/2, 1/2)$; from there it's reasonably easy to convert the above calculations.

if I want the estimated average time to occurrence to be within ±20% of the actual mean based on a confidence level, how many samples do I need to record?

This is basically the inverse of the above problem, where we specify information about the interval and solve for $n$

Edit: the approach I took here before didn't quite solve the problem as stated, but a related problem (I solved the true mean being within 20% of the estimates). This version should be right:

An interval that's like this: $0.8 \mu <\hat \mu < 1.2 \mu$

is equivalent to (after dividing by $\mu\,\hat\mu$ and inverting):

$0.833 \hat\mu < \mu < 1.25 \hat\mu$

which is in the right form to play with.

So our interval is $l_n \hat{\mu}<\mu<u_n \hat{\mu}$.

You want to choose $n$ that $u_n>1.25$ and that $l_n<0.833$.

Unfortunately, the distribution of a $\text{gamma}(n,1/n)$ (and, correspondingly, and inverse gamma) is different for every $n$.

One could of course proceed by trial and error and this will work perfectly well, if sometimes being a little tedious.

Half a dozen guesses for the lower limit gets us to:

> n=86;1/qgamma(.025,n,scale=1/n)
[1] 1.250202
> n=87;1/qgamma(.025,n,scale=1/n)
[1] 1.248503

While for the other end:

> n=105;1/qgamma(.975,n,scale=1/n)
[1] 0.8332447
> n=106;1/qgamma(.975,n,scale=1/n)
[1] 0.8339306

Which means we need 106 observations to be within the desired limit on both ends and get the desired coverage.

Note that for the lower end we needed a larger sample size to achieve the desired width. [If you are willing to use an asymmetric interval (asymmetric in probability, I mean), in order to get a symmetric-in-percentage-error interval -- i.e. to be substantially more likely to be wrong out one end of the interval than the other -- then you can achieve a smaller sample size while still getting $1-\alpha$ coverage.]

Computing an asymmetric interval

If we try different values of $n$, we can find how much tail area is cut off either side of the +/-20%, and then find the smallest $n$ that totals to no more than $\alpha$ in the tails. In our example, the first place to try is to 'split the difference' and try $n=96$:

 n=96;pgamma(.8,n,scale=1/n);pgamma(1.2,n,scale=1/n,lower.tail=FALSE)
[1] 0.01907025
[1] 0.03033813

Those two tail areas give just under 0.05. Trying n=95 gives more than 0.05 in the tails combined, so n=96 is as small as you can go.

(Note that the 0.03033 value corresponds to the lower tail for the inverse gamma, but it's the probability of the estimate being at least 20% too high)

So if you are happy with more probability of being too high than too low, you can do n=96 for a 95% interval. [A higher coverage ($1-\alpha$) would give a larger $n$, and a lower coverage would give a smaller $n$.]


Simulation example in R to show that this approach is doing what it should. In this example, we calculate the limits ($l$ and $u$) for a symmetric-in-probability (2.5% in each tail) interval based on the 106 just computed:

n=106;l=1/qgamma(.975,n,scale=1/n);u=1/qgamma(.025,n,scale=1/n)
l;u
[1] 0.8339306
[1] 1.221422

Since $l$ and $u$ are within $(0.8333, 1.25)$ if the true mean lies inside that interval the estimate can't be further than 20% from the true mean and by construction, 95% of the time the generated interval includes the true mean, then at least 95% of the time we generate such an interval, the estimate must be within 20% of the true mean.

Now we simulate ten thousand intervals to show it does have the required coverage:

sum(replicate(10000,
    isTRUE(all.equal(c(l,u)*mean(rexp(106,rate=1/3.21))<3.21,c(TRUE,FALSE)))))
[1] 9507

- looks like it's right.


To avoid (or at least greatly reduce) the trial and error (not that it takes more than a few moments), we can make use of an approximation that will get us very close to the right answer - the Wilson-Hilferty transformation:

If $X \sim χ²(k)$ then $\scriptstyle\sqrt[3]{X/k}$ is approximately normally distributed with mean $\scriptstyle 1-2/(9k)$ and variance $\scriptstyle 2/(9k)$.

$\endgroup$
  • $\begingroup$ Thank you very much! This is exactly the kind of explanation I was looking for. I have a much better idea what I should do after reading your post. $\endgroup$ – user2675918 Dec 15 '13 at 22:41
  • $\begingroup$ If you'd like an example of the shorter, asymmetric-in-probability interval on the $\pm$20% case, let me know. I may end up doing one anyway, for completeness' sake. $\endgroup$ – Glen_b -Reinstate Monica Dec 15 '13 at 23:07
  • $\begingroup$ Sure, more examples are always better, if its not too much trouble. $\endgroup$ – user2675918 Dec 16 '13 at 1:07
  • $\begingroup$ Done, at least in outline. Considering the second question is more easily tackled with a pivot of $\hat\mu/\mu$, I may have been better using that from the start and dealing with a slightly more fiddly calculation at the start, to save more fiddling later. I may go back and rework the post with that in mind. To double check: I am correct in thinking that your $\pm$20% is asking that $\hat\mu/\mu$ lies between 0.8 and 1.2, right? $\endgroup$ – Glen_b -Reinstate Monica Dec 16 '13 at 2:09
  • $\begingroup$ Hi Glen you are right, that is what I had in mind. But really your original post is already very informative and I think we can easily do the calculation on our own after following your post. $\endgroup$ – user2675918 Dec 17 '13 at 1:06

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.