How many new kinds the next sample will show

Counting a population by its repeats ran the birthday problem backwards. Instead of a known year and the chance of a shared birthday, it took the shared birthdays and inferred the year: the number of matching pairs among the draws measures how concentrated the population is, and the number of kinds seen exactly once measures, by a rule Alan Turing and I. J. Good worked out at Bletchley Park, how much of the population’s probability lies on kinds not yet seen at all. It ended on a question that rule cannot answer — how many kinds there are in total — and on the observation that without a model of the population’s long tail, nothing can.

There is a question between the two that can be answered, exactly, and it is the one a collector actually asks. A naturalist who has spent two years trapping butterflies in Malaya and found some hundreds of species wants to know how many new species another year of trapping would add. The question is not how many species exist in Malaya; it is how many will turn up in a sample of a given size. In 1956 Good and Geoffrey Toulmin gave an answer built from nothing but the counts of species seen once, twice, three times and so on. It is unbiased, it is simple, and it has a cliff in it at exactly the point where the future sample becomes larger than the past one.

An alternating sum of the counts #

Write $f_k$ for the number of kinds that appeared exactly $k$ times in the sample: $f_1$ singletons, $f_2$ doubletons, and so on. Suppose the further sample is $t$ times as large as the first. Good and Toulmin’s estimate of the number of kinds that will appear in it and did not appear in the first is

$$U(t) = f_1 t - f_2 t^2 + f_3 t^3 - f_4 t^4 + \cdots.$$

The first term is the obvious guess: the kinds seen once are the ones that only just made it into the sample, and a similar number of kinds only just failed to, so a sample of the same size should turn up about $f_1$ new kinds. The later terms correct that guess, and they alternate because the correction alternates — the doubletons say the singletons overstated the rare tail, the tripletons that the doubletons overcorrected, and so on.

The reason the sum is exactly right is a short calculation. Suppose a kind appears in the first sample a Poisson-distributed number of times with mean $\lambda$, and in the second with mean $t\lambda$. It contributes to $f_k$ with probability $e^{-\lambda}\lambda^k/k!$, so its expected contribution to $U(t)$ is

$$\sum_{k \ge 1} (-1)^{k+1} t^k e^{-\lambda} \frac{\lambda^k}{k!} = e^{-\lambda}\bigl(1 - e^{-t\lambda}\bigr),$$

which is exactly the probability that the kind is missed by the first sample and caught by the second. Adding over all kinds, the expected value of $U(t)$ is the expected number of new kinds, whatever the population is. No assumption about the shape of the population’s tail was needed: the alternating sum knows nothing about the population, and is unbiased for every population at every $t$.

The hero figure takes one sample of about two thousand draws from a population in which the $i$-th commonest kind has a share proportional to $1/i$ — the shape that word frequencies and many species counts roughly follow, which the one that hardly ever comes up met as Zipf’s law. The sample saw 756 kinds, 561 of them once. At $t = 1$, a second sample as large as the first, the sum predicts 492 new kinds, and the true expected number is 489. For every $t$ below 1 the prediction lies on the truth. Then, at $t = 1.04$, it is off by more than the whole answer and leaves the chart.

Why the samples are Poisson #

The calculation above assumed that each kind’s count is a Poisson variable, and the figures sample that way: each kind is drawn its own Poisson number of times, independently of the others, so that the sample’s size is itself random with mean two thousand. A real naturalist’s sample has a fixed size, or a fixed duration, and then the counts of different kinds are not independent — if one kind takes more of the sample, the others must take less. For a sample in which no single kind is a large share, the difference is negligible, for the same reason that how many get their own hat found the number of people who receive their own hat settling on a Poisson law even though the hats are not handed out independently. Good and Toulmin worked with fixed-size samples and found the same sum, with corrections that vanish as the sample grows; the Poisson version is the one in which the sum is exactly unbiased and the algebra is a single line.

