Stats for

Generalized and Multilevel Regression

Poisson Regression

Poisson regression predicts a count, like how often something happens.

Poisson regression takes a count outcome, such as how many panic attacks someone has, and tests how one or more predictors relate to the rate at which the events happen.

Example

Do people randomized to a brief CBT program have panic attacks less often than people randomized to a waitlist? Two hundred adults with panic disorder were randomized to CBT or a waitlist, then kept a daily panic diary for 21 to 56 days; we counted their attacks, allowing for baseline anxiety.

Count Outcome + One or More Predictors + Time Observed (Offset) → Poisson Regression
Tests
Population incidence rate ratios (how many times higher the rate of events is per unit of each predictor)
Observations
Independent
Parametric
Yes
Headline Statistic
\text{IRR}
Effect Size
Incidence Rate Ratio
Groups
1

The Basics

When to Use It

Use it when your outcome is a count of events (how many times something happened to each person) and you want to know how one or more predictors relate to how often it happens. If people were observed for different lengths of time, an offset lets the model compare rates. For example:

  • Do people randomized to a brief CBT program have panic attacks less often than people randomized to a waitlist?
  • Does loneliness predict how many alcoholic drinks older adults report in a week?
  • Do children with higher ADHD symptom ratings get more classroom disciplinary referrals per school term?

Vocabulary You Need

Definition Count Outcome

An outcome that records how many times something happened to or was done by each person over some stretch of time or space: panic attacks in a month, cigarettes smoked in a day, errors on a task. Counts are whole numbers that can't go below zero, so when they are small they pile up near zero with a long tail to the right, and their spread grows with their average. Those features break ordinary regression's assumptions, which is what Poisson and negative binomial regression are built for. With large counts (averages well above 10 or so) they often behave well enough for ordinary regression. Glossary →

Definition Offset (\ln(t_i))

People observed for longer have more chance to rack up events. An offset puts each person's exposure (days in a diary, hours of observation, the size of a population) into the model on the log scale with its coefficient fixed at 1, so the model predicts a count of t × rate and the coefficients describe the rate per unit of exposure. It's not a predictor: nothing is estimated for it. Leaving it out when exposures differ mixes up having more events with having been watched longer; adding log exposure as an ordinary predictor instead lets the data choose a coefficient other than 1, which answers a different question. Glossary →

Definition Incidence Rate Ratio (\text{IRR})

The exponentiated coefficient from a Poisson or negative binomial regression, e^b. It multiplies the expected rate of events: an IRR of 0.70 for a treatment means the treated group's rate is 70% of the comparison group's (30% lower), and an IRR of 1.06 for a one-point rise in a predictor means 6% more events per point. 1 means no difference. Because it multiplies, effects compound: two points is 1.06 × 1.06, not 1.12. It's a ratio of rates (events per unit of time), so it only reads as a ratio of counts when everyone is observed for the same length of time. Also called the rate ratio; some software labels it Exp(B). Glossary →

Definition Dispersion Statistic (\hat{\phi})

The Pearson chi-square of a fitted Poisson model divided by its residual degrees of freedom, φ̂ = Σ (yᵢ − μ̂ᵢ)²/μ̂ᵢ ÷ (n − k − 1). Each term compares a person's squared miss with the variance the model expects, which for a Poisson count is the mean itself, so φ̂ near 1 means the counts vary about as much as the model assumes. Above 1 is overdispersion: the counts are more spread out than the model allows, usually because people differ in ways the predictors don't capture, and every standard error is too small by a factor of about √φ̂. Below 1 (underdispersion) is rare. There's no agreed cutoff; when φ̂ is clearly above 1, fit a negative binomial model and compare. Glossary →

Assumptions

1. Independent Observations The people (or other units) in your study don't influence each other.
What It Means
Each person's data comes from a separate, unconnected source, so knowing one person's scores tells you nothing about anyone else's. In designs that measure people more than once, their own scores are meant to be related; it's different people who must be independent. Methods built for clustered data, such as multilevel models and GEE, also allow members of the same group to be related, as long as the analysis names the grouping; then it's the groups that must be independent.
Applies To
Observations
How to Check
By design, not by a test. Did each person contribute their own data, with no clusters such as roommates, couples, therapy groups or classrooms whose members influence each other? If there are clusters, does the analysis model every one of them?
How Much It Matters
A lot. Events that come in clusters, such as several people from the same clinic or group, or repeated diaries from the same person treated as separate cases, make the standard errors too small and every p too small.

If It's Violated

2. Variance Equal to the Mean Your counts vary about as much as the Poisson model says they should.
What It Means
A Poisson model has no separate spread parameter: it assumes that among people with the same predicted count, the variance of the counts equals that count. Real counts usually vary more (overdispersion): people may differ in ways the predictors don't capture, one event may make the next more likely, or the model may be missing a predictor or a curve. If the model's form is right, the coefficients stay roughly right, but every standard error is too small, so intervals are too narrow and p-values too small.
Applies To
Outcome
How to Check
Fit the Poisson model and divide its Pearson χ² by its residual degrees of freedom (the dispersion statistic). Near 1 is fine; clearly above 1 means overdispersion. Compare the fitted Poisson and negative binomial models' log-likelihoods, and compare the number of zeros you observed with the number the model expects: many more zeros than expected points to excess zeros. Plot the Pearson residuals against the fitted counts: their spread should stay level, and a fan that widens with the fitted count means the variance grows faster than the mean, which favors negative binomial over quasi-Poisson (Ver Hoef & Boveng, 2007). A large dispersion statistic can also come from a missing predictor or curve, so check the model's form before blaming the variance.
How Much It Matters
A lot, for the standard errors. The coefficients stay roughly right, but each standard error is too small by about √φ̂, so intervals are too narrow and p too small, and the problem doesn't shrink with more people. In the panic diary study φ̂ = 1.42: the Poisson standard errors are about 16% too small, the treatment's p of .004 becomes .019 in the negative binomial model, and the interval widens. 23 people had no attacks, where the Poisson model expects 14 and the negative binomial model 19, a gap well within chance. A likelihood-ratio test of the Poisson model against the negative binomial tests θ at the edge of its range, so halve its p (here p = .0002).

