HOW MANY TRIALS IT TAKES TO FIND ALL DIFFERENT SPECIES OF A KIND - - PDF document

HOW MANY TRIALS IT TAKES TO FIND ALL DIFFERENT SPECIES OF A KIND

Vassilis G. Papanicolaou Department of Mathematics National Technical University of Athens Zografou Campus 157 80, Athens, GREECE

e-mail: papanico@math.ntua.gr 1

Consider a population whose members are of N different types (e.g. colors). For 1 ≤ j ≤ N we denote by pj the probability that a member of the population is of type j. The members of the population are sampled independently with replacement and their types are recorded. Our main object of study is the number TN of trials it takes until all N types are detected (at least once). Of course, P {TN ≥ k} = 1, if 1 ≤ k ≤ N. 2

It is convenient to introduce the events Ak

j, 1 ≤ j ≤ N, that the type j is not detected

until trial k (included). Then P {TN ≥ k} = P

• Ak−1

1

∪ · · · ∪ Ak−1

N

• ,

k = 1, 2, ... . By invoking the inclusion-exclusion principle one gets P {TN ≥ k} = P

• Ak−1

1

• + · · · + P
• Ak−1

N

• −P
• Ak−1

1

Ak−1

2

• − · · · − P
• Ak−1

N−1Ak−1 N

• .

. . +(−1)N−1P

• Ak−1

1

· · · Ak−1

N

• ,
• r, in a more compact notation

P {TN ≥ k} =

• J⊂{1,...,N}

J=∅

(−1)|J|−1 P

• j∈J

Ak−1

j

• ,

(1) where the sum extends over all 2N − 1 nonempty subsets J of {1, ..., N}, while |J| denotes the cardinality of J. 3

Now P

• Ak−1

j

• = (1 − pj)k−1 ,

P

• Ak−1

j

Ak−1

i

• = [1 − (pj + pi)]k−1 ,

and, in general, if J ⊂ {1, ..., N}, then P

• j∈J

Ak−1

j

• =
• 1 −
• j∈J

pj k−1 . (2) Combining formulas (2) and (1) yields P {TN ≥ k} =

• J⊂{1,...,N}

J=∅

(−1)|J|−1

• 1 −
• j∈J

pj k−1 , k = 1, 2, ... (3) (valid also for the trivial case k = 1 under the convention 00 = 1). 4

• Remark. A side result of the above analysis is the (somehow nontrivial) algebraic formula
• J⊂{1,...,N}

(−1)|J|

• 1 −
• j∈J

pj n = 0, for n = 0, 1, ..., (N − 1). 5

Lemma 0. Let X be a random variable taking values in N = {0, 1, 2, ...}. If g : N → R is a function such that E [g(X)] < ∞, then E [g(X)] = g(0) +

∞

• k=1

[g(k) − g(k − 1)] P {X ≥ k} . (4) In particular (g(k) = k), E [X] =

∞

• k=1

P {X ≥ k} . (5) 6

• Proof. We have

E [g(X)] =

∞

• k=0

g(k)P {X = k} = g(0)P {X = 0} + g(1)P {X = 1} + g(2)P {X = 2} + g(3)P {X = 3} + · · ·. The above sum can be rewriten as g(0) [P {X = 0} + P {X = 1} + P {X = 2} + P {X = 3} + · · ·] + [g(1) − g(0)] [P {X = 1} + P {X = 2} + P {X = 3} + · · ·] + [g(2) − g(1)] [P {X = 2} + P {X = 3} + · · ·] + [g(3) − g(2)] [P {X = 3} + · · ·] . . . which is the right-hand side of (4).

• Remark. If g : N → R is increasing, then (4) is valid even if E [g(X)] = ∞ (similarly if

g is decreasing). 7

Combination of (3) in (5) yields E [TN] =

• J⊂{1,...,N}

J=∅

(−1)|J|−1

∞

• k=1
• 1 −
• j∈J

pj k−1 . Summation of the geometric series gives E [TN] =

• J⊂{1,...,N}

J=∅

(−1)|J|−1

• j∈J

pj −1 , (6)

• r

E [TN] =

N

• m=1

(−1)m−1

• 1≤j1<···<jm≤N

1 pj1 + · · · + pjm . 8

We proceed by noticing that 1 −

N

• j=1
• 1 − e−pjt

=

• J⊂{1,...,N}

J=∅

(−1)|J|−1 exp

• −t
• j∈J

