The coupon collector’s problem is undoubtedly one of the most popular classical problems in combinatorial probability. It is not only of purely theoretical interest, but also finds numerous applications (see, e.g., [8, 25, 1, 4, 5]). More recently, an extended version of this problem, known in the literature as the collector’s brotherhood problem, has attracted considerable attention; see, e.g., [15, 14, 2, 11].

Following [14], we give its statement as follows. A person collects coupons, each of which belongs to one of n ∈ N different types. The coupons arrive one by one at discrete times, the type of each coupon being equiprobable and independent of types of preceding ones. Each time the person receives a coupon which he does not yet have, he puts it in his album. Otherwise, he gives it to his younger brother. In his turn, the latter puts it in his own album, if the coupon is new for him. Otherwise, he gives it to the next younger brother, and so on. All the brothers try to complete their own collections, using the same policy. The main collector is labeled 0, the younger brother 1, the next younger brother 2, etc.

For r ∈ N 0 = N ∪ { 0 }, let T r ( n ) stand for the time the rth person completes his collection. By a classical result due to Erdős and Rényi [12], (1) E T r ( n ) = n ln n + r n ln ln n + ( γ − ln r ! ) n + O ( n ) , lim n → ∞ P { T r ( n ) n − ln n − r ln ln n < x } = exp ( − e − x r ! ) , x ∈ R , with γ = − Γ ′ ( 1 ) standing for the Euler–Mascheroni constant. So, the limiting distribution is of Gumbel type.

Denote by U r ( n ) , r ∈ N, the number of empty spots in the album of the rth brother at time T 0 ( n ) , that is, when the main collector completed his album. In [26], by using the optional stopping theorem for a specially constructed martingale, it was shown that E U 1 ( n ) = H 1 ( n ) , where H 1 ( n ) = ∑ k = 1 n 1 k stands for the nth harmonic number. Subsequently, it turned out that this can be further generalized to (2) E U r ( n ) = H r ( n ) , n , r ∈ N , where H r ( n ) are hyperharmonic numbers defined by the recursive formula H r ( n ) = ∑ k = 1 n H r − 1 ( k ) k , r ≥ 2 . This was proved in [15] and [14] by using a non-elementary generating function argument, and in [2] by a simpler and more direct probabilistic reasoning. Moreover, the latter argument made it possible to obtain another, more explicit representation: H r ( n ) = ∑ k = 1 n n k ( − 1 ) k + 1 k r , r ≥ 1 . By (10.1) in [15], for the above hyperharmonic numbers the following asymptotic formula holds: (3) r ! H r ( n ) ln r n → 1 as n → ∞ . It should be noted that the term “hyperharmonic numbers” often refers to a different numerical sequence (see, e.g., the corresponding article in Wikipedia).

Note that this problem allows for different, though equivalent, formulations. For instance, we give another one in terms of a balls-into-bins model. Suppose we have an infinite sequence of balls and n bins of unlimited capacity. Each time a single ball is placed into one of the bins, which is chosen equiprobably and independently of the previous history. In this setting, T r ( n ) means the first time when each bin contains at least r + 1 balls. Similarly, U r ( n ) stands for the number of bins containing at most r balls at time T 0 ( n ) , that is, when each bin first contains at least one ball.

Until now, the question of limiting distributions for U r ( n ) as n → ∞ remained open. This problem was explicitly stated on p. 446 in [11]. In order to determine these limiting distributions it would be natural to use the generating function of U r ( n ) . The latter can be easily deduced from Proposition 4.1 in [15]: E u U r ( n ) = n u ∑ a + b + c 1 + + ⋯ + c r = n − 1 n − 1 a , b , c 1 , … , c r ( − 1 ) b ( u − 1 ) ∑ k = 1 r c k ( ∑ k = 1 r k c k ) ! ∏ k = 1 r ( k ! ) c k ( n − a ) − 1 − ∑ k = 1 r k c k . It is, however, clear that the generating function of such an ominous form can hardly be used for this purpose. Therefore, as in our previous papers [19, 20], we take a different approach based on the use of point processes. This allows not only to determine the limiting distributions for U r ( n ) (looking ahead, they will turn out to be exponential) but also to obtain the limit theorem in an ultimate infinite-dimensional form. In other words, we will prove the distributional convergence of the sequence ( U r ( n ) , r ∈ N ), normalized in a proper way, to a random element ( E , E , … ) of R ∞ as n → ∞, where the same Exp ( 1 )-distributed random variable E appears in every coordinate. This means that U r ( n ) exhibit a remarkable asymptotic stability. The precise statement is given in Theorem 1 below.