If It's Violated

  • Alternative Test Fit a negative binomial model, which estimates the extra variation (θ) and reads the same way. The usual choice in psychology.
  • Correction Keep the Poisson coefficients and multiply every standard error by √φ̂ (a quasi-Poisson model). Ver Hoef and Boveng (2007) explain when it and the negative binomial disagree.
  • Robust Estimation Sandwich (robust) standard errors for the Poisson coefficients, which don't assume any particular variance.
  • Alternative Test If there are far more zeros than even a negative binomial model expects, because some people could never have the event, a zero-inflated or hurdle model separates the never-group from the rest (no page yet).
3. Linearity on the Log Scale Each continuous predictor relates to the log of the expected count or rate in a straight line.
What It Means
The model assumes that each one-unit rise in a continuous predictor multiplies the expected rate by the same factor (adds the same amount to its log), whether it happens low or high on the predictor's range. On the count scale that's a curve, so the counts themselves aren't expected to rise in a straight line, and a raw scatterplot of counts, bunched at 0 and fanning out as the mean grows, can't show it. A threshold, a leveling off or a U shape breaks it. Yes/no predictors meet it automatically.
Applies To
Relationship
How to Check
Plot the Pearson or deviance residuals against each continuous predictor and look for a curve or trend (a flat band is fine), or add a centered squared term or a spline for the predictor and compare the two models with a likelihood-ratio test.
How Much It Matters
Moderate to a lot, for what the rate ratios mean. The model assumes each predictor changes the log of the rate by a constant amount per unit, so each extra GAD-7 point multiplies the rate by the same 1.06. If anxiety only matters above some level, one rate ratio misdescribes everyone. In the panic diary study, a centered GAD-7² term in the negative binomial model shows no sign of a curve, χ²(1) = 0.44, p = .51.

If It's Violated

  • Transformation Log-transform a predictor whose effect levels off, or split it into a few bands.
  • Correction Add a centered squared term, or a spline, to let the effect on the log rate bend.

Limitations

This Is Not Your Test If…

What It Doesn't Tell You

  • Why. Randomization lets us read the CBT difference as an effect of being offered the program rather than waiting: everyone is analyzed in the group they were randomized to (intention to treat), whether or not they finished. It isn't the effect of completing CBT. Comparing only those who finished with the waitlist would lose randomization's protection, since people who finish a program differ from those who drop out. The anxiety rate ratio is only an association: anything that goes with baseline anxiety could be behind it.
  • Whether it matters. A rate 20% lower is about 0.6 fewer attacks a month for a typical waitlisted person. Whether that's worth the program is a clinical judgment, not a statistical one.
  • Who else. The rate ratios describe people like these 200, over diaries of three to eight weeks. They also assume each person's rate stays steady through the diary; attacks that bunch up early, or fade as the diary goes on, break that.
  • Whether there are two kinds of zeros. Some people may never have an attack in a given month while others are simply lucky. A Poisson or negative binomial model treats every zero alike. If many more people report zero than even the negative binomial expects (here 23 against its 19, a gap well within chance, so these data give no sign of it), a zero-inflated or hurdle model separates the two groups (Coxe et al., 2009); there's no page for it yet.