pj

• .

Thus, ∞

• 1 −

N

• j=1
• 1 − e−pjt
• dt =
• J⊂{1,...,N}

J=∅

(−1)|J|−1

• j∈J pj

. and hence E [TN] = ∞

• 1 −

N

• j=1
• 1 − e−pjt
• dt,

(7)

• r, by substituting x = e−t in the integral,

E [TN] = 1

• 1 −

N

• j=1

(1 − xpj)

• dx

x . (8) Formulas (6), (7), and (8) are well-known. 9

Likewise, for the generating function of TN F(z) := E

• z−TN

, we have the formulas F(z) = 1 − (z − 1)

• J⊂{1,...,N}

J=∅

(−1)|J|−1 z − 1 +

• j∈J pj

, F(z) = 1 − (z − 1) ∞

• 1 −

N

• j=1
• 1 − e−pjt
• e−(z−1)tdt,

and F(z) = 1 − (z − 1) 1

• 1 −

N

• j=1

(1 − xpj)

• xz−2dx.

10

And for the second moment of TN we have E [TN (TN + 1)] = 2

• J⊂{1,...,N}

J=∅

(−1)|J|−1

• j∈J

pj −2 E [TN (TN + 1)] = 2 ∞

• 1 −

N

• j=1
• 1 − e−pjt
• t dt,

and E [TN (TN + 1)] = −2 1

• 1 −

N

• j=1

(1 − xpj)

• ln x

x dx. 11

Naturally, the simplest case regarding the previous formulas occurs when one takes p1 = · · · = pN = 1 N . It is easy to check (e.g. by taking logarithms) that, for any fixed t > 0, the maximun of the quantity

N

• j=1
• 1 − e−pjt

, subject to the constraint p1 + · · · + pN = 1,

• ccurs when all pj’s are equal.

It follows that E [TN] attains its minimum value when all pj’s are equal (see also M.V. Hildebrand [11]). The same is true for E [TN (TN + 1)]. As for F(z) = E

• z−TN

, where z > 1, it follows that it attains its maximum value, when all pj’s are equal. 12

Let p1 = · · · = pN = 1/N. Then E [TN] = N

N

• m=1

N m (−1)m−1 m = 1

• 1 −
• 1 − x1/NN dx

x . The substituting u = 1 − x1/N in the integral yields E [TN] = N 1 1 − uN 1 − u du = N 1 N

• m=1

um−1

• du = NHN,

where HN is “the N-th harmonic number” HN =

N

• m=1

1 m. 13

In a similar way we get F(z) = 1 − (z − 1)

N

• m=1

N m

• (−1)m−1

z − 1 + (m/N) = 1 − (z − 1) 1

• 1 −
• 1 − x1/NN

xz−2dx, and E [TN (TN + 1)] = 2N 2

N

• m=1

Hm m = N 2

• H2

N + N

• m=1

1 m2

• ,

which also implies V [TN] = N 2

N

• m=1

1 m2 − NHN. 14

Reminder (The Euler-Maclaurin sum formula). If S(N) =

N

• m=0

f(m), then, as N → ∞, S(N) ∼ 1 2f(N) + N f(x)dx + C +

∞

• k=1

(−1)k+1 Bk+1 (k + 1)!f (k)(N), where C is a constant and Bk is the k-th Bernoulli number defined by the formula z ez − 1 =

∞

• k=0

Bk k! zk (e.g. B0 = 0, B1 = −1/2, B2 = 1/6, B4 = −1/30, Bk = 0, for all odd k ≥ 3). For example HN =

N

• m=1

1 m ∼ ln N + γ + 1 2N −

∞

• k=2

Bk kN k , where γ = 0.5772... is Euler’s constant. Also

N

• m=1

1 m2 ∼ π2 6 − 1 N + 1 2N 2 −

∞

• k=2

Bk N k+1. 15

If p1 = · · · = pN = 1/N, then by using the Euler-Maclaurin sum formula one obtains E [TN] ∼ N ln N + γN + 1 2 −

∞

• k=2

Bk kN k−1, E [TN (TN + 1)] = N 2

• (ln N)2 + 2γ ln N + γ2 + π2

6 + O ln N N

• ,

and V [TN] = π2N 2 6 − N ln N − (γ + 1)N + O ln N N

• ,

as N → ∞. 16

• QUIZ. The town F has population 1825 (= 5 × 365) while the town S has population

