I simulated a simple poisson regression model as follows

$x \sim Normal(2, 1)$

$\mu = \exp(0.6 + 0.9 x)$

$y \sim \text{Poisson}(\mu)$

n = 1001
x = rnorm(n, 2)
mu = exp(0.6 + 0.9 * x)
sim <- rpois(n, lambda = mu)

However, the mean and variance are not equal

# 16.68731
# 349.5331

but the mean is quite close to the standard deviation:

# 18.69581

I expected the mean to equal the standard deviation. Why is that not occurring here? How can I fix it?

I eventually want to be able to simulate a Poisson regression with multiple continuous covariates and at least 2 categorical covariates, Possibly with hierarchical effects. So any code to extend to that case would be great as well.

  • 4
    $\begingroup$ This is confused on several levels. The main problem is already explained in the first answer. Here is one more: Even for an unconditional Poisson, the mean does not equal the SD unless exceptionally the SD and the variance coincide at 1. $\endgroup$
    – Nick Cox
    Commented Jul 27, 2022 at 16:13

1 Answer 1


Your simulated data sim are not Poisson. That is, they are not unconditionally Poisson, but only conditionally on the predictor. Thus, equidispersion does not hold - or, more precisely, equidispersion will hold conditionally on the predictor: if you simulate many times for a fixed value of x (or equivalently, mu), the result will be equidispersed. But then, that will not come as a surprise.

What you have is a case of a compound distribution, specifically, a Poisson-lognormal distribution: the parameter to a Poisson distribution is itself lognormally distributed.

  • 3
    $\begingroup$ Great explanation to something I completely overlooked there. Makes total sense. Thanks! $\endgroup$ Commented Jul 27, 2022 at 16:34

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

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