Common Mistakes

  • Don'tLog-transform the counts and run ordinary regression.

    InsteadFit a Poisson or negative binomial model to the raw counts. Logging needs a fudge for zeros and biases the estimates (O'Hara & Kotze, 2010).

  • Don'tCompare raw counts when people were observed for different lengths of time.

    InsteadInclude log(time observed) as an offset, so the model compares rates.

  • Don'tReport Poisson standard errors without checking the dispersion statistic.

    InsteadReport φ̂, and switch to a negative binomial or quasi-Poisson model when it's clearly above 1.

  • Don'tEnter log(time observed) as an ordinary predictor when you want rates.

    InsteadUse it as an offset, which fixes its coefficient at 1. As a predictor it answers a different question.

  • Don'tRead an incidence rate ratio of 0.70 as 0.70 fewer events.

    InsteadRead it as a rate 30% lower, and translate it into events per week or month for a typical person.

  • Don'tDescribe a regression coefficient as the effect of raising the predictor.

    InsteadSay the prediction is that much higher per unit of the predictor, on the model's scale (the outcome itself, or its log-odds or log rate in a logistic, ordinal or Poisson model). Causal language needs a design that supports it, such as an experiment.

  • Don'tRead a bigger χ² as a stronger effect.

    InsteadReport an effect size suited to the design, such as Cramér's V, an odds ratio or a model's coefficients. χ² grows with sample size even when the effect doesn't.

  • Don'tRead a p-value as the chance that the null hypothesis is true.

    InsteadRead it as: if the null were true, how often data at least this extreme would turn up.

Result Decoder

Significance Guidelines

Benchmarks Significance (p Value)
  • p < .001

    Very strong evidence against the null hypothesis.

  • .001 ≤ p < .01

    Strong evidence against the null hypothesis.

  • .01 ≤ p < .05

    Statistically significant at α = .05. Moderate evidence against the null hypothesis.

  • .05 ≤ p < .10

    Not significant. Weak evidence; sometimes called a “trend.” Treat it as inconclusive.

  • p ≥ .10

    Not significant. Little evidence against the null, which is not the same as evidence of no effect.

These evidence labels are one common convention for describing p in graded language rather than as significant or not; the boundaries are rough guides, not rules, and other authors draw them differently. Muff et al. (2022) argue for this kind of language of evidence over a single cutoff. α is chosen before the analysis, and p never says how big or important an effect is.

More about statistical significance ↓

How It's Reported

Equation Form

Example

negative binomial: IRR = 0.80, 95% CI [0.67, 0.96], z = −2.34, p = .019; θ = 8.38; model χ²(2) = 33.32, p < .001; Poisson: φ̂ = 1.42, IRR = 0.79, 95% CI [0.68, 0.93], z = −2.90, p = .004

Table 1. Poisson Regression Predicting the Rate of Panic Attacks

Example
bSEzpIRR95% CI (IRR)
Intercept −3.1050.126−24.63< .001
CBT (vs waitlist) −0.2310.079−2.90.0040.79[0.68, 0.93]
Baseline GAD-7 (per point) 0.0630.0096.63< .0011.06[1.05, 1.08]

Note. N = 200. Offset: ln(days of diary). IRR = incidence rate ratio, eb; the 95% CIs are Wald intervals, eb ± 1.96 SE. Likelihood-ratio χ²(2) = 53.25, p < .001. Dispersion φ̂ = 1.42, so these standard errors are too small; see Table 2. CBT is the coefficient worked through in the Napkin Math.

Table 2. Negative Binomial Regression, the Same Model Allowing Extra Variation

Example
bSEzpIRR95% CI (IRR)
Intercept −3.1220.149−20.94< .001
CBT (vs waitlist) −0.2220.095−2.34.0190.80[0.67, 0.96]
Baseline GAD-7 (per point) 0.0630.0115.54< .0011.07[1.04, 1.09]

Note. N = 200. Offset: ln(days of diary). θ = 8.38, estimated by maximum likelihood (variance = μ + μ²/θ). Likelihood-ratio χ²(2) = 33.32, p < .001, against a negative binomial model with no predictors. The coefficients barely move; the standard errors grow by about a fifth, which is what the dispersion statistic predicted.

Full Text (APA Style)

Example

Over diaries of 21 to 56 days, participants randomized to CBT recorded 2.17 panic attacks per 30 days and those randomized to the waitlist 2.78 (M = 3.23 attacks per person, SD = 2.58). A Poisson regression with log diary length as an offset and baseline GAD-7 as a covariate was mildly overdispersed (Pearson χ²/df = 1.42), so we report a negative binomial regression with the same terms (θ = 8.38). Holding baseline anxiety constant, assignment to CBT was associated with a 20% lower rate of panic attacks, IRR = 0.80, 95% CI [0.67, 0.96], z = −2.34, p = .019, and each point of baseline anxiety with a 7% higher rate, IRR = 1.07, 95% CI [1.04, 1.09]. The Poisson model gave a similar treatment estimate, IRR = 0.79, 95% CI [0.68, 0.93], but understated its uncertainty (p = .004). In the negative binomial model, treatment and baseline anxiety together improved on a model with neither, likelihood-ratio χ²(2) = 33.32, p < .001.

  • Say that the model has an offset and what the exposure is, so readers know the ratios compare rates, not raw counts.
  • Report incidence rate ratios with their intervals, not just the log-scale b's, and translate one into events per week or month.
  • Report the dispersion statistic. If it's clearly above 1, report the negative binomial (or quasi-Poisson) results as the main ones, as here, or at least alongside.
  • Name the interval type. Software gives Wald intervals (as here) or profile-likelihood intervals, which differ slightly.

Reading Each Piece

The equation form above, piece by piece, in this order: how big, how certain, then how surprising.

Negative binomial: IRR = 0.80, 95% CI [0.67, 0.96], p = .019The answer Allowing for the extra variation, CBT still goes with a rate about 20% lower, plausibly from 4% to 33% lower. If CBT did nothing, a result at least this far from 1, in either direction, would turn up about 2% of the time, so it's statistically significant, but less decisively than the Poisson model claimed. This is the result to lead with.
IRR = 0.79 Holding baseline anxiety constant, people assigned to CBT had panic attacks at 0.79 times the waitlist's rate, 21% less often (Poisson; the negative binomial's 0.80 in Step 9 is the one to report). For a typical waitlisted person with about 2.8 attacks a month, that's roughly 0.6 fewer a month. There are no conventional benchmarks for a rate ratio; whether that matters is a clinical question.
95% CI [0.68, 0.93] If the Poisson model's assumptions held, CBT's true rate ratio is plausibly anywhere from 0.68 (32% fewer attacks) to 0.93 (7% fewer). The whole range is below 1. Step 8 shows the interval is too narrow for these data.
z = −2.90, p = .004 CBT's coefficient is 2.90 standard errors below zero. If CBT made no difference to the rate, a result at least this far out, in either direction, would turn up about 4 times in 1,000, if the Poisson model's standard errors were right. They aren't quite (Step 8).
χ²(2) = 53.25, p < .001 Together, treatment and baseline anxiety fit the counts far better than a single shared rate for everyone. Most of the drop is anxiety's; adding CBT to a model that already has anxiety lowers the deviance by 8.50. The 2 is the degrees of freedom: the two predictors added. Overdispersion inflates this χ² too (Step 8); the same test in the negative binomial model gives χ²(2) = 33.32, still p < .001.
φ̂ = 1.42 The attack counts' variance is 42% larger than a Poisson model allows (their spread about 19% larger): mild overdispersion, most likely because people differ in how prone they are to attacks in ways CBT and baseline anxiety don't capture. The Poisson interval and p above are too optimistic.