2190 (= 6 × 365). Let f and s be the probabilities that all 365 birthdays are represented in F and S respectively. Estimate s f . 17

Answer. s f ≈ 4.8. In fact, f ≈ 0.085 and s ≈ 0.4051

• Hint. It can be shown (see R. Durrett [8]) that, as N → ∞,

P {TN − N ln N ≤ Nx} → exp

• −e−x

. 18

Some Applications. The above formulas are associated to what is usually called the “Coupon Collector Problem” (CCP), where N different coupons (of arbitrary occurring fre- quencies) are sampled independently with replacement. We now mention some applications. The first three examples introduce probabilistic computational algorithms which can be modeled by the CCP.

• 1. Constraint classification in mathematical programming. In 1983, Karwan et al. [14]

described a class of randomized algorithms for classifying all the constraints in a mathemat- ical programming problem as necessary or redundant. The basic algorithm, also known as PREDUCE (Probabilistic REDUCE), can be briefly described as follows: Given an interior feasible point, each iteration consists of generating a ray in a random direction, and recording the first constraint it intersects. Such a constraint is a necessary one. The algorithm gener- ates rays until a stopping rule is satisfied. Then, all the constraints which were not hit at all are classified as redundant—possibly erroneously. Each iteration corresponds to drawing one coupon, with N being the number of necessary constraints. Thus, the CCP model can help to determine an efficient stopping rule. 19

• 2. Multistart methods in global optimization. In an environment with multiple local op-

tima, this method chooses a random initial point and uses any optimization algorithm to find a local optimum of the objective function. We repeatedly pick different random initial points, which in turn lead to different local optima, according to the partition of the domain to at- traction regions. Here, the local optima are the coupons, and the probability of finding any local optimum is proportional to the size of its attraction region. The global optimizer wants to know how long it is expected to take until all the local optima are found, thus yielding the global optimum. For more information we refer to A. Torn and A. Zallinskas [18]. 20

• 3. Determining the convex closure of a finite subset S in Rd. To determine the subset of

S that spans its convex closure, we generate a random (d − 1)-dimensional hyperplane and measure the orthogonal distance of each point from that hyperplane. It can be shown that the point with the greatest distance belongs to the convex closure. This is done repeatedly until a stopping rule indicates that all the spanning points were detected. 21

The next two examples illustrate applications of the CCP in engineering and natural sci- ences.

• 4. Engineering application. Fault detection in electronic hardware. Many fault detection

procedures in electronic hardware are based on repeated tests of the hardware. In each such test, a certain number of faults are detected, while some other faults are not detected. Here, each detection of a fault corresponds to a drawing of a coupon. It is of interest to determine the expected number of tests necessary until all the faults are detected. 22

• 5. Biological application. Estimating the number of species. Consider a complex ecosys-

tem inhabited by a great variety of species. Biologists who explore that environment want to detect all the different species of a certain type. For example, all the different species of birds in the Galapagos Islands near Equador. Observing an animal from a given species constitutes a detection of that species. The reader should easily see that this scenario falls into the category

• f the CCP.

23

Large N Asymptotics. When N is large it is not obvious at all what information one can

• btain from the formulas for E [TN] and for E [TN (TN + 1)].

For this reason there is a need to develop efficient ways for deriving asymptotics of E [TN] and E [TN (TN + 1)], as N → ∞ (we have already analyzed the very special case of equal probabilities). Reminder (see, e.g., R. Durrett [8]). If the quantity bN grows fast enough so that V [TN] =

• (b2

N), as N → ∞, then

TN − E [TN] bN → 0 in probability. The asymptotics of E [TN] and E [TN (TN + 1)] are necessary for a good choice of bN. 24

How to formulate the asymptotics problem? Let α = {aj}∞

j=1 be a sequence of strictly positive numbers. Then, for each integer

N > 0, one can create a probability measure πN = {p1, ..., pN} on the set {1, ..., N} by taking pj = aj AN , where AN =

N

• j=1

aj. Notice that pj depends on α and N, thus, given α, it makes sense to consider the asymptotic behavior of E [TN], as N → ∞. 25

Motivated by (6) we introduce the notation EN(α) :=

• J⊂{1,...,N}

J=∅

(−1)|J|−1

• j∈J

aj −1 (9) =

N

• m=1

(−1)m−1

• 1≤j1<···<jm≤N

1 aj1 + · · · + ajm . Then, as in (7) and (8), one has EN(α) = ∞

