# How large does a Poisson distribution's mean need to be to use normal distribution statistics?

As the mean of a Poisson distribution increases, the Poisson distribution approximates a normal distribution. I assume that once the Poisson mean becomes large enough, we can use normal distribution statistics. Therefore we can start saying things '68% of the distribution will lie within 1 standard deviation of the Poisson mean' once the mean of a Poisson distribution becomes large enough.

How large does the mean of a Poisson distribution need to be before we can start using normal distribution statistics?

• It all depends on which statistics you want to use. Do you have anything in mind besides the one you quoted? And how accurate do you want the approximation(s) to be? – whuber Jun 23 '15 at 22:02
• Depends on not only what kinds of things you're doing but also the region of the distribution you're doing it in (tails tend to take longer) and your degree of tolerance for inaccuracy. The question is likely too broad as presently framed. – Glen_b Jun 24 '15 at 2:37

Update: As @whuber has pointed out with his comments, a better way to look at this is by computing the true coverage probabilities for the Poisson. The simulation, while also with its uses, does not reveal the interesting pattern seen in the plot below.

This was based on @whuber's code (see his first comment on this answer):

f <- function(mu) ppois(mu + sqrt(mu), mu) - ppois(mu - sqrt(mu), mu) ## computing coverage probabilities as a function of the true mean
curve(f, from = 1, to = 2000, n = 300, main = "True Coverage Probabilities for Poisson", xlab = expression(lambda), ylab = "Probability")
abline(h = 0.6827, col = "red") ## coverage prob. for a Normal RV


What this does: If $X \sim \textrm{Pois}(\lambda)$ then $E(X) = \lambda$ and $V(X) = \lambda$. This means that the interval in question is $I := (\lambda - \sqrt \lambda, \lambda + \sqrt \lambda)$. The function f computes $$\mathbb P_\lambda(X \in I) = F_X(\lambda + \sqrt \lambda; \lambda) - F_X(\lambda - \sqrt \lambda; \lambda)$$ where in R $F_X(t; \lambda)$ is obtained via the ppois function.

Original answer: This is in no way a categorical answer but I thought you might like to see a simulation. Note that I'm using samples of size $n = 20000$ because you didn't mention that you cared about the sample size, so I wanted each sample to reflect asymptotic properties.

The simulation shows that the Poisson random variables (RVs) do not behave indistinguishably from the Normal RVs until around $\lambda \approx 1000$ with respect to the coverage rate and this choice of $n$. We can also see the variation in the coverage of random samples of Normal RVs even though they all exactly have the property that we are investigating at the population level. Note that in this simulation I compared a random sample to its sample mean and sample standard deviation rather than the population mean and population standard deviation. I chose to do so because I felt this to be more interesting for a discussion about the distribution of a statistic calculated from a sample.

Here's the code to make the plot.

set.seed(1)
lambda.seq <- round(seq(10, 2000, length = 300)) ## lambdas to try
res.norm <- res.pois <- numeric(length(lambda.seq)) ## these get the results
nsim <- 20000 ## number of observations at each iteration

for(i in 1:length(lambda.seq))
{
sims <- rpois(nsim, lambda.seq[i]) ## simulating Poissons
res.pois[i] <- mean(sims > mean(sims) - sd(sims) & sims < mean(sims) + sd(sims)) ## computing proportion of Poisson RVs within 1 SE of the mean

sims <- rnorm(nsim, lambda.seq[i], lambda.seq[i]) ## simulating normals
res.norm[i] <- mean(sims > mean(sims) - sd(sims) & sims < mean(sims) + sd(sims))
}

plot(res.pois ~ lambda.seq, pch = 19, col = "red", main = "1 SD Coverage Convergence of Poisson to Normal", xlab = expression(lambda), ylab = expression(paste("Proportion  within  ", mu, "  \u00b1  ", sigma, sep = "")))
points(res.norm ~ lambda.seq, pch = 19, col = "blue")
abline(h = .6827) ## this is the true coverage for a normal distribution
legend("topright", c("Poisson", "Normal"), pch = 19, col = c("red", "blue"), bty = "n")

• Why did you elect to perform a simulation, which is computationally expensive and gives only approximate results, when you could produce exact results with much less time and programming effort? f <- function(mu) ppois(mu + sqrt(mu), mu) - ppois(mu - sqrt(mu), mu); curve(f(x), 1, 100, 1001). I'm sure readers would appreciate a brief description of what your code is doing and, if possible, an interpretation of the output and an explanation. Also consider clarifying what you mean by "reliably". (I happen to think the 68% approximation is pretty reliable even for $\lambda=1$.) – whuber Jun 23 '15 at 23:38
• whuber I personally think that's a nice answer – luciano Jun 24 '15 at 5:23
• @whuber: Of course you are correct that this could have been done without approximation and more efficiently. I went this route because the question was about statistics so I felt that looking at sampling distributions was more interesting. I like seeing the randomness in coverage proportions for truly Normal samples along with the Poisson samples. Also I'm not particularly concerned about the computational aspect for such a tiny simulation. – jld Jun 24 '15 at 18:08
• I appreciate that explanation. Normally I would find it more compelling, but the problem here is that a very interesting and detailed pattern present in the exact, correct calculations is obscured by the variation in the simulations, to the point of misleading readers. – whuber Jun 24 '15 at 18:43
• @whuber Good point. I've updated the answer to reflect your comments. Please let me know if you feel that I've misrepresented/failed to sufficiently credit you. – jld Jun 24 '15 at 20:14

Here's another perspective. Let's say you want to approximate cdf values of a Poisson with a normal fit. You didn't specify a desired accuracy, so I'll just give an example. Suppose you want to be within 0.01 of the true cdf at all points.

Empirically, if $\lambda$ is the mean of the Poisson, I observe that $\lambda \geq 44$ is enough to get the stated accuracy by matching the mean and variance of the normal to the Poisson and using a continuity correction.

To really answer we need to know what you want to do! If you want a confidence interval, say, more might be needed than for standard errors. But why do you want a normal approximation? You can construct confidence intervals directly from the likelihood function, which probably is better.

And, note that for the poisson distribution, the variance is equal to the expectation. So with increasing $$\lambda$$ (the expectation) the variance is increasing as well. For confidence intervals you really want a pivotal quantity ( a pivot is a function of data and parameters with a distribution not depending on the parameters), or, if that is impossible (as in most discrete distribution cases), at least with constant expectation and variance.

One way to obtain that is to use a variance-stabilizing transformation, which in the Poisson case is the square root transformation. So, if $$X$$ is Poisson ($$\lambda$$) then $$\sqrt{X}$$ have variance approximately 1/4. A normal approximation for the transformed variable usually works better, as it is closer to pivotal.

Also have a look at the post GLM vs square root data transformation