Calculations

The Running Example

Example

The panic diary study: 200 adults, randomized to CBT or a waitlist, each with a count of panic attacks over their own diary period (21–56 days) and a baseline GAD-7 anxiety score.

Null Hypothesis:
\text{IRR}_{\text{CBT}} = 1
Alternative (Two-Sided):
\text{IRR}_{\text{CBT}} \neq 1

IRR (the incidence rate ratio) is the CBT group's rate of panic attacks divided by the waitlist's, among people with the same baseline anxiety, in the population these 200 represent. 1 means being assigned to CBT makes no difference to the rate. The whole model has its own null too, that neither CBT nor baseline anxiety predicts the rate, tested by the likelihood-ratio χ².

What Your Data Should Look Like

One row per person: a column with their count (a whole number, 0 or more), one for each predictor, and, if people were observed for different lengths of time, one for that exposure (days, hours, sessions). The model uses the log of the exposure as an offset. Don't divide the counts by the exposure yourself: the model needs the counts, not rates.

Example
participantcbtgad7daysattacks
118260
2110270
3114220
4121366
5114282
6012533
… 194 more rows
MSD12
1. Attacks3.232.58
2. CBT (1) vs Waitlist (0)0.500.50−.16
3. Baseline GAD-711.004.06.30−.02

N = 200.

Correlations are Pearson r.

Simulated data.

Download all 200 rows (data.csv)

Napkin Math