• 1 −

N

• j=1
• 1 − e−ajt
• dt

(10) and EN(α) = 1

• 1 −

N

• j=1

(1 − xaj)

• dx

x . (11) 26

Since EN(sα) = 1 sEN(α), we have E [TN] = EN(A−1

N α) = ANEN(α).

(12) Thus, the problem of estimating E[TN], as N → ∞, can be treated as two separate problems, namely estimating AN and estimating EN(α). Our analysis focuses on estimating EN(α). The estimation of the sum AN will be considered an external matter which can be handled by existing powerful methods, such as the Euler-Maclaurin sum formula, the Laplace method for sums (see, e.g., C.M. Bender and S.A. Orszag [3]), or even summation by parts. 27

Likewise, in order to analyze E [TN (TN + 1)] we should introduce QN(α) := 2

• J⊂{1,...,N}

J=∅

(−1)|J|−1

• j∈J

aj −2 Then one has E [TN(TN + 1)] = QN(A−1

N α) = A2 NQN(α).

28

The Dichotomy. For convenience, we denote f α

N(x) = N

• j=1

(1 − xaj), 0 ≤ x ≤ 1. The following properties of f α

N are immediate:

(i) f α

N(0) = 1,

f α

N(1) = 0;

(ii) f α

N is monotone decreasing in x;

(iii) f α

N+1(x) ≤ f α N(x), in particular

lim

N f α N(x) = ∞

• j=1

(1 − xaj)

• exists. Thus, by applying the MCT we get the existence of

L1(α) := lim

N EN(α) =

1

• 1 −

∞

• j=1

(1 − xaj)

• dx

x (13) Notice that L1(α) > 0, for any α (since, for every x ∈ [0, 1), f α

N(x) < 1 and decreases with

N). However, we may have L1(α) = ∞. Therefore, we have the following dichotomy: (i) 0 < L1(α) < ∞

• r

(ii) L1(α) = ∞. (14) 29

L2(α) := lim

N QN(α) = −2

1

• 1 −

∞

• j=1

(1 − xaj)

• ln x

x dx. 30

The Case L1(α) < ∞. Obviously, as N → ∞, E [TN] = ANL1(α) [1 + o(1)] . (15) It would be desirable to give a simple necessary and sufficient condition for L1(α) < ∞. Theorem 1. L1(α) < ∞ if and only if there exists a ξ ∈ (0, 1) such that

∞

• j=1

ξaj < ∞. (16) This condition is easily applied. 31

We need the following lemma: Lemma 1. Let {bj}∞

j=1 be a sequence of real numbers such that 0 ≤ bj ≤ 1, for all j. If

∞

j=1 bj < ∞, then ∞

• j=1

bj −

• 1≤l<j

blbj ≤ 1 −

∞

• j=1

(1 − bj) ≤

∞

• j=1

bj.

• Proof. Induction and limit.

32

Proof of Th. 1. Assume first that there is a ξ ∈ (0, 1) such that (16) is true. Then, by (13) and Lemma 1 we have L1(α) ≤ ξ ∞

• j=1

xaj

• dx

x + 1

ξ

• 1 −

∞

• j=1

(1 − xaj)

• dx

x ≤ ξ ∞

• j=1

xaj−1

• dx−ln ξ. (17)

Now, (16) implies that ξaj → 0, hence aj → ∞. Therefore minj {aj} = aj0 > 0.Evaluating the integral in (17) gives L1(α) ≤

∞

• j=1

ξaj aj − ln ξ ≤ 1 aj0

∞

• j=1

ξaj − ln ξ < ∞. Conversely, if

∞

• j=1

ξaj = ∞, for all ξ ∈ (0, 1], then, by a well-known property of infinite products (see, e.g. [17])

∞

• j=1

(1 − xaj) = 0, for all x ∈ (0, 1], and hence (13) yields L1(α) = 1 dx x = ∞. 33

• Remarks. (a) If lim inf aj < ∞, then L1(α) = ∞.

On the other hand, lim aj = ∞ does not imply that L1(α) < ∞, since we may have ξaj → 0, but ∞

j=1 ξaj = ∞. For example, if aj = (ln j)p, 0 < p < 1, then lim aj = ∞ and

L1(α) = ∞. However, if

∞

• j=1

1/aj < ∞, then by Tonelli’s theorem 1 ∞

• j=1

xaj

• dx

x =

∞

• j=1

1 xaj dx x

