Confidence interval for Poisson distribution coefficient

Question

This is an exam question, testing if water is bad - that is if a sample has more than 2000 E.coli in 100ml. We have taken $n$ samples denoted $X_i$, and model the samples as a Poisson distribution with parameter $\theta$, and as mutually independent.

1 Determine the maximum likelihood estimator for $\theta$, denote it $\hat\theta_n$ from now on.

Under the Poisson distribution: $$ P(X=k)=\frac{\theta^k}{k!}e^{-\theta}\text{ so the likelihood is } L(\theta)=\frac{\theta^{\sum x_i}\, e^{-n\theta}}{x_1!\ldots x_n!}$$

$$ \ln (L(\theta)=\ln \left(\frac{\theta^{\sum X_i}\, e^{-n\theta}}{x_1!\ldots x_n!}\right)= \sum_{i=1}^n x_i\ln(\theta)-n\theta+\text{ term with no theta} $$ Now differentiating w.r.t. $\theta$: $$ \frac{\partial}{\partial\theta}\ln(\theta)=\frac{1}{\theta}\sum x_i-n\text{ and a second time }\frac{\partial^2}{\partial\theta^2}\ln(\theta) = \frac{-1}{\theta^2}<0$$ So the maximum likelihood estimator $\hat\theta_n=n^{-1}\sum{x_i}$

2 Calculate the expected value and variance: I got $\theta$ and $\theta/n$ respectively - which I think are right.

3 Show for $\epsilon>0$ $$ P(|\hat\theta_n-\theta|>\epsilon\sqrt{\theta})\leq\frac{1}{n\epsilon^2}$$ This follows from Chebyshev's inequality with the constant as $\epsilon\sqrt{\theta}$.

Now the bit I can't match up with my notes and complete: construct a confidence interval of level $\alpha$. We have covered Bienayme-Tchebychev and Hoeffding for confidence intervals of proportions, but I don't think that is right here. Can I use an interval of the form: $$\left[\hat\theta_n\pm\frac{1}{\sqrt{n}}\sqrt{\hat\theta_n (1-\hat\theta_n)}\,\phi^{-1}(1-\frac{\alpha}{2}) \right]$$ where $\phi$ is the inverse cdf of the standard normal?

The final part is to construct a test of level $\alpha$ for the null hypothesis the water is bad with the data $n=20$, $\sum_{i=1}^{20} X_i=38600$, $\alpha=0.05$ What decision is taken?

Answer 1

You are right, for this, you cannot use a CI based on proportions. Actually, the method of constructing CI based on proportions is a consequence of the Central Limit Theorem for Binomial Parameter Estimation. It isn't something arbitrary.

We use a result from the Central Limit Theorem.

Suppose $X_1,X_2,...,X_n$ form an i.i.d. sample from a relatively well behaved distribution. Let $Y_n$ be the random sample mean and let $\mu=E(X_i)$. Let $g$ be a differentiable function so that $g'(\mu)\neq 0$. Then, $\sqrt{n}(g(Y_n)-g(\mu))\to N(0,(g'(\mu))^2Var(X_i))$ in distribution.

This is the Delta Method.

What we will do is we will try to find a confidence interval for $g(\mu)$ rather than $\mu$ and hence get a confidence interval for $\mu$. To do that, our variance term for the CI for $g(\mu)$ must be a constant, say $1$. So you choose $g$ such that $(g'(\mu))^2Var(X_i)=1$, your $\mu$ being the argument for $g$.

In our case, $\mu=\theta$ and $Var(X_i)=\theta$ so we want to find $g$ so that $g'(\theta)=\dfrac{1}{\sqrt{\theta}}$, implying $g(\theta)=2\sqrt{\theta}$ as one of the solutions (no need to worry about the constant term as it is irrelevent; your job is to find $g$ which satisfies this, not its properties).

Thus, frame a CI for $g(\theta)=2\sqrt{\theta}$ using the CLT, which you should be able to do. You want to look at:

$$P[-c\leq \sqrt{n}(2\sqrt{Y_n}-2\sqrt{\theta})\leq c]=1-\alpha$$

You get your $c$ based on $\alpha$ and hence you get a confidence interval for $\theta$.

Please take some time and try to use this method to find a suitable level $\alpha$ test. Try for 1 day, after that if you still cannot do it, please comment below the answer. All the technicalities have been covered, now you should be able to peform the main inference based on your CI.

Answer 2

I have been following this Question for several days, and it seems to me that some clarification is in order.

MLE for Poisson mean. Suppose $X_i,$ for $i = 1,\dots,n$ are iid $Pois(\lambda).$ You have correctly shown that the MLE of $\lambda$ is $$\hat \lambda = \bar X = \frac{1}{n}\sum_{i=1}^n = T/n,$$ where the total $T \sim Pois(\Lambda)$ and $\Lambda = n\lambda.$ For large $T$ (and hence presumably large $\Lambda$), $T$ is approximately $Norm(\mu = \Lambda, \sigma^2 = \Lambda).$

Two-sided test for Poisson mean. For $n = 20,$ a two-sided test of $H_0: \Lambda = 40000$ against $H_A: \Lambda \ne 40000$ rejects at level $\alpha = .05$ when $$|Z| = |T - 40000|/40000 > 1.96,$$ where 1.96 is the 97.5 percentile of standard normal.

CI based on test. "Inverting the test" to get a 95% confidence interval (CI) for $\Lambda$ means to find an interval of values of $\Lambda$ for which $H_0$ is not rejected at level $\alpha.$ Simplifying the results of solving a quadratic equation by conflating 1.96 and 2, one obtains the 95% CI $T + 2 \pm 2\sqrt{T + 1}.$

This result assumes the approximate normality of $T$, If $T$ is very large, as in your example, the further simplification to $T \pm 2\sqrt{T}$ is harmless. This gives the 95% CI $(38207.06, 38992.94)$ for $\Lambda.$ And, upon dividing by $n = 20$ the CI for $\lambda$ is $(1910.4, 1949.6)$ or $1930 \pm 19.6$ for $\lambda.$

Note: This differs from the result in your comment; I suspect you have forgotten to multiply the standard error by 1.96 (or 2) to get the margin of error.

CI not assuming normality. For small $T$ it is not reasonable to assume normality. In general, a very good 95% CI for $\lambda$ is based on percentiles 2.5 and 97.5 of a gamma distribution with shape parameter $T$ and rate parameter $n$ (scale $1/n$). In R the computation for your data is

 qgamma(c(.025,.975), 38600, 20)
 ## 1910.794 1949.301

For your data with large $T,$ the relevant gamma distribution is nearly normal, and the CI is essentially the same as we got above from the normal approximation. Owing to the skewness of Poisson and gamma distributions, this type of interval is not symmetrical about $\bar X.$ (This type of interval comes from a Bayesian context with Poisson likelihood and a noninformative 'improper' gamma prior, yielding a gamma posterior distribution.)

One-sided hypothesis test. For your test of hypothesis, which is presumably one-sided, you need to be specific about the null and alternative hypotheses. For example, if you are assuming normality and testing $H_0: \Lambda = \Lambda_0 = 40000$ against $H_A: \Lambda \ne \Lambda_0,$ then you would reject for $Z = (T - \Lambda_0)/\sqrt{\Lambda_0} > 1.645.$ For your data $Z = -0.035.$ There is no question that you will conclude that the 'water must be safe' because $T < 40000.$ Notice, however, water with $\Lambda = 40000$ might sometimes produce $T$ just a little above $40000.$ Then you would not reject the null hypothesis--and have the task of explaining to non-statisticians why the water is likely OK anyhow.