The whole test, one idea at a time. Under each heading, a line says why the step matters, and the green note reads the equation in plain English with the panic diary study's numbers. Collapse all to read just the headings: they tell the whole story on their own.

  1. 1. Put Everyone's Count on the Same Clock

    The idea: a Poisson distribution describes counts of events that happen independently at a steady rate. Someone who averages 3 attacks a month has some months with none and some with 6, and under a Poisson distribution the variance of those monthly counts equals their average, 3. Poisson regression predicts each person's expected count, μ, from the predictors, so it's really a model of rates. People kept their diaries for different numbers of days, and a longer diary catches more attacks. Dividing each group's attacks by its total days turns counts into rates we can compare.

    \text{rate} = \dfrac{\sum y_i}{\sum t_i} \times 30 =\dfrac{282}{3900} \times 30

    \text{CBT} = 2.17 per 30 days

    \text{Waitlist} = 2.78 per 30 days

    \text{Ratio} = 0.78

    \sum y_i is the CBT group's 282 attacks and \sum t_i their 3,900 diary days, so 282 ÷ 3,900 × 30 = 2.17 attacks per 30 days. The waitlist's 363 attacks over 3,917 days is 2.78, so CBT's rate is 2.17 ÷ 2.78 = 0.78 of the waitlist's, before allowing for anxiety.

  2. 2. Find the Equation That Predicts Each Person's Rate

    The model predicts the natural log (ln) of each person's expected count, μᵢ. The log undoes raising e to a power, so a predictor that adds to ln μ multiplies μ itself, and predictions can never go below zero. Log diary length enters with its coefficient fixed at 1 (the offset), which makes the other coefficients about the rate per day. Software climbs to the coefficients that make the observed counts most likely.

    \ln(\mu_i) = \ln(t_i) + b_0 + b_1 x_{1i} + b_2 x_{2i} =\ln(t_i) - 3.105 - 0.231\,x_1 + 0.063\,x_2

    b_0 = −3.105

    b_1 = −0.231

    b_2 = 0.063

    x_1 is 1 for CBT and 0 for the waitlist, x_2 is baseline GAD-7 and t_i is days of diary. There's no formula for the b's, so software finds them by iteratively reweighted least squares, a few rounds of weighted regression, and we take its numbers. b_1 = −0.231 says being assigned to CBT goes with a log rate 0.231 lower, holding anxiety constant.

  3. 3. Predict How Many Attacks One Person Should Have

    Undoing the log shows what the equation means for one person. The exponent gives their expected rate per day; multiplying by their days gives an expected count.

    \hat{\mu}_i = t_i \times e^{\,b_0 + b_1 x_{1i} + b_2 x_{2i}} =26 \times e^{\,-3.105 - 0.231 + 0.063 \times 8}

    = 1.53 attacks

    Participant 1 was assigned to CBT, had a baseline GAD-7 of 8 and kept the diary for t = 26 days. The exponent is −3.105 − 0.231 + 0.063 × 8 = −2.83, and e^{-2.83} = 0.059 attacks a day (using the unrounded b's). Times 26 days, \hat{\mu} = 1.53 attacks. They recorded none, which a count this low makes unsurprising.

  4. 4. Turn the Treatment Coefficient into a Rate Ratio

    A coefficient on the log scale is hard to picture. Raising e to it turns an added amount on the log scale into a multiplier on the rate.

    \text{IRR} = e^{\,b} =e^{-0.231}

    = 0.79

    \text{IRR}_{\text{GAD-7}} = 1.06

    b is CBT's −0.231, and e^{-0.231} = 0.79: holding anxiety constant, the CBT group's rate of attacks is 0.79 times the waitlist's, 21% lower. Baseline anxiety's e^{0.063} = 1.06 means each GAD-7 point goes with a 6% higher rate.

  5. 5. Give a Range of Plausible Rate Ratios

    The interval is built on the log scale, where the estimate's sampling distribution is close to normal, then turned back into rate ratios. That's why it isn't symmetric around 0.79.

    \left[\, e^{\,b - 1.96\,SE_b},\; e^{\,b + 1.96\,SE_b} \,\right] =e^{\,-0.231 \pm 1.96 \times 0.0794}

    = [0.68, 0.93]

    SE_b = 0.0794 is the standard error of CBT's coefficient. −0.231 ± 1.96 × 0.0794 runs from −0.386 to −0.075 on the log scale; e to those powers gives 0.68 to 0.93.

  6. 6. Ask How Surprising the Treatment Effect Would Be If CBT Did Nothing

    As in any regression, a coefficient is judged by how many standard errors it sits from zero. With maximum likelihood the yardstick, in large samples, is the normal curve rather than a t distribution.

    z = \dfrac{b_j}{SE_{b_j}} =\dfrac{-0.231}{0.0794}

    = −2.90

    p = .004

    b_j is CBT's −0.231 and SE_{b_j} its 0.0794, so z = −2.90: 2.90 standard errors below zero. Twice the chance of a standard normal landing below −2.90 is p = .004.

  7. 7. Ask Whether the Predictors Together Beat a Model with None

    Deviance measures how far a model falls short of fitting every count exactly, like the leftover sum of squares in ordinary regression. How much it drops when the predictors go in is the likelihood-ratio test.

    \chi^2_{LR} = D_{\text{null}} - D_{\text{model}}, \quad df = k =355.74 - 302.49

    = 53.25

    D_{\text{null}} = 355.74

    D_{\text{model}} = 302.49

    p < .001

    \chi^2_{\text{NB}} = 33.32

    p_{\text{NB}} < .001

    D_{\text{null}} = 355.74 is the deviance of a model with only an intercept (and the offset); D_{\text{model}} = 302.49 adds CBT and anxiety. The drop, 53.25, is judged against a χ² with k = 2 degrees of freedom, for p < .001.

  8. 8. Check Whether the Counts Vary More Than a Poisson Model Allows

    A Poisson model assumes the variance of a count equals its mean. Comparing each person's squared miss with that variance, and averaging, shows whether the counts really are that tidy. If they vary more, every standard error so far is too small.

    \hat{\phi} = \dfrac{1}{n - k - 1}\sum \dfrac{(y_i - \hat{\mu}_i)^2}{\hat{\mu}_i} =\dfrac{280.29}{197}

    = 1.42

    SE_{\text{quasi}} = 0.095

    \text{zeros} = 23

    \text{expected} = 13.66

    The Pearson χ², each (y_i - \hat{\mu}_i)^2 / \hat{\mu}_i added up, is 280.29, and n - k - 1 = 200 − 2 − 1 = 197. 280.29 ÷ 197 = 1.42: the counts' variance is 42% larger than the model allows. Their spread, and with it every standard error, is about √1.42 = 1.19 times what the model says, so CBT's should be about 0.095. There are 23 zeros where the model expects 13.66.

  9. 9. Refit, Letting People's Rates Vary More Than the Predictors Explain

    The negative binomial model keeps the same equation but adds a parameter, θ, for the extra variation: a count's variance becomes \mu + \mu^2/\theta instead of \mu, so the standard errors allow for the extra variation, provided it grows with the mean that way. If the conclusion survives, it doesn't depend on the Poisson model's tidy variance.

    \text{IRR} = e^{\,b} =e^{-0.222}

    = 0.80

    95\%\ \text{CI} = [0.67, 0.96]

    p = .019

    \theta = 8.38

    Maximum likelihood puts \theta at 8.38, so a person expected to have 3 attacks has a variance of 3 + 9 ÷ 8.38 = 4.07 rather than 3. Refit that way, CBT's coefficient is −0.222 and e^{-0.222} = 0.80, 95% CI [0.67, 0.96], p = .019: nearly the same effect, a wider interval and a larger p, still significant.

Run It in Software

The code below reads the running example's data from data.csv; save it in the same folder you run the code from.

SPSS (version 30, as documented by IBM)
  1. Transform → Compute Variable: make log_days = LN(days).
  2. Analyze → Generalized Linear Models → Generalized Linear Models.
  3. Type of Model: choose Poisson loglinear.
  4. Response: move attacks into Dependent Variable.
  5. Predictors: move cbt and gad7 into Covariates, and log_days into Offset Variable.
  6. Model: move cbt and gad7 into Model (Main effects).
  7. Statistics: check Include exponential parameter estimates, and click OK.
  8. For the negative binomial refit, go back to Type of Model and choose Custom: set Distribution to Negative binomial and Link function to Log, and under Parameter choose Estimate value. (The ready-made Negative binomial with log link fixes the parameter at 1.)

Or paste this syntax:

COMPUTE log_days = LN(days).
EXECUTE.

GENLIN attacks WITH cbt gad7
  /MODEL cbt gad7 INTERCEPT=YES OFFSET=log_days
    DISTRIBUTION=POISSON LINK=LOG
  /CRITERIA METHOD=FISHER(1) SCALE=1 COVB=MODEL
  /PRINT CPS DESCRIPTIVES MODELINFO FIT SUMMARY SOLUTION (EXPONENTIATED).

* Negative binomial refit, estimating the ancillary parameter (1/theta).
GENLIN attacks WITH cbt gad7
  /MODEL cbt gad7 INTERCEPT=YES OFFSET=log_days
    DISTRIBUTION=NEGBIN(MLE) LINK=LOG
  /PRINT FIT SUMMARY SOLUTION (EXPONENTIATED).

Reading the Output

Read Exp(B), the incidence rate ratio, with its 95% Wald Confidence Interval for Exp(B), then Sig. Check Value/df for the Pearson Chi-Square in Goodness of Fit before trusting any of it. SPSS tests each coefficient with a Wald χ², which is z² (2.905² = 8.44).

Goodness of Fit

Value df Value/df
Deviance 302.490 197 1.535
Pearson Chi-Square 280.287 197 1.423

Omnibus Test

Likelihood Ratio Chi-Square df Sig.
53.248 2 <.001

Parameter Estimates

Parameter B Std. Error Wald Chi-Square df Sig. Exp(B) 95% Wald CI for Exp(B) Lower Upper
(Intercept) -3.105 .1261 606.611 1 <.001 .045 .035 .057
cbt -.231 .0794 8.438 1 .004 .794 .680 .928
gad7 .063 .0095 43.901 1 <.001 1.065 1.045 1.085
(Scale) 1

Pearson χ²/df = 1.42, so run the negative binomial syntax too. SPSS estimates the negative binomial's ancillary parameter, 1/θ, as about .119 (R's θ = 8.38), and its Exp(B) for cbt is about .80; SPSS's standard errors for it can differ from R's in the third decimal, because SPSS estimates the ancillary parameter and the coefficients together.