• =

∞

• j=1

1 aj < ∞, hence ∞

j=1 xaj < ∞, for x ∈ (0, 1). Therefore L1(α) < ∞.

34

(b) Consider the error term, defined by δN = L1(α) − EN(α). Then δN = 1

N

• j=1

(1 − xaj)

• 1 −

∞

• j=N+1

(1 − xaj)

• dx

x ≤ 1

• ∞
• j=N+1

xaj

• dx

x =

∞

• j=N+1

1 aj . (18) 35

The Case L1(α) = ∞. The problem of determining the asymptotic behavior of EN(α) becomes much more challenging in this case. We begin the analysis by writing aj in the form aj = 1 f(j), where f(x) > 0 is a real function. For our main result we need to impose the following conditions on the limiting behavior

• f f(x) as x → ∞:

(i) f(x) ր ∞, (ii) f ′(x) f(x) ց 0 (iii) f ′′(x)/f ′(x) [f ′(x)/f(x)] ln [f ′(x)/f(x)] → 0. 36

Condition (iii) is merely technical. If we introduce the function g(x) = f ′(x) f(x) , namely the logarithmic derivative of f, then condition (iii), given (ii), can be expressed in the simpler form g′(x) g(x)2 ln g(x) → 0. The above three conditions are satisfied by a variety of commonly used functions. For example, (i) f(x) = xp, p > 0, (ii) f(x) = exp(xr), 0 < r < 1,

• r convex combinations of products of such functions.

37

Theorem 2. If α = {1/f(j)}∞

j=1, where f satisfies (i)–(iii), then, as N → ∞,

EN(α) ∼ F(N) := −f(N) ln f ′(N) f(N)

• .

. Proof of Th. 2. Since F(N) > 0, (10) can be written as EN(α) = F(N) ∞

• 1 − exp

N

• j=1

ln

• 1 − e− F (N)

f(j) s

ds. 38

Lemma 2. Suppose f satisfies (i)–(iii) and F(x) := −f(x) ln f ′(x) f(x)

• .

(19) Then lim

N N

• j=1

ln

• 1 − e− F (N)

f(j) s

= −∞, if s < 1; 0, if s ≥ 1. (20) 39

Proof of the lemma. Since F(N) f(j) → ∞, as N → ∞, and ln(1 − x) ∼ −x, as x → 0, we can replace

N

• j=1

ln

• 1 − e− F (N)

f(j) s

by −

N

• j=1

e− F (N)

f(j) s.

Furthermore, f is increasing, N

1

e− F (N)

f(x) sdx ≤

N

• j=1

e− F (N)

f(j) s ≤

N+1

1

e− F (N)

f(x) sdx,

hence we can work with N

1

e− F (N)

f(x) sdx.

40

First, integration by parts gives N

1

e− F (N)

f(x) sdx =

• f(x)2e− F (N)

f(x) s

sF(N)f ′(x) N

x=1

− N

1

e− F (N)

f(x) s

sF(N) f(x)2 f ′(x) ′ dx. (21) Now, N

1

e− F (N)

f(x) s

sF(N) f(x)2 f ′(x) ′ dx = 2 s N

1

f(x) F(N)e− F (N)

f(x) sdx

−1 s N

1

f ′′(x)/f ′(x) f ′(x)/f(x) f(x) F(N)e− F (N)

f(x) sdx,

thus the conditions (i), (ii) and (iii) and the definition of F give that, as N → ∞, 41

N

1

e− F (N)

f(x) s

sF(N) f(x)2 f ′(x) ′ dx = o N

1

e− F (N)

f(x) sdx

• .

The role of condition (iii) is to control the term in which f ′′(x) appears. Then, (21) implies N

1

e− F (N)

f(x) sdx ∼ f(N)2e− F (N) f(N) s

sF(N)f ′(N) . Using again (19), i.e. the definition of F, we get N

1

e− F (N)

f(x) sdx ∼

1 s ln [f(N)/f ′(N)] f(N) f ′(N) 1−s . (22) Since we assumed f(N)/f ′(N) ր ∞, (22) implies that lim

N

N

1

e− F (N)

f(x) sdx =

∞, if s < 1; 0, if s ≥ 1. (23) 42

Using Lemma 2 we obtain that as N → ∞, EN(α) = F(N)[1 + o(1)]+ F(N) ∞

1

• 1 − exp

N

• j=1

ln

• 1 − e− F (N)

f(j) s

