Poisson regression: modelling how many times customers do things

The customer who ordered seven times

Every December the same thing happens at The Gift Bow, my Solidus demo hamper shop. Someone in the team notices a customer who has ordered again and again, gets excited, and wants to build a VIP programme around people like her. Last time it was a woman in Salzburg who had placed seven orders in two years. Her first hamper was a £75 gift for a colleague, she opens most of the emails we send, and she has never bought anything for herself. Is she remarkable? Is she the tip of a segment we should chase? Or is she simply what happens when you have a few thousand customers and count for long enough?

That is a counting question, and counting questions have their own kind of model. Orders per customer per year, support tickets per account per month, hotel stays per guest per season, quotes raised per Werkbank workspace per week: none of these is a yes or no, and none of them is a smooth continuous quantity like revenue. They are whole numbers, they start at zero, they cannot go below it, and most of the mass sits at the bottom with a long thin tail to the right. Treating them like revenue and pushing them through ordinary least squares is one of the most common quiet errors I see in marketing analytics. This chapter is about the right tool for the job: Poisson regression, its more forgiving cousin the negative binomial, and what to do when your zeros come from two different places. The numbers throughout are illustrative, taken from a demo shop and rounded so the reasoning is easy to follow.

Where this sits in the series

We are in part two of the series, the chapters on dependent variable techniques, where a single outcome is explained by a set of drivers. The previous chapter, Polynomial distributed lags: how long does marketing keep working?, dealt with a continuous outcome spread over time. This one moves to outcomes that are counts. The next chapter, Logistic regression, lift charts and market basket analysis, takes the other common non continuous case, an outcome that is simply yes or no. Read the three together and you have most of the regression toolbox a small business will ever need.

Why ordinary regression stumbles on counts

Ordinary least squares assumes the outcome can wander freely in both directions around the line, with roughly the same spread everywhere. Counts break both assumptions.

First, the floor. A straight line fitted to orders per year will happily predict minus 0.3 orders for a low value gift buyer in a region with weak engagement. There is no such customer. You can clip the prediction at zero, but then your model is no longer the model you fitted, and the coefficients were estimated under the assumption that minus 0.3 was a legitimate place to be.

Second, the spread. For counts, the variance grows with the mean. Customers who order once a year on average scatter between zero and three. Corporate accounts who order five times a year scatter between one and twelve. OLS assumes constant variance, so it gives the noisy high frequency customers the same weight as the quiet ones, and its standard errors are wrong, usually too optimistic. You end up confident about effects that are not there.

Third, the shape. Counts at low means are skewed. The average customer at The Gift Bow orders 1.4 times a year, but the most common count is zero and a handful of accounts order fifteen times. A symmetric error term does not describe this. I argued in Counting and summing beats machine learning that most business questions are answered by counting well before they need anything clever. Poisson regression is what counting well looks like once you want to know why the counts differ.

The Poisson distribution and the log link

The Poisson distribution describes the number of events in a fixed window when events arrive independently at a steady average rate. If a customer's underlying rate is $\lambda$ orders per year, the probability of seeing exactly $k$ orders is

$$P(Y = k) = \frac{\lambda^{k} e^{-\lambda}}{k!}, \qquad k = 0, 1, 2, \ldots$$

Here $\lambda$ (lambda) is the expected number of events in the window, $k$ is the count you actually observe, $e$ is Euler's number and $k!$ is k factorial. The distribution has one parameter, and that parameter is both its mean and its variance. Remember that, because it is the property that will bite us later.

Poisson regression lets the rate depend on the customer. We do not model the rate directly, because a linear combination of drivers can be negative and a rate cannot. Instead we model its logarithm:

$$\log(\lambda_i) = \beta_0 + \beta_1 x_{1i} + \beta_2 x_{2i} + \cdots + \beta_p x_{pi}$$

$\lambda_i$ is the expected count for customer $i$, the $x$ values are that customer's characteristics (first order value, whether the first purchase was a gift, whether they open emails, region) and the $\beta$ values are the coefficients we estimate. This is the log link. Whatever the right hand side adds up to, the exponential of it is positive, so the predicted rate is always a legitimate rate.