Note that in the proof of this limit theorem we use a rather unusual trick which, hopefully, can be applied elsewhere. Namely, appealing to a certain analogy between normalization and thinning, we prove the distributional convergence for thinned random variables rather than for normalized ones. This convergence can be deduced from that of specially constructed point processes. The way we then return to normalized random variables is described in Section 4.

We also prove another limit theorem which can be called delayed one. The idea behind is as follows. Due to their definition, the random variables U r ( n ) increase unboundedly as n → ∞. Thus, in order to obtain a distributional limit, we had to use some normalization. If one still wants to obtain a limit result without any normalization, one should, instead of U r ( n ) , consider the number of empty spots in the album of the rth brother not at time T 0 ( n ) but after some delay, that is, at time T 0 ( n ) + d r ( n ) with some non-random d r ( n ) → ∞ as n → ∞. With the right choice of d r ( n ) , the corresponding numbers of empty spots will not go to infinity as n does, but will instead converge in distribution to some limit law. In Theorem 2, we find this right form of d r ( n ) and establish convergence to a geometric distribution. Like the previous result, we state and prove this one in an infinite-dimensional form, that is, in the sense of distributional convergence in R ∞ .

Our third main result, Theorem 3 below, describes the asymptotic behavior of the sequence of completion times, that is, of the (centered and normalized in a proper way) random elements ( T r ( n ) , r ∈ N 0 ) in R ∞ . Actually, this is an infinite-dimensional extension of (1), the limit theorem due to Erdős and Rényi. Surprisingly, the limiting random element turns out to consist of independent entries, which means that the completion times T r ( n ) , after appropriate centering and normalizing, become asymptotically independent. This may seem all the more unexpected since T r ( n ) increase in r, which excludes any kind of independence for non-centered completion times. So, this asymptotic independence is achieved only at the expense of centering and normalizing. As a direct application of the latter result, we derive a limit theorem for inter-completion times with limiting logistic distribution.

Denote by Y i , r ( n ) , n ∈ N, i ≤ n, r ∈ N 0 , the arrival time of the ( r + 1 )th coupon of type i. So, at time Y i , r ( n ) the collector labeled r receives such a coupon into his collection. Hence, P { Y i , r ( n ) = k } = k − 1 r ( 1 n ) r + 1 ( 1 − 1 n ) k − r − 1 , k ≥ r + 1 , which is one version of the negative binomial distribution, namely, the one that counts trials up to the ( r + 1 )th success. Note that (4) T r ( n ) = max i ≤ n Y i , r ( n ) , (5) U r ( n ) = card { i : Y i , r ( n ) > max i ≤ n Y i , 0 ( n ) } .

For different i, the random variables Y i , r ( n ) are identically distributed but not independent. There is a standard way to get rid of this dependence, namely, poissonization. This technique goes back at least to Karlin and was systematically exploited by Holst in the context of random allocations (see [18] or [19] for details). Consider a poissonized scheme, that is, assume that coupons arrive at random times with independent Exp ( 1 )-distributed intervals E j , j ∈ N. Similarly to the above, denote by Z i , r ( n ) the arrival time of the ( r + 1 )th coupon of type i in the poissonized scheme. Then, for different i, Z i , r ( n ) are independent gamma-distributed random variables, Z i , r ( n ) ∼ Γ ( r + 1 , 1 n ), which is a consequence of the thinning theorem for Poisson point processes (see, e.g., Theorem 5.8 in [24]). For any fixed n, the sequences ( Y i , r ( n ) ) and ( Z i , r ( n ) ) can be given on a common probability space and coupled by (6) Z i , r ( n ) = ∑ j = 1 Y i , r ( n ) E j , i ≤ n , r ∈ N 0 . Moreover, ( Y i , r ( n ) , i ≤ n , r ∈ N 0 ) is independent of ( E j , j ∈ N ), because the inter-arrival times E j do not affect the order in which the coupons of different types arrive. Finally, note that, along with (5), a similar formula in terms of Z i , r ( n ) holds: (7) U r ( n ) = card { i : Z i , r ( n ) > max i ≤ n Z i , 0 ( n ) } .