ds. (24) Now, by the DCT lim

N

∞

1

• 1 − exp
• − C1−s

N

s ln CN

• ds = 0,

where CN = f(N)/f ′(N) → ∞. Using (22) this implies that lim

N

∞

1

• 1 − exp
• −

N

1

e− F (N)

f(x) sdx

• ds = 0.

By similar inequalities as in the proof of the lemma we once again can argue that the limit is valid for the sum as well, namely lim

N

∞

1

• 1 − exp

N

• j=1

ln

• 1 − e− F (N)

f(j) s

ds = 0, and the assertion of the theorem follows immediately from (24).

• 43

It is plausible that the result of Theorem 2 also applies to sequences {aj}∞

j=1 other than

those satisfying the assumptions (i)–(iii). For sequences to which the result of the theorem cannot be shown to apply, we can obtain upper and lower bounds of EN(α) (often sharp) by using Lemma 1. In particular, if ǫj = e−ωj, where ωj ր ∞, then by Lemma 1 EN(α) ≤ ǫN N

• j=1

xaj

• dx

x − ln ǫN =

N

• j=1

e−ωNaj aj + ωN. (25) In order to best utilize the estimate given in (25), we should find a sequence {ωj}∞

j=1 with a

growth rate such that

N

• n=1

e−ωNaj aj =

N

• j=1

e− ln aj−ωNaj = O(ωN). Such a sequence {ωj}∞

j=1 would yield

EN(α) = O (ωN) . (26) In many cases, we will even have

N

• j=1

e− ln aj−ωNaj ∼ CωN, for some C = 0. (27) 44

To illustrate this idea consider a sequence of the type aj = e−jcj, where {cj}∞

j=1 is any

nondecreasing sequence of positive numbers. It is easy to see that such sequences {aj}∞

j=1 do

not satisfy all the conditions (i)–(iii). However, it follows immediately that ωj = 1 aj = ejcj satisfies (27). By using Lemma 1, it may be concluded that EN(α) ∼ CωN = CeNcN, for some C = 0. (28) A specific situation where (28) holds will be illustrated in Example 4 below. 45

The next theorem helps us obtain asymptotic estimates by comparison with sequences α for which the asymptotic estimates of EN(α) are known. First, we recall the following notation. Suppose that {sj}∞

j=1 and {tj}∞ j=1 are two sequences of nonnegative terms. The symbol

sj ≍ tj means that there are two constants C > c > 0 and an integer n0 > 0 such that ctj ≤ sj ≤ Ctj, for all j ≥ j0, i.e. sj = O(tj) and tj = O(sj). 46

Theorem 3. Let α = {aj}∞

j=1 and β = {bj}∞ j=1 be sequences of strictly positive terms

such that limN EN(α) = limN EN(β) = ∞. (i) If there exists an j0 such that aj = bj, for all j ≥ j0, then EN(β) − EN(α) is bounded; (ii) if aj = O(bj), then EN(β) = O (EN(α)), as N → ∞; (iii) if aj = o(bj), then EN(β) = o (EN(α)), as N → ∞; (iv) if aj ≍ bj, then EN(β) ≍ EN(α), as N → ∞; (v) if aj ∼ bj, then EN(β) ∼ EN(α), as N → ∞. 47

• Proof. (i) follows easily from (10), whereas (ii)–(v) follow from (i) and (12) . For instance,

to prove (v), we first fix an ǫ > 0. Then (1 − ǫ)bj ≤ aj ≤ (1 + ǫ)bj, for all j ≥ j0(ǫ). Thus, by part (i) there is an M = M(ǫ) such that 1 1 + ǫEN(β) − M ≤ EN(α) ≤ 1 1 − ǫEN(β) + M, for all N ≥ N0(ǫ). If we divide by EN(β) and then let N → ∞, we obtain (v) since ǫ is arbitrary and limN EN(β) = ∞.

• 48

EXAMPLES Example 1. The case aj = 1, for all j, has been already analyzed in detail. This case can also provide us with an application of Theorem 3: If β = {bn}∞

n=1 is a sequence such that 0 < limbn ≤ limbn < ∞, then by part (iv) of

Theorem 4, EN(β) ≍ ln N. If, in addition, lim bn = b exists, then by part (v) of Theorem 4, EN(β) ∼ b−1 ln N. 49

Example 2. aj = jp, p > 0. In this case Lp

def

= L1(α) = 1

• 1 −