Output laid out from this page's computed results, not captured from a live run; labels and layout may differ slightly in your version.

R
library(MASS)
panic <- read.csv("data.csv")
fit <- glm(attacks ~ cbt + gad7 + offset(log(days)), family = poisson, data = panic)
summary(fit)
round(exp(cbind(IRR = coef(fit), confint.default(fit))), 3)

# Likelihood-ratio test of the whole model
null <- glm(attacks ~ 1 + offset(log(days)), family = poisson, data = panic)
anova(null, fit, test = "Chisq")

# Overdispersion: Pearson chi-square / residual df
sum(residuals(fit, type = "pearson")^2) / df.residual(fit)

# Negative binomial refit
nb <- glm.nb(attacks ~ cbt + gad7 + offset(log(days)), data = panic)
summary(nb)
round(exp(cbind(IRR = coef(nb), confint.default(nb))), 3)
Call:
glm(formula = attacks ~ cbt + gad7 + offset(log(days)), family = poisson,
    data = panic)

Coefficients:
             Estimate Std. Error z value Pr(>|z|)
(Intercept) -3.104982   0.126068 -24.629  < 2e-16 ***
cbt         -0.230693   0.079417  -2.905  0.00367 **
gad7         0.062823   0.009482   6.626 3.45e-11 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for poisson family taken to be 1)

    Null deviance: 355.74  on 199  degrees of freedom
Residual deviance: 302.49  on 197  degrees of freedom
AIC: 836.98

Number of Fisher Scoring iterations: 5

              IRR 2.5 % 97.5 %
(Intercept) 0.045 0.035  0.057
cbt         0.794 0.680  0.928
gad7        1.065 1.045  1.085
Analysis of Deviance Table

Model 1: attacks ~ 1 + offset(log(days))
Model 2: attacks ~ cbt + gad7 + offset(log(days))
  Resid. Df Resid. Dev Df Deviance  Pr(>Chi)
1       199     355.74
2       197     302.49  2   53.248 2.737e-12 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
[1] 1.422778

Call:
glm.nb(formula = attacks ~ cbt + gad7 + offset(log(days)), data = panic,
    init.theta = 8.376933413, link = log)

Coefficients:
            Estimate Std. Error z value Pr(>|z|)
(Intercept) -3.12170    0.14907 -20.942  < 2e-16 ***
cbt         -0.22166    0.09456  -2.344   0.0191 *
gad7         0.06349    0.01145   5.545 2.94e-08 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for Negative Binomial(8.3769) family taken to be 1)

    Null deviance: 265.64  on 199  degrees of freedom
Residual deviance: 228.52  on 197  degrees of freedom
AIC: 826.7

Number of Fisher Scoring iterations: 1


              Theta:  8.38
          Std. Err.:  3.16

 2 x log-likelihood:  -818.702
              IRR 2.5 % 97.5 %
(Intercept) 0.044 0.033  0.059
cbt         0.801 0.666  0.964
gad7        1.066 1.042  1.090

Read the IRR table first (the exponentiated coefficients and their Wald intervals), then z value and Pr(>|z|). The deviance lines give the likelihood-ratio test; the [1] 1.422778 is the dispersion statistic. confint(fit) without .default gives profile-likelihood intervals instead, which differ slightly. The negative binomial's null deviance is on its own scale, so don't compare it with the Poisson one; compare the two models' AIC (836.98 vs 826.7) or log-likelihoods instead.