The log link also changes how you read coefficients, and this is the part people either love or resist. In ordinary regression a coefficient is an amount added. Here it is a multiplier. Exponentiate the coefficient and you get the incidence rate ratio:

$$\text{IRR}_j = e^{\beta_j} = \frac{\lambda(x_j + 1)}{\lambda(x_j)}$$

The IRR for driver $j$ is the factor by which the expected rate changes when $x_j$ goes up by one unit, holding everything else fixed. An IRR of 1.70 means 70 percent more orders per year. An IRR of 0.70 means 30 percent fewer. Multiplicative effects are how most marketing actually behaves: a good email programme does not add 0.4 orders to everyone, it lifts everyone's rate by a proportion, so the corporate account gains more absolute orders than the occasional gift buyer. The model says so without being told.

Exposure: not everyone has been a customer for the same time

There is a trap waiting in the data before you fit anything. A customer who joined 23 months ago has had 23 months to place orders; one who joined 4 months ago has had 4. If you simply count orders per customer, the old customers look like heavy buyers and the model learns that tenure causes loyalty, which is backwards.

The fix is an exposure offset. You model the rate per unit of time and include the log of each customer's observed time as a term whose coefficient is fixed at one:

$$\log(\lambda_i) = \log(t_i) + \beta_0 + \beta_1 x_{1i} + \cdots + \beta_p x_{pi}$$

$t_i$ is the exposure, here the number of years customer $i$ has been observable. Moving $\log(t_i)$ to the left gives $\log(\lambda_i / t_i)$, the log of the yearly rate, which is what we actually want to explain. Every count model I build for a client has an offset, and I check the offset before I check anything else. For Seeblick Hotels the exposure is the number of seasons a guest has been in the database; for Werkbank it is the number of months a workspace has been live.

The Gift Bow model: reading incidence rate ratios

Here is what a fitted model looks like on 4,180 customers from the demo shop, observed for between six months and two years. The outcome is orders, the offset is log years observed, and the drivers are the ones the team argues about at Christmas. Austrian orders are converted to sterling before modelling so the first order value sits on one scale.

Driver	Coefficient	Standard error	IRR	95 percent interval	Plain reading
Intercept	minus 0.15	0.08	0.86	0.74 to 1.01	The baseline: a UK self buyer, £40 first order, no email opens, not corporate, orders 0.86 times a year
First order value, per £10 above £40	0.045	0.009	1.046	1.03 to 1.06	Every extra £10 in the first basket lifts the yearly rate by about 4.6 percent
First purchase was a gift	minus 0.36	0.06	0.70	0.62 to 0.79	Gift buyers reorder at 70 percent of the self buyer rate
Opened an email in the first 90 days	0.53	0.06	1.70	1.51 to 1.91	Early email engagement goes with a 70 percent higher rate
Austria rather than UK	minus 0.15	0.07	0.86	0.75 to 0.99	Austrian customers order 14 percent less often, and the interval nearly touches 1
Corporate account	0.88	0.10	2.41	1.98 to 2.93	Corporate accounts order at roughly two and a half times the rate

Read it from the top. The intercept is the rate for a customer with every driver at its reference value; exponentiated it is 0.86 orders a year, the anchor everything else multiplies. The first order value coefficient is per £10, so a £75 first order is 3.5 steps above the £40 reference and contributes $1.046^{3.5} \approx 1.17$, a 17 percent higher rate. The gift coefficient is negative, so its IRR is below one: gift buyers come back less, which every hamper seller suspects and few have quantified. Email engagement is the largest effect we can influence directly, at 1.70. The Austria effect is real but small and its interval runs up to 0.99, so I would not build a country strategy on it. The corporate effect is large and expected.

Now the Salzburg customer. She is an Austrian gift buyer with a £75 first order who opens emails and is not corporate. Multiply it through:

$$\lambda = 0.86 \times 1.046^{3.5} \times 0.70 \times 1.70 \times 0.86 \approx 1.03 \text{ orders per year}$$

Over the two years we have watched her, the model expects about 2.06 orders. She placed seven, roughly three and a half times the expected rate. Whether that is remarkable depends entirely on the next section.

Does the Poisson fit? The overdispersion check

The Poisson distribution insists that variance equals mean. Real customers rarely oblige. Some are quietly loyal, some bought once by accident, some are small businesses hiding behind a personal email address. That hidden heterogeneity makes the observed counts spread wider than a Poisson with the same mean would allow. The condition is called overdispersion, and it is the rule rather than the exception in customer data.

The first check is a picture. Here is the observed distribution of yearly order rates against what a Poisson with the same overall mean of 1.4 would predict:

Bars are the 4,180 observed customers, the line is the Poisson expectation with the same mean. Too many zeros, too few ones, and a tail the Poisson says should not exist.

The pattern is the classic signature: more zeros than expected, fewer customers in the middle, and a fat right tail. The Poisson predicts two customers with seven orders and essentially none with eight or more. The shop has 80 of them.

The second check is a number. After fitting the model, compute the Pearson dispersion statistic:

$$\hat{\phi} = \frac{1}{n - p} \sum_{i=1}^{n} \frac{(y_i - \hat{\lambda}_i)^2}{\hat{\lambda}_i}$$

$y_i$ is the observed count for customer $i$, $\hat{\lambda}_i$ the model's predicted count, $n$ the number of customers and $p$ the number of estimated coefficients. Each term asks how far the customer sits from the prediction, in units of the variance the Poisson would allow. If the Poisson is right, $\hat{\phi}$ is close to 1. For the raw counts above it is about 2.1: the variance is twice the mean. After the covariates, which soak up some of the spread (the corporate flag alone explains a lot of the tail), it is still 1.9. Anything above about 1.3 makes me stop trusting the Poisson standard errors, and above 1.5 I switch model.

What goes wrong if you do not? The coefficients are usually still roughly right; Poisson regression is quite robust for the point estimates. The standard errors are not. They are too small by a factor of roughly $\sqrt{\hat{\phi}}$, so with $\hat{\phi} = 1.9$ every interval in the table above should be about 38 percent wider. The Austria effect, whose interval already touched 0.99, would stop being significant. That is exactly the kind of effect people build campaigns around and then wonder why the campaign did nothing.

Negative binomial and zero inflation: choosing the model

The negative binomial model is the standard remedy. It keeps the log link and the IRR interpretation, so nothing you learned above is wasted, but it lets the variance grow faster than the mean:

$$\operatorname{Var}(Y_i) = \lambda_i + \alpha \lambda_i^{2}$$

$\alpha$ (alpha) is a dispersion parameter estimated from the data. When $\alpha = 0$ you are back at the Poisson. For The Gift Bow the fitted $\alpha$ is about 0.48 once the drivers are in the model. Refit with the negative binomial and the IRRs barely move, but the intervals widen honestly and the tail probabilities change dramatically.

Back to the Salzburg customer. Under the Poisson with an expected 2.06 orders, the probability of seven or more is about 0.5 percent, one customer in 190. She looks like a star. Under the negative binomial with $\alpha = 0.48$, the same probability is about 3.7 percent, one in 27. Among 4,180 customers we should expect roughly 150 people like her purely from ordinary variation in loyalty. She is a good customer, she deserves a thank you, but she is not a segment, and a VIP programme designed around her would be designed around noise. This is the practical difference between the two models: the negative binomial tells you how surprised to be.

Zero inflation is the other structural issue, and it is a different one. Overdispersion is about spread. Zero inflation is about the zeros coming from two populations. At a hamper shop, some zeros are customers who might well order again and simply have not yet; they sit on the Poisson or negative binomial curve. Other zeros are people who bought one gift for a leaving colleague, will never think of us again, and could not be moved by any email. They are structural zeros. A zero inflated model fits two parts: a logistic part for the probability of being a structural never again customer, and a count part for everyone else. It is worth it when your zero share is far above what even the negative binomial predicts and when you have a plausible driver for the never again group, such as a corporate gift shipped to a third party address.