∞

• j=1

(1 − xjp)

• dx

x . (Notice that Lp decreases with p). By Theorem 2 Lp < ∞. Now, in accordance with (12), we also need to estimate AN. In this case AN =

N

• n=1

np = N p+1 p + 1

• 1 + O

1 N

• by the Euler-Maclaurin sum formula. In fact this formula gives the full asymptotic expansion
• f AN. It even gives AN exactly as a polynomial in N of degree p+1, if p is a positive integer.

Therefore (12) implies E[TN] = ANEN(α) = N p+1 p + 1Lp [1 + o(1)] . 50

The case p = 1 is of particular interest. From Euler’s pentagonal-number formula, namely

∞

• j=1

(1 − xj) = 1 +

∞

• k=1

(−1)k xω(k) + xω(−k) , where ω(l) = 3l2 − l 2 , l = 0, ±1, ±2, ... (see, e.g., T.M. Apostol [1], Th. 14.3) we obtain L1 = 1 ∞

• k=1

(−1)k+1 xω(k) + xω(−k)

• dx

x =

∞

• k=1

(−1)k+1

• 1

ω(k) + 1 ω(−k)

• =

∞

• k=1

12(−1)k+1 9k2 − 1 = 4π √ 3 3 − 6. 51

Furthermore we can derive an upper bound for the error δN = L1 − EN(α), as follows. For any λ ∈ (0, 1), we have δN = 1−N−λ N

• j=1

(1 − xj) −

∞

• j=1

(1 − xj)

• dx

x + 1

1−N−λ N

• j=1

(1 − xj)

• 1 −

∞

• j=N+1

(1 − xj)

• dx

x ≤ 1−N−λ N

• j=1

(1 − xj) − 1 −

N

• k=1

(−1)k xω(k) + xω(−k)

• dx

x + 1−N−λ

• ∞
• k=N+1

(−1)k xω(k) + xω(−k)

• dx

x + 1

1−N−λ

dx x = 1−N−λ N

• j=1

(1 − xj) − 1 −

N

• k=1

(−1)k xω(k) + xω(−k)

• dx

x + O(N −λ). The integrand in the last integral is bounded by NxN+1 (see the proof of Th. 14.3 in [1]). Hence 52

δN ≤ N N + 1

• 1 − 1

N λ N + O(N −λ) = O(N −λ), for any λ ∈ (0, 1). Notice that the general error bound given in (18) does not cover this case. Furthermore, nu- merical evidence suggests that the order of the error δN is even smaller than N −λ, λ < 1. In conclusion we have E[TN] =

• 2π

√ 3 3 − 3

• N(N + 1)
• 1 + O(N −λ)
• ,

for any λ < 1. (29) 53

Example 3. aj = 1/jp, p > 0. In this case Theorem 1 implies L1(α) = ∞. If f(x) = xp, then f satisfies the assumptions (i)–(iii) and hence Theorem 2 applies. Therefore EN(α) ∼ N p ln N. On the other hand AN ∼ N 1−p 1 − p, if 0 < p < 1, AN = HN ∼ ln N, if p = 1 and AN ∼ ζ(p), if p > 1, where ζ(p) denotes the Riemann zeta function (the above asymptotics follow from the Euler- Maclaurin sum formula). Thus E[TN] ∼ N 1 − p ln N, if 0 < p < 1, E[TN] ∼ N(ln N)2, if p = 1 and E[TN] ∼ ζ(p)N p ln N, if p > 1. Note that p = 1 corresponds to the so-called Zipf distribution, and indeed the asymptotic estimate here is in agreement with the known result for the Zipf distribution (see P. Flajolet et al, [10]). 54

Example 4. aj = epj, bj = e−pj, p > 0. For the sequence α = {aj}∞

n=0 we have, by

Theorem 2, that L1(α) < ∞. Furthermore, by (18) δN = EN(α) − L1(α) ≤

∞

• j=N+1

e−pj = e−p(N+1) 1 − e−p . Also, AN = ep(N+1) − 1 ep − 1 . Therefore, by (12) E[TN] = ep(N+1) ep − 1 L1(α) + O(1). (30) In the special case of aj = 2j, we have (see J.M. Borwein and P.B. Borwein [5])

∞

• j=0

(1 − x2j) =

∞

• k=0

(−1)δ(k)xk, where δ(k) is the number of ones in the binary expansion of k. Therefore, by (13) L1(α) = 1

• 1 −