We now state our first main result, an infinite-dimensional distributional limit theorem for ( U r ( n ) , r ∈ N ) in the space R ∞ equipped with the product topology.

Let us now turn to the delayed limit theorem. By analogy with U r ( n ) , for fixed d r ( n ) ∈ N denote by U ˆ r ( n ) the number of empty spots in the album of the rth brother at time T 0 ( n ) + d r ( n ) : (9) U ˆ r ( n ) = card { i : Y i , r ( n ) > max i ≤ n Y i , 0 ( n ) + d r ( n ) } . The following result asserts that, with the right choice of d r ( n ) , the random elements ( U ˆ r ( n ) , r ∈ N ) of R ∞ converge in distribution to some limiting random element with geometrically distributed marginals.

As our final result, we give an infinite-dimensional extension of (1), the limit theorem by Erdős and Rényi.

In particular, this result allows obtaining limit distributions for times between completions of different collections. Denote by Δ r 1 , r 2 ( n ) = T r 2 ( n ) − T r 1 ( n ) , r 2 > r 1 , such an inter-completion time.

It is well known (and may be easily checked by means of characteristic functions) that the difference of two independent Gumbel-distributed random variables has a logistic distribution. Particularly, B r 2 − B r 1 = d L r 1 , r 2 . Thus, (17) follows from (16).

One of the key points in the proof of Theorem 1 is convergence of specially constructed point processes to a Poisson one. We now proceed to the corresponding construction.

For fixed n ∈ N and r ∈ N 0 , let (18) ψ ( n ) ( x ) = x n − ln n , x ∈ R , (19) η r ( n ) = ∑ i = 1 n δ ψ ( n ) ( Z i , r ( n ) ) , where δ u stands for the Dirac measure 𝟙 { u ∈ · }. The processes η r ( n ) describe ( r + 1 )th arrivals of different coupon types in the poissonized scheme. Note that in a similar case in [19], in order to provide the necessary convergence, we used an r-dependent centering/normalizing function (20) ψ r ( n ) ( x ) = x n − ln n − r ln ln n , x ∈ R . In our case, however, we want to consider the numbers of empty spots in different albums at the same point in time, namely when the 0th album is completed. Hence, we have to deal with the same centering for different r. So, we will achieve the desired convergence in a different way — by means of thinnings.

For a proper point process η, let us denote by T p η, p ∈ [ 0 , 1 ], its p-thinning, i.e. the point process which independently keeps the points of η with probability p and removes otherwise (see, e.g., [24], Section 5.3 for details). On the space (21) X = N 0 × ( R ∪ { + ∞ } ) , endowed with some relevant metric, say, (22) d ( ( r 1 , x ) , ( r 2 , y ) ) = | r 2 − r 1 | + | e − y − e − x | , consider the point processes H ( n ) , n ≥ 3, given as follows: for B r ∈ B ( R ∪ { + ∞ } ), r ∈ N 0 , define (23) H ( n ) ( ⋃ r = 0 ∞ ( { r } × B r ) ) = ∑ r = 0 ∞ ( T ln − r n η r ( n ) ) ( B r ) , where all thinnings are performed independently across levels and of the underlying coupon process. (We require n ≥ 3 in order to have ln − r n ∈ [ 0 , 1 ] for all r ∈ N 0 .) In other words, we are thinning out the processes η r ( n ) and glue them into one “multilevel” point process H ( n ) . Note that, by construction in (19), (24) ( T ln − r n η r ( n ) ) ( { + ∞ } ) = η r ( n ) ( { + ∞ } ) = 0 anyway. The reason we, nevertheless, consider the semi-compactified real axis R ∪ { + ∞ } instead of just R will become clear from what follows (see Section 5).

