In probability theory, the coupon collector's problem refers to mathematical analysis of "collect all coupons and win" contests. It asks the following question: if each box of a given product (e.g., breakfast cereals) contains a coupon, and there are n different types of coupons, what is the probability that more than t boxes need to be bought to collect all n coupons? An alternative statement is: given n coupons, how many coupons do you expect you need to draw with replacement before having drawn each coupon at least once? The mathematical analysis of the problem reveals that the expected number of trials needed grows as Θ ( n log ( n ) ) {\displaystyle \Theta (n\log(n))} . For example, when n = 50 it takes about 225 trials on average to collect all 50 coupons. Sometimes the problem is instead expressed in terms of an n-sided die.
Solution
Calculating the expectation Let time T be the number of draws needed to collect all n coupons, and let ti be the time to collect the i-th coupon after i − 1 coupons have been collected. Then T = t 1 + ⋯ + t n {\displaystyle T=t_{1}+\cdots +t_{n}} . Think of T and ti as random variables. Observe that the probability of collecting the i-th new coupon is p i = n − ( i − 1 ) n = n − i + 1 n {\displaystyle p_{i}={\frac {n-(i-1)}{n}}={\frac {n-i+1}{n}}} . Therefore, t i {\displaystyle t_{i}} has geometric distribution with expectation 1 p i = n n − i + 1 {\displaystyle {\frac {1}{p_{i}}}={\frac {n}{n-i+1}}} . By the linearity of expectations we have:
E ( T )
= E ( t 1 + t 2 + ⋯ + t n )
= E ( t 1 ) + E ( t 2 ) + ⋯ + E ( t n )
= 1 p 1 + 1 p 2 + ⋯ + 1 p n
= n n + n n − 1 + ⋯ + n 1
= n ⋅ ( 1 1 + 1 2 + ⋯ + 1 n )
= n ⋅ H n . {\displaystyle {\begin{aligned}\operatorname {E} (T)&{}=\operatorname {E} (t_{1}+t_{2}+\cdots +t_{n})\\&{}=\operatorname {E} (t_{1})+\operatorname {E} (t_{2})+\cdots +\operatorname {E} (t_{n})\\&{}={\frac {1}{p_{1}}}+{\frac {1}{p_{2}}}+\cdots +{\frac {1}{p_{n}}}\\&{}={\frac {n}{n}}+{\frac {n}{n-1}}+\cdots +{\frac {n}{1}}\\&{}=n\cdot \left({\frac {1}{1}}+{\frac {1}{2}}+\cdots +{\frac {1}{n}}\right)\\&{}=n\cdot H_{n}.\end{aligned}}}
Here Hn is the n-th harmonic number. Using the asymptotics of the harmonic numbers, we obtain:
E ( T ) = n ⋅ H n = n log n + γ n + 1 2 + O ( 1 / n ) , {\displaystyle \operatorname {E} (T)=n\cdot H_{n}=n\log n+\gamma n+{\frac {1}{2}}+O(1/n),}
where γ ≈ 0.5772156649 {\displaystyle \gamma \approx 0.5772156649} is the Euler–Mascheroni constant. Using the Markov inequality to bound the desired probability:
P ( T ≥ c n H n ) ≤ 1 c . {\displaystyle \operatorname {P} (T\geq cnH_{n})\leq {\frac {1}{c}}.}
The above can be modified slightly to handle the case when we've already collected some of the coupons. Let k be the number of coupons already collected, then:
E ( T k )
= E ( t k + 1 + t k + 2 + ⋯ + t n )
= n ⋅ ( 1 1 + 1 2 + ⋯ + 1 n − k )
= n ⋅ H n − k {\displaystyle {\begin{aligned}\operatorname {E} (T_{k})&{}=\operatorname {E} (t_{k+1}+t_{k+2}+\cdots +t_{n})\\&{}=n\cdot \left({\frac {1}{1}}+{\frac {1}{2}}+\cdots +{\frac {1}{n-k}}\right)\\&{}=n\cdot H_{n-k}\end{aligned}}}
And when k = 0 {\displaystyle k=0} then we get the original result.
Calculating the variance Using the independence of random variables ti, we obtain:
Var ( T )
= Var ( t 1 + ⋯ + t n )
= Var ( t 1 ) + Var ( t 2 ) + ⋯ + Var ( t n )
= 1 − p 1 p 1 2 + 1 − p 2 p 2 2 + ⋯ + 1 − p n p n 2
= ( n 2 n 2 + n 2 ( n − 1 ) 2 + ⋯ + n 2 1 2 ) − ( n n + n n − 1 + ⋯ + n 1 )
= n 2 ⋅ ( 1 1 2 + 1 2 2 + ⋯ + 1 n 2 ) − n ⋅ ( 1 1 + 1 2 + ⋯ + 1 n )
< π 2 6 n 2 {\displaystyle {\begin{aligned}\operatorname {Var} (T)&{}=\operatorname {Var} (t_{1}+\cdots +t_{n})\\&{}=\operatorname {Var} (t_{1})+\operatorname {Var} (t_{2})+\cdots +\operatorname {Var} (t_{n})\\&{}={\frac {1-p_{1}}{p_{1}^{2}}}+{\frac {1-p_{2}}{p_{2}^{2}}}+\cdots +{\frac {1-p_{n}}{p_{n}^{2}}}\\&{}=\left({\frac {n^{2}}{n^{2}}}+{\frac {n^{2}}{(n-1)^{2}}}+\cdots +{\frac {n^{2}}{1^{2}}}\right)-\left({\frac {n}{n}}+{\frac {n}{n-1}}+\cdots +{\frac {n}{1}}\right)\\&{}=n^{2}\cdot \left({\frac {1}{1^{2}}}+{\frac {1}{2^{2}}}+\cdots +{\frac {1}{n^{2}}}\right)-n\cdot \left({\frac {1}{1}}+{\frac {1}{2}}+\cdots +{\frac {1}{n}}\right)\\&{}<{\frac {\pi ^{2}}{6}}n^{2}\end{aligned}}}
since π 2 6 = 1 1 2 + 1 2 2 + ⋯ + 1 n 2 + ⋯ {\displaystyle {\frac {\pi ^{2}}{6}}={\frac {1}{1^{2}}}+{\frac {1}{2^{2}}}+\cdots +{\frac {1}{n^{2}}}+\cdots } (see Basel problem). Bound the desired probability using the Chebyshev inequality:
P ( | T − n H n | ≥ c n ) ≤ π 2 6 c 2 . {\displaystyle \operatorname {P} \left(|T-nH_{n}|\geq cn\right)\leq {\frac {\pi ^{2}}{6c^{2}}}.}
Stirling numbers Let the random variable X be the number of dice rolls performed before all faces have occurred. The subpower is defined a { b } = a ! { b a } {\displaystyle a^{\{b\}}=a!\left\{{b \atop a}\right\}} , where { b a } {\displaystyle \left\{{b \atop a}\right\}} is a Stirling number of the second kind. Sequences of k {\displaystyle k} die rolls are functions k → n {\displaystyle k\rightarrow n} counted by n k {\displaystyle n^{k}} , while surjections (that land on each face at least once) are counted by n { k } {\displaystyle n^{\{k\}}} , so the probability that all faces were landed on within the k-th throw is P ( X ≤ k ) = n { k } n k {\displaystyle P(X\leq k)={\frac {n^{\{k\}}}{n^{k}}}} . By the recurrence relation of the Stirling numbers, the probability that exactly k rolls are needed is P ( X = k ) = n { k } n k − n { k − 1 } n k − 1 = ( n − 1 ) { k − 1 } n k − 1 {\displaystyle P(X=k)={\frac {n^{\{k\}}}{n^{k}}}-{\frac {n^{\{k-1\}}}{n^{k-1}}}={\frac {(n-1)^{\{k-1\}}}{n^{k-1}}}}
Generating functions Replacing z {\displaystyle z} with 1 + z {\displaystyle 1+z} in the probability generating function produces the o.g.f. for E [ ( X k ) ] {\displaystyle E\left[{X \choose k}\right]} . Using the partial fraction decomposition ( 1 x − 1 n ) − 1 = ∑ k = 0 n ( n k ) ( − 1 ) n − k 1 − k x {\displaystyle {{\frac {1}{x}}-1 \choose n}^{-1}=\sum _{k=0}^{n}{n \choose k}{\frac {(-1)^{n-k}}{1-kx}}} , we can take the expansion
( n x + 1 n ) − 1 = ∑ i = 0 n ( n i ) ( − 1 ) n − i 1 − i ( 1 − n x + 1 + n ) = ∑ i = 0 n ( n i ) ( − 1 ) n − i ( 1 + n 1 + n − i + i n ∑ k = 1 ∞ ( i − 1 ) k − 1 ( n + 1 − i ) k + 1 x k ) {\displaystyle {\begin{aligned}&{{\frac {n}{x+1}} \choose n}^{-1}\\=&\sum _{i=0}^{n}{n \choose i}{\frac {(-1)^{n-i}}{1-i(1-{\frac {n}{x+1+n}})}}\\=&\sum _{i=0}^{n}{n \choose i}(-1)^{n-i}\left({\frac {1+n}{1+n-i}}+in\sum _{k=1}^{\infty }{\frac {(i-1)^{k-1}}{(n+1-i)^{k+1}}}x^{k}\right)\end{aligned}}}
revealing that for k > 0 {\displaystyle k>0} ,
E [ ( X k ) ] = n ∑ i = 0 n ( n i ) ( − 1 ) n − i i ( i − 1 ) k − 1 ( n + 1 − i ) k + 1 {\displaystyle E\left[{X \choose k}\right]=n\sum _{i=0}^{n}{n \choose i}(-1)^{n-i}i{\frac {(i-1)^{k-1}}{(n+1-i)^{k+1}}}}
Given an o.g.f. f, since ( x 1 − x ) i = ∑ n = 0 ∞ ( k − 1 i − 1 ) x k {\displaystyle \left({\frac {x}{1-x}}\right)^{i}=\sum _{n=0}^{\infty }{k-1 \choose i-1}x^{k}} , a variation of the binomial transform is [ x k ] f ( x 1 + x ) = ∑ i = 0 k ( k − 1 i − 1 ) ( − 1 ) k − i [ x i ] f ( x ) {\displaystyle [x^{k}]f\left({\frac {x}{1+x}}\right)=\sum _{i=0}^{k}{k-1 \choose i-1}(-1)^{k-i}[x^{i}]f(x)} . (Specifically, if ( n x + 1 n ) − 1 = f ( x 1 + x ) {\displaystyle {{\frac {n}{x+1}} \choose n}^{-1}=f\left({\frac {x}{1+x}}\right)} , f ( x ) = ( n − n x n ) − 1 {\displaystyle f(x)={n-nx \choose n}^{-1}} .) Rewriting the binomial coefficient via the gamma function and expanding as the exp {\displaystyle \exp } of the polygamma series (in terms of generalised harmonic numbers), we find [ x i i ! ] ( n − x n ) − 1 = ∑ P ∈ p e r m s ( i ) ∏ c ∈ P H n ( | c | ) {\displaystyle \left[{\frac {x^{i}}{i!}}\right]{n-x \choose n}^{-1}=\sum _{P\in \mathrm {perms} (i)}\prod _{c\in P}H_{n}^{(|c|)}} , so
E [ ( X k ) ] = ∑ i = 0 k ( k − 1 i − 1 ) ( − 1 ) k − i n i