Here is how I decide:

The decision path I follow for every count outcome. Most client work ends at the negative binomial.

For the demo shop the negative binomial predicts about 1,610 zeros against 1,560 observed. That is close enough. The zeros are explained by heterogeneity, not by a hidden never again population, so I stop there. At Seeblick Hotels the story was different: guests who booked once through a platform and never gave an email address were almost all structural zeros for direct bookings, and the zero inflated model was clearly better.

From model to targeting list

A fitted count model is only useful if it changes who you talk to. Three lists come out of it almost for free.

The first is the expected rate list. Every customer gets a predicted $\hat{\lambda}_i$ per year. Multiply by average contribution per order and you have expected yearly contribution, which is the honest basis for how much to spend keeping each person. This connects directly to the retention maths: a small lift in the rate for the middle of the list is worth more than a large lift for a handful at the top.

The second is the surprise list. Compare what each customer did with what the model expected. The Pearson residual $(y_i - \hat{\lambda}_i)/\sqrt{\hat{\lambda}_i + \alpha \hat{\lambda}_i^2}$ standardises the gap. Customers far above expectation are worth a human look, not a campaign: a good share of them turn out to be small companies ordering under a personal address, which means a phone call and a corporate account rather than a discount code.

The third is the lever list. Because email engagement carries an IRR of 1.70, the customers where an engagement nudge is worth most are those with a high baseline rate and no engagement yet. The model gives you the expected gain for each person directly: $\hat{\lambda}_i \times (1.70 - 1)$ additional orders per year if you could move them, which you cannot fully, so I discount it.

Customer	Profile	Expected per year	Observed per year	Standardised residual	Action
A	UK corporate, £120 first order, engaged	5.1	5.5	0.1	Keep; account manager already in place
B	Salzburg gift buyer, £75, engaged	1.03	3.5	2.0	Thank you note; check for a hidden business
C	UK self buyer, £95, not engaged	1.10	1.0	minus 0.1	Engagement nudge: highest expected gain
D	UK gift buyer, £35, not engaged	0.58	0	minus 0.7	Leave alone; likely structural zero
E	Vienna self buyer, £60, engaged	1.38	0	minus 0.9	Win back email; expected but silent

Five rows are enough to show the logic. Customer B is the Salzburg buyer: a residual of 2.0 is notable but not extraordinary, which is what the tail probability told us. Customer C is the one most marketers would ignore and the model would prioritise: a decent baseline with the biggest lever untouched. Customer D would get a discount from a naive campaign and it would be wasted.

Running it yourself

You need one row per customer with the count, the exposure in years and the drivers as they were known at the start of the window. That last point matters: use the first order value and the first 90 day email behaviour, not lifetime values, or you leak the outcome into the inputs. In SQL it is a single aggregation:

select c.id,
       count(o.id) as orders,
       extract(epoch from now() - c.first_order_at) / 31557600.0 as years_observed,
       (c.first_order_total_gbp - 40) / 10.0 as first_value_10,
       c.first_order_is_gift::int as gift,
       c.opened_email_first_90d::int as email_open,
       (c.country = 'AT')::int as austria,
       c.is_corporate::int as corporate
from customers c
left join orders o on o.customer_id = c.id and o.completed_at > c.first_order_at
where c.first_order_at < now() - interval '6 months'
group by c.id;

Fitting takes seconds in any statistics package. In Python with statsmodels:

import numpy as np
import statsmodels.api as sm
import statsmodels.formula.api as smf

formula = "orders ~ first_value_10 + gift + email_open + austria + corporate"
pois = smf.glm(formula, data=df, family=sm.families.Poisson(),
               offset=np.log(df["years_observed"])).fit()
phi = pois.pearson_chi2 / pois.df_resid
print(f"dispersion {phi:.2f}")          # near 1.0 is fine, above 1.5 switch