The independence is also what makes the population’s shape irrelevant to the mean. Each kind contributes its own term to the expectation, and those terms add whatever the kinds’ shares are. A population with a few dominant kinds and a long tail of rare ones, like the Zipf population here, and a population of equally common kinds, give the same unbiasedness. What they do not give is the same variance, and that difference is where the rest of this essay goes. It is a familiar asymmetry from any unevenness brings the match sooner, where an uneven distribution of birthdays always produced a shared birthday sooner: unevenness concentrates the sample on the common kinds, produces more repeats, and leaves fewer new kinds for a later sample to find — and here it also produces the large counts that make the extrapolation unstable.

The terms shrink, then grow #

The cliff is visible in the terms themselves.

At $t = 0.5$ each term is the count $f_k$ multiplied by a half to the power $k$, and the terms fall away so fast that the first four or five settle the sum. At $t = 1$ the terms are just the counts $f_k$, which fall slowly in a population like this one because it has a long head as well as a long tail: there are kinds seen dozens of times, and the commonest was seen 222 times. At $t = 1.5$ every term is multiplied by $1.5^k$, and for the commonest kinds that is an astronomically large number. The kind seen 222 times contributes a term of about $10^{39}$, with a sign depending on whether 222 is odd or even, to a sum whose true value is about 690.

In expectation those enormous terms cancel exactly, which is why the estimator is still unbiased. In any one sample they do not, because the counts of the common kinds fluctuate: a kind expected 220 times might appear 210 times or 235, and moving it from $f_{221}$ to $f_{222}$ changes the sum by more than $1.5^{222}$. The estimate is correct on average and has a spread that dwarfs its average.

Right on average and useless beyond one #

The figure below draws two hundred independent samples from the same population and measures, at each $t$, how much the estimates vary from sample to sample.

Up to $t = 1$ the plain sum is excellent. Averaged over two hundred samples it is within a fraction of a per cent of the truth at every $t$, and its spread from sample to sample is about 5% of the answer at $t = 1$ — most of which is the genuine randomness of the first sample, not a fault of the method. At $t = 1.25$ the spread is $10^{24}$ times the answer; at $t = 1.5$, $10^{42}$; at $t = 3$ it passes $10^{120}$, and beyond that it is too large for double-precision arithmetic to represent. The point at which it fails is not approximate. Below $t = 1$ the variance is a convergent series; above it, for any population with kinds of large mean, it diverges.

That $t = 1$ is the boundary has a clean interpretation. A sample can say a great deal about a future sample of its own size or smaller: every kind that will turn up in a smaller sample has, in expectation, already left a trace in the counts. A larger future sample will reach kinds rarer than anything the first sample could detect, and the alternating sum tries to estimate their number by extrapolating the counts to a region they do not cover. The extrapolation is unbiased, but extrapolating a power series past its radius of convergence is the problem that turning a slow series geometric faced from the other side, and it fails here for the same reason.

The common kinds are what break the sum #

Since the trouble comes from large $k$, a population without common kinds should not suffer from it. That can be tested.

With ten thousand equally common kinds and two thousand draws, each kind is expected only 0.2 times. No kind appears more than five times in any of the two hundred samples, so the sum has at most five terms, each of moderate size, and the plain estimate works far beyond $t = 1$: its typical error is 3.7% at $t = 2$ and 5.5% at $t = 3$. The same estimator on the Zipf population, with the same number of draws, explodes past $t = 1$. The difference is entirely in the head of the distribution. The kinds the estimator is trying to count are the rare ones, and the kinds that break it are the common ones, which contribute nothing to the answer because they are certain to have been seen already.

That is what makes the failure avoidable rather than fundamental. The common kinds carry no information about new kinds, so their terms can be suppressed without losing anything the sum needs, and every repair of the estimator is a way of doing that.

Truncating the sum at a random point #

The obvious repair, cutting the sum off after a fixed number of terms, introduces a bias that depends sharply on where the cut is made. Bradley Efron and Ronald Thisted, estimating in 1976 how many new words a further body of Shakespeare’s writing would contain, smoothed the sum instead, damping its later terms by factors that fall from one towards zero, and their estimate was used again in 1985 to test the attribution of a newly found poem. In 2016 Alon Orlitsky, Ananda Theertha Suresh and Yihong Wu gave a version with a proof attached. Choose a random cut-off $L$ from a Poisson distribution with mean

$$r = \frac{1}{2t} \log \frac{n (t+1)^2}{t - 1},$$