∞

• k=0

(−1)δ(k)xk

• dx

x =

∞

• k=1

(−1)δ(k)−1 k . 55

Now, for the sequence β = {bj}∞

j=0 we have L(β) = ∞. Furthermore f(x) = e−px does

not satisfy all conditions (i)–(iii), thus we cannot use Theorem 3. However BN =

N

• j=0

bj = 1 − e−p(N+1) 1 − e−p and EN(β) = EN(e−pNα) = epNEN(α) = epN L1(α) + O(e−pN)

• = epNL1(α) + O(1).

Therefore, for the sequence β, the estimate of E[TN] will be the same as in (30). In fact, the sequences α and β produce the same coupon probabilities. This follows from the fact that for each N, if we let cN = e−pN, then {aj : 0 ≤ j ≤ N} = {cNbj : 0 ≤ j ≤ N}, i.e. the elements of the two truncated sequences are proportional to each others. Therefore, α and β lead to the same asymptotics. The sequence β was mentioned in order to illustrate (28) and to exemplify a situation not covered by Theorem 2. 56

Example 5. If α is of the type aj = ejcj, where cj ր ∞, then L1(α) < ∞ and AN =

N

• j=1

aj ∼ aN. (31) Examples of such sequences are aj = ejr with r > 1, aj = jj and aj = j!. For instance, if aj = j!, then AN ∼ N! and E[TN] ∼ N!L1(α). Example 6. aj = 1/(j!). Notice that f(x) = Γ(x + 1) does not satisfy the conditions (i)–(iii). However, (26) is valid with ωj = j!, thus EN(α) = O(N!). If (28) is true, then EN(α) ∼ CN!, in which case E[TN] ∼ CeN!. 57

References

[1] T.M. Apostol, Introduction to Analytic Number Theory (Springer-Verlag, New York, 1976). [2] L.E. Baum and P. Billingsley, Asymptotic distributions for the Coupon Collector’s prob- lem, Annals of Mathematical Statistics 36 (1965) 1835–1839. [3] C.M. Bender and S.A. Orszag, Advanced Mathematical Methods for Scientists and En- gineers (McGraw-Hill, New York, 1978). [4] S. Boneh and V.G. Papanicolaou, General Asymptotic Estimates for the Coupon Col- lector Problem, Journal of Computational and Applied Mathematics 67 (2) (Mar. 1996) 277–289. [5] J.M. Borwein and P.B. Borwein, Strange series and high precision fraud, American Mathematical Monthly 99 (7) (1992) 622–640. [6] R.K. Brayton, On the asymptotic behavior of the number of trials necessary to complete a set with random selection, Journal of Mathematical Analysis and Applications 7 (1963) 31–61. [7] F.N. David and D.E. Barton, Combinatorial Chance (Charles Griffin & Co., London, 1962). [8] R. Durrett, Probability: Theory and Examples (Third Edition, Duxbury Advanced Series, Brooks/Cole—Thomson Learning. Belmont, CA, USA, 2005). 58

[9] W. Feller, An Introduction to Probability Theory and Its Applications, Vol. II, John Wiley & Sons, Inc., New York, 1966. [10] P. Flajolet, D. Gardy and L. Thimonier, Birthday paradox, coupon collectors, caching algorithms and self-organizing search, Discrete Applied Mathematics 39 (1992) 207– 229. [11] M.V. Hildebrand, The Birthday problem, a notice appeared in the American Mathemati- cal Monthly 100 (7) (1993) p. 643. [12] L. Holst, On Birthday, Collectors’, Occupancy and other classical Urn problems, Inter- national Statistical Review 54 (1986) 15–27. [13] S. Janson, Limit theorems for some sequential occupancy problems, Journal of Applied Probability 20 (1983) 545–553. [14] M. Karwan, V. Lotfi, J. Telgen, and S. Zionts, Redundancy in Mathematical Program- ming, Springer-Verlag, Berlin, 1983. [15] D.J. Newman and L. Shepp, The Double Dixie Cup problem, American Mathematical Monthly 67 (1960) 58–61. [16] S. Ross, A First Course in Probability, Seventh Edition (Pearson Prentice Hall, 2006). [17] W. Rudin, Real and Complex Analysis (McGraw-Hill, New York, 1987). [18] A. Torn and A. Zallinskas, Global Optimization, Lecture Notes in Computer Science,

• No. 350, Springer-Verlag, New York, 1989.

59