Output from a real run: R 4.6.1 (MASS 7.3-65), run 2026-09-29.

Python
import numpy as np
import pandas as pd
import statsmodels.api as sm
import statsmodels.formula.api as smf

panic = pd.read_csv("data.csv")
fit = smf.glm("attacks ~ cbt + gad7", data=panic, family=sm.families.Poisson(),
              offset=np.log(panic["days"])).fit()
print(fit.summary())
print(np.exp(fit.conf_int()).assign(IRR=np.exp(fit.params)).round(3))
print(fit.pearson_chi2 / fit.df_resid)

nb = smf.negativebinomial("attacks ~ cbt + gad7", data=panic,
                          offset=np.log(panic["days"])).fit()
print(nb.summary())
                 Generalized Linear Model Regression Results
==============================================================================
Dep. Variable:                attacks   No. Observations:                  200
Model:                            GLM   Df Residuals:                      197
Model Family:                 Poisson   Df Model:                            2
Link Function:                    Log   Scale:                          1.0000
Method:                          IRLS   Log-Likelihood:                -415.49
Date:                Tue, 29 Sep 2026   Deviance:                       302.49
Time:                        13:32:48   Pearson chi2:                     280.
No. Iterations:                     4   Pseudo R-squ. (CS):             0.2337
Covariance Type:            nonrobust
==============================================================================
                 coef    std err          z      P>|z|      [0.025      0.975]
------------------------------------------------------------------------------
Intercept     -3.1050      0.126    -24.629      0.000      -3.352      -2.858
cbt           -0.2307      0.079     -2.905      0.004      -0.386      -0.075
gad7           0.0628      0.009      6.626      0.000       0.044       0.081
==============================================================================
               0      1    IRR
Intercept  0.035  0.057  0.045
cbt        0.680  0.928  0.794
gad7       1.045  1.085  1.065
1.4227777539086204
Optimization terminated successfully.
         Current function value: 2.046756
         Iterations: 13
         Function evaluations: 17
         Gradient evaluations: 17
                     NegativeBinomial Regression Results
==============================================================================
Dep. Variable:                attacks   No. Observations:                  200
Model:               NegativeBinomial   Df Residuals:                      197
Method:                           MLE   Df Model:                            2
Date:                Tue, 29 Sep 2026   Pseudo R-squ.:                 0.03910
Time:                        13:32:48   Log-Likelihood:                -409.35
converged:                       True   LL-Null:                       -426.01
Covariance Type:            nonrobust   LLR p-value:                 5.825e-08
==============================================================================
                 coef    std err          z      P>|z|      [0.025      0.975]
------------------------------------------------------------------------------
Intercept     -3.1217      0.150    -20.841      0.000      -3.415      -2.828
cbt           -0.2217      0.094     -2.346      0.019      -0.407      -0.037
gad7           0.0635      0.012      5.509      0.000       0.041       0.086
alpha          0.1194      0.045      2.648      0.008       0.031       0.208
==============================================================================

The coefficient table is on the log scale; exponentiate to get rate ratios, as the second block does. Deviance and Pearson chi2 are in the header. statsmodels' negative binomial reports α = 1/θ (0.1194 = 1/8.38) and estimates it together with the coefficients, so its standard errors (0.094 for cbt) differ from R's glm.nb (0.0946) in the third decimal; the page follows R. Its LLR p-value tests the predictors against a negative binomial model with none.

Output from a real run: Python 3.9.6 (statsmodels 0.14.6), run 2026-09-29.

Impact

Effect Size

In a count model the effect size is the incidence rate ratio: e raised to a coefficient, which says how many times higher the rate of events is per unit of the predictor, holding the others constant. It's a ratio, so studies that counted over different periods can be compared if each person's rate is steady over time, though not studies of different kinds of events. For a continuous predictor, check that both studies use the same unit (per point, or per SD).

\text{IRR} = e^{\,b} = e^{-0.222} = 0.80

People randomized to CBT had panic attacks at 0.80 times the waitlist's rate in the negative binomial model that allows for the extra variation, and 0.79 times in the Poisson model: about 20% less often, or about 0.6 fewer attacks a month for a typical waitlisted person.

A continuous predictor's ratio is per unit, and it compounds: each GAD-7 point multiplies the rate by 1.065 in the Poisson model, so 5 points multiply it by 1.065⁵ ≈ 1.37, not the 1.33 you get by adding 5 × 0.065 to 1.

There are no conventional small, medium and large cutoffs for an incidence rate ratio, and borrowing the odds ratio's would mislead. Judge it in the outcome's own terms: how many fewer events, for whom, over what period. As with an odds ratio, 0.50 and 2.00 are equally strong in opposite directions.

Statistical Significance

A Poisson regression answers two kinds of significance question. Each coefficient's z test asks: if this predictor were unrelated to the rate once the others were known, how often would its coefficient land at least this far from zero, in either direction? The likelihood-ratio χ² asks the same about all the predictors together, from how much the deviance drops when they go in. Both lean on the Poisson model's assumption that the counts vary no more than their mean: the z test through its standard errors, the χ² through the likelihood.