where $n$ is the size of the first sample, and keep the term $t^k f_k$ only if $k \le L$ — or, averaging over $L$, multiply the term by the probability that $L \ge k$. The weights fall off faster than $t^k$ grows, so the common kinds are suppressed, and the choice of $r$ balances the bias this introduces against the variance it removes.

The smoothed sum behaves exactly as the theory says. Up to $t = 1$ it is the plain sum. Beyond 1 its spread grows only slowly, to 8% of the answer at $t = 6$, and it acquires a bias that is always in one direction: it falls short, by 5% at $t = 2$, 18% at $t = 4$ and 26% at $t = 6$. The shortfall has the right sign to be honest. A sample of two thousand draws contains no evidence about kinds so rare that they would be expected less than once in five thousand, and a sample six times as large would turn up many of them. The smoothed estimate leaves those kinds out rather than guessing at them, and reports a lower bound that is close to the truth for small $t$ and increasingly conservative for large.

How far a sample can see #

The last figure asks how the reach of the smoothed estimate depends on the size of the first sample.

The population here is large — two hundred thousand kinds, so that even the largest sample sees only a small part of it — and the first sample is quadrupled three times, from 500 draws to 32,000. The point at which the typical error passes 10% moves from $t = 1.7$ to $2.3$, $2.6$ and $3.2$: about half a unit further for each factor of four. That is growth like the logarithm of the sample size, and Orlitsky, Suresh and Wu proved that it is the best possible. For any estimator whatsoever, there are populations on which the error stays bounded only if $t$ grows no faster than a constant times $\log n$. At the rate the figure measures, a sample a thousand times larger buys only two or three more units of $t$, and no amount of statistical ingenuity will let any sample predict a future one a thousand times its size.

The birthday problem, run forwards from the sample #

The essays on the birthday problem have moved from a known population to the chance of a repeat, and then from observed repeats back to the population. The Good–Toulmin sum sits between those directions. It uses the repeats in a sample — the counts $f_2, f_3, \ldots$ are precisely the birthday problem’s matches, sorted by how many people share each day — to predict the matches and non-matches of a future sample, without passing through any estimate of the population at all. Twenty-three people showed how quickly repeats appear in a known population; here the repeats that have appeared are used to forecast the ones that will not. And a room where nobody is alone asked how long it takes until every kind present has been seen twice, which is the question of when $f_1$ falls to zero — the moment at which the sum’s first term, and with it the whole prediction of new kinds, runs out.

The method’s history runs back to the same place as the Turing–Good rule for the missing mass. Good described in later accounts how he and Turing reasoned about frequencies of frequencies while breaking the naval Enigma, where the “kinds” were the naval cipher’s code groups and the question was how much of its codebook remained unseen; evidence measured in decibans followed the other half of that work, the scoring of evidence for and against a guess. The Good–Toulmin paper of 1956 was the peacetime sequel, and its first application was to the naturalist’s question: R. A. Fisher had been asked in the 1940s, by the entomologist Steven Corbet, how many new species of butterfly another stretch of trapping in Malaya would produce, and the alternating sum is the general answer to that question.

Still open: the tail itself #

The smoothed estimator’s reach is limited because the first sample says nothing about kinds rarer than it can detect. Every extension past $t$ of order $\log n$ needs an assumption about those kinds — a parametric form for the tail of the population, or a bound on how many kinds it can contain — and which assumptions are safe, in which applications, is the open part of the subject. Ecologists estimating species richness, linguists estimating vocabularies, and computer scientists estimating the number of distinct values in a database column all use estimators that extrapolate beyond the Good–Toulmin range, and the theory says that each of them is right only on the populations its assumptions describe.

There is a sharper question inside that one. The bound of Orlitsky, Suresh and Wu is a worst case over all populations. For the populations that occur — word frequencies, species abundances, which follow power laws with particular exponents — it is not known how much further than $\log n$ a sample can see, or whether the smoothing that is optimal in the worst case is close to optimal on them. A theory that measured the extrapolation range from properties of the observed counts themselves, and told a naturalist with a given set of $f_k$ how far ahead it was safe to predict, would be the natural completion of what Good and Toulmin began, and it does not yet exist.