Before stating the main theorem of this section, we reformulate the basic definitions related to convergence of point processes as applied to our problem (for a detailed exposition in the abstract setting, see [28, 29], or [22]). Let M p ( X ) denote the space of all locally finite (with respect to d in (22)) point measures on X given by (21). For μ , μ 1 , μ 2 , … ∈ M p ( X ), μ n are said to converge vaguely to μ (denoted by μ n → v μ) if ∫ X f d μ n → ∫ X f d μ for each continuous compactly supported non-negative test function f defined on X. In our case, this means that ∫ R ∪ { + ∞ } g ( t ) μ n ( { r } × d t ) → ∫ R ∪ { + ∞ } g ( t ) μ ( { r } × d t ) for each r ∈ N 0 and each continuous compactly supported g : R ∪ { + ∞ } → [ 0 , + ∞ ). As usual, the set M p ( X ), equipped with the topology of the above convergence, can be metrized as a complete separable metric space. This setting allows us to consider the distributional convergence of point processes H ( 3 ) , H ( 4 ) … , denoted as H ( n ) → v d H. The main result of this section, Theorem 4 below, asserts that the point processes H ( n ) converge in this sense toward a non-homogeneous Poisson process. Theorem 4.

Let H be a Poisson point process on X with intensity measure λ given by (25) λ ( ⋃ r = 0 ∞ ( { r } × B r ) ) = ∑ r = 0 ∞ 1 r ! ∫ B r e − x d x , B r ∈ B ( R ∪ { + ∞ } ) . Then H ( n ) → v d H as n → ∞.

In this section, we consider a somewhat unusual approach to proving limit theorems, which is the basic tool in the proof of Theorem 1. This approach does not seem to have been used in this context before, and may prove useful elsewhere as well.

For an N 0 -valued random variable X, define its p-thinning, p ∈ [ 0 , 1 ], by (34) p ⊙ X = ∑ i = 1 X ε i , where ε i are Bin ( 1 , p )-distributed and independent of each other and of X. Note that, in the notation of the previous section, ( T p η ) ( B ) = d p ⊙ η ( B ) for a proper point process η and a Borel set B. This operation, going back to [27], was then used in [32] to introduce the concepts of discrete self-decomposability and stability.

The idea behind the proposed technique is to replace the normalization by thinning. In other words, instead of proving that p ( n ) X ( n ) → d Y with some constants p ( n ) → 0 and a limiting random variable Y, we will prove that p ( n ) ⊙ X ( n ) → d Z with some new limiting random variable Z. The connection between the distributions of Y and Z can be guessed from the following examples. For X ( n ) = n a.s., we have 1 n X ( n ) = 1 → d 1 and 1 n ⊙ X ( n ) → d Z ∼ Pois ( 1 ) by the Poisson limit theorem. More generally, consider an i.i.d. sequence ( ξ i , i ∈ N ) with E ξ i = λ, and denote X ( n ) = ∑ i = 1 n ξ i . Due to the law of large numbers, 1 n X ( n ) → d λ, and, by Theorem 3.1 in [17], 1 n ⊙ X ( n ) → d Z ∼ Pois ( λ ). All this suggests that, in general, Z must have the mixed Poisson distribution with mixing distribution of Y (see, e.g., Chapter 2 in [16]). This means that (35) P { Z = k } = ∫ [ 0 , + ∞ ) e − y y k k ! F Y ( d y ) , k ∈ N 0 , where F Y stands for the distribution function of Y.

The main result of this section, Theorem 5 below, justifies this approach. With an eye to the future, we will state and prove it in multidimensional form. Let Y = ( Y 1 , … , Y s ) be a random vector with a.s. non-negative components. By analogy with (35), the random vector Z = ( Z 1 , … , Z s ) is said to have a multivariate mixed Poisson distribution with mixing distribution of Y (see [13] or [23]) if (36) P { Z 1 = k 1 , … , Z s = k s } = ∫ · · · ∫ [ 0 , + ∞ ) s e − y 1 y 1 k 1 k 1 ! · … · e − y s y s k s k s ! F Y ( d y 1 , … , d y s ) , k 1 , … , k s ∈ N 0 . Remark 5.

In what follows, we will need another equivalent interpretation of the multivariate mixed Poisson distribution. Namely, we may define Z by Z r = N r ( Y r ), r = 1 , … , s, where N r stand for unit-rate Poisson counting processes, independent of each other and of Y.