If CBT made no difference to the rate, a coefficient at least this far from zero, in either direction, would turn up about 0.4% of the time by the Poisson model (p = .004). But the counts' variance is 42% larger than that model allows, so its standard errors are too small; the negative binomial model, which allows for the extra variation, puts p at .019. Still significant, but less decisively. Together, CBT and baseline anxiety beat a model with neither, χ²(2) = 33.32, p < .001, in the negative binomial model (53.25 in the Poisson model, inflated the same way).

p is not the probability that CBT does nothing, and it says nothing about how many attacks CBT prevents; that's what the rate ratio and its interval are for. A Poisson p computed without checking for overdispersion can be far too small.

Benchmarks p at a Glance
Very Strong
p < .001
Strong
.001 ≤ p < .01
Significant
.01 ≤ p < .05
Weak
.05 ≤ p < .10
Little
p ≥ .10

Strength of evidence against the null hypothesis.

Full benchmarks and caveats ↑

Power and Sample Size

Power is the chance that a study detects an effect that's really there. A real effect doesn't guarantee a significant result: with too few people, it can easily come out non-significant. For a count model, power depends on the rate ratio you expect, how many events people have (their base rate times how long they're observed), how many people you have, α, and how much more the counts vary than a Poisson model allows. Overdispersion costs power just as it inflates standard errors, so a plan that assumes Poisson variation promises more than it will deliver.

Suppose we expected CBT to lower the rate of attacks by about 20% (IRR = 0.80), with waitlisted people of average anxiety having 0.09 attacks a day, diaries of 21 to 56 days, baseline GAD-7 in the model as a covariate, and counts overdispersed as ours were (θ ≈ 8.4). We simulated 4,000 studies at each sample size and fitted a negative binomial model to each, estimating θ afresh each time. 200 people, 100 per group, would detect the effect about 65% of the time; about 280 would give 80%. Had the counts truly had Poisson variation, 200 would have given about 80%, which is what a plan that ignores overdispersion promises. A larger effect is easier: a 30% reduction (IRR = 0.70) needs about 115 people for 80%.

Researchers usually aim for at least 80% power. Plan it before collecting data, for the smallest rate ratio that would matter, and allow for overdispersion. There's no calculator for count models on this page. Simulation, as above, is the most flexible way to plan: generate data like the study you expect many times, fit the model you'll use, and count how often it comes out significant. Longer observation raises power too, since each person contributes more events, but with overdispersion the gain levels off: differences between people don't average out over time, only over more people.

Other Info

Follow-Up and Post Hoc Tests

None needed. This test answers a single question, so there are no planned contrasts or post hoc tests to run.

References

  • Nelder, J. A., & Wedderburn, R. W. M. (1972). Generalized linear models. Journal of the Royal Statistical Society. Series A (General), 135(3), 370–384. https://doi.org/10.2307/2344614
    “These generalized linear models are illustrated by examples relating to four distributions; the Normal, Binomial (probit analysis, etc.), Poisson (contingency tables) and gamma (variance components).” (Abstract)
  • Coxe, S., West, S. G., & Aiken, L. S. (2009). The analysis of count data: A gentle introduction to Poisson regression and its alternatives. Journal of Personality Assessment, 91(2), 121–136. https://doi.org/10.1080/00223890802634175
    “We also discuss the problems of excess zeros in which a subgroup of respondents who would never display the behavior are included in the sample and truncated zeros in which respondents who have a zero count are excluded by the sampling plan.” (Abstract)
    “In cases in which the outcome variable is a count with a low arithmetic mean (typically < 10), standard ordinary least squares regression may produce biased results.” (Abstract)
    “Two variants of Poisson regression, overdispersed Poisson regression and negative binomial regression, are introduced that may provide more optimal results when a key assumption of standard Poisson regression is violated.” (Abstract)
  • Gardner, W., Mulvey, E. P., & Shaw, E. C. (1995). Regression analyses of counts and rates: Poisson, overdispersed Poisson, and negative binomial models. Psychological Bulletin, 118(3), 392–404. https://doi.org/10.1037/0033-2909.118.3.392
    “In the overdispersed model, a factor can be estimated that corrects the regression model's inferential statistics.” (Abstract)
    “The simplest, the Poisson regression model, is likely to be misleading unless restrictive assumptions are met because individual counts are usually more variable ("overdispersed") than is implied by the model.” (Abstract)
  • O'Hara, R. B., & Kotze, D. J. (2010). Do not log-transform count data. Methods in Ecology and Evolution, 1(2), 118–122. https://doi.org/10.1111/j.2041-210X.2010.00021.x
    “We found that the transformations performed poorly, except when the dispersion was small and the mean counts were large.” (Summary, point 4)
  • Ver Hoef, J. M., & Boveng, P. L. (2007). Quasi-Poisson vs. negative binomial regression: How should we model overdispersed count data? Ecology, 88(11), 2766–2772. https://doi.org/10.1890/07-0043.1
    “While they often give similar results, there can be striking differences in estimating the effects of covariates.” (Abstract)
  • Muff, S., Nilsen, E. B., O'Hara, R. B., & Nater, C. R. (2022). Rewriting results sections in the language of evidence. Trends in Ecology & Evolution, 37(3), 203–210. https://doi.org/10.1016/j.tree.2021.10.009
    “we suggest language of evidence that allows for a more nuanced approach to communicate scientific findings as a simple and intuitive alternative to statistical significance testing” (Abstract)

Spot an error? Something still confusing? Tell us at support@statsforpsychology.com.