nb = smf.glm(formula, data=df, family=sm.families.NegativeBinomial(alpha=0.48),
             offset=np.log(df["years_observed"])).fit()
print(np.exp(nb.params))                 # incidence rate ratios

In R the same is glm(..., family = poisson, offset = log(years_observed)) followed by MASS::glm.nb. Estimating alpha properly rather than fixing it takes one more line in either language.

Before you trust it, run four checks. Plot observed against predicted counts by decile of predicted rate; the points should sit on the diagonal. Compare the predicted share of zeros with the observed share. Hold out the most recent six months and see whether the model's ranking of customers holds up on orders it never saw. And look at the largest residuals by hand; if the top ten are all one kind of customer, you are missing a driver, not looking at noise.

Time budget: a day for the data, half a day for the modelling, and a day to turn it into lists people will actually use. The modelling is the smallest part, which is true of nearly everything in this series.

Pitfalls

Counting the wrong window. If the outcome window overlaps with the period the drivers were measured in, the model explains the past with the past. Freeze the drivers at the start and count from there.

Forgetting the offset. I have seen a model that proudly showed tenure as the strongest driver of orders. Of course it was: the customers had been counted for longer. The offset removes this in one line and the whole story changes.

Reading IRRs as additive. An IRR of 1.70 on top of a baseline of 0.3 orders is 0.21 extra orders a year. The same IRR on a baseline of 5 is 3.5 extra orders. Always multiply through to the customer's own baseline before promising anyone a result.

Trusting Poisson intervals when the dispersion is 2. The coefficients look fine, the significance stars look fine, and half the significant effects are not significant. Check the dispersion every single time; it is one number.

Chasing the tail. Heavy buyers exist in every negative binomial, by construction. Before you announce a segment, ask how many customers at that level the model already expects. Usually it expects about as many as you found.

Adding zero inflation because it fits better. A zero inflated model will almost always fit slightly better, because it has more parameters. Use it when you can name the never again population and act on it differently; otherwise it is complexity without a decision attached.

How I do this for clients

The deliverable is a scored customer list and a one page decision memo, not a model file. It starts with the free workshop, where we agree what to count (orders, stays, tickets, quotes), what the exposure is, and the four or five drivers that were actually known at the start of each customer's window. Then I take two weeks on my own account before you commit to anything.

The data I need is unglamorous: an export of customers with their first order details and a table of subsequent events, or read access to the shop or CRM database. For a Solidus shop like The Gift Bow this is a single query against the orders and users tables; for a Rails SaaS like Werkbank it is workspaces and their usage events. I build the exposure correctly, fit the Poisson, run the dispersion check, move to the negative binomial when the data demands it, and only go zero inflated when we can name the group.

What you get: every customer scored with an expected yearly rate and an error band, the three lists above (expected value, surprises, levers) as CSV or straight into your email tool, the IRR table in plain language, and a memo that says what the model recommends, what it would be worth if the lever works half as well as the coefficient suggests, and where it is uncertain. You own all of it, including the code.

On cost: the count model is a well scoped piece of data science work, priced as a fixed project after the two weeks, and the follow on activity can run as performance marketing on a commission basis if you prefer paying for outcomes. The pricing page sets out how the two fit together. I am one accountable person from data to campaign, which is the point.

Questions to put to your own team

• What exactly is the count, and over what window per customer? If nobody can name the exposure, the model has no offset.
• Were the drivers measured before the counting window started, or during it?
• What is the dispersion statistic of the fitted Poisson model, and what did you do when it came out above 1.5?
• Are the effects reported as multipliers, and have they been multiplied through to a real customer's baseline?
• How many customers with this many orders does the model already expect? Is the "segment" larger than that?
• Does the predicted share of zeros match the observed share, and if not, can you name who the extra zeros are?
• Has the customer ranking been checked on a period the model never saw?
• Which list changed because of this work, and who is working through it?

This series is inspired by Mike Grigsby's Marketing Analytics (Kogan Page). The explanations, examples and numbers here are my own.