What are the chances of finding all the balls?

Question

Given $x$ distinguishable balls (say they have different colors), sample with replacement repeatedly until all the balls that have been sampled have been sampled at least twice. I am interested in

$$P(\text{number of distinct balls sampled} = x).$$

That is, the probability that you would have seen all the balls when you stop. You can of course stop before then, if for example you sample the same ball the first two times and $x > 1$.

@leonbloy gives the different values of $x$ from $2$ to $20$:

2   0.5
3   0.4444444
4   0.4479167
5   0.4689333
6   0.4958611
7   0.5240295
8   0.5512568
9   0.5765074
10  0.5993544
11  0.6197125
12  0.6376863
13  0.6534778
14  0.6673292
15  0.6794889
16  0.6901919
17  0.6996503
18  0.7080495
19  0.7155483
20  0.7222809

Calculating by hand, the correct probability for $x = 3$ is in fact $4/9$.

Bounty

Is it possible to compute what the probabilities are in general, either exactly or at least what they are asymptotic to, and do they tend to $1$ as $x$ tends to infinity?

Answer 1

Let $n$ be the total number of balls, and $P_{i,j}$ be the probability of success, given that we have $i$ balls that have been picked once, and $j$ balls that have not been picked.

Then we get the recursion.

$$P_{i,j} = \frac{i}{n} P_{i-1,j}+ \left(1 - \frac{i+j}{n}\right) P_{i,j} + \frac{j}{n} P_{i+1,j-1} \tag 1$$

or

$$P_{i,j} = \frac{i}{i+j} P_{i-1,j}+ \frac{j}{i+j} P_{i+1,j-1} \tag 2$$

with the boundary conditions $P_{0,j}=0$ for $j>0$ and $P_{i,0}=1$.

We are interested in $g(x)=P_{1,x-1}$

I'm not sure if $(2)$ can be solved analytically, but it certainly can be solved numerically.

Python code:

import functools
@functools.lru_cache()
def calcp(i,j):
    if j == 0:
        return 1
    if i == 0: 
        return 0
    return calcp(i-1,j)i/(i+j) + calcp(i+1,j-1)j/(i+j)
for x in range(2,41):
    p =  calcp(1,x-1)
    print(f"x={x} {p:.7f}")

First values:

x   g(x)
2   0.5
3   0.4444444
4   0.4479167
5   0.4689333
6   0.4958611
7   0.5240295
8   0.5512568
9   0.5765074
10  0.5993544
11  0.6197125
12  0.6376863
13  0.6534778
14  0.6673292
15  0.6794889
16  0.6901919
17  0.6996503
18  0.7080495
19  0.7155483
20  0.7222809

It's not obvious for me that the probability tends to one.

Edit: As @Simd points out in the comments, there is no need to use recursion to compute numerically. Here are results for greater values

1000 0.889054162308888
2000 0.898252010471831
3000 0.902866350902487
4000 0.905863651086147
5000 0.908049682762985
6000 0.909753752952075
7000 0.911140967547163
8000 0.912305225768371
9000 0.913304743087483
10000 0.91417796089266
20000 0.91947919109262
30000 0.92226664355406

Empirically, it might seem that $g(x) \sim 1 - \frac{a}{\log(x)}$

Of course, all this still does not answer the question.

Edit 2: Perhaps this paper "On Solutions to a General Combinatorial Recurrence" (Spivey 2011) can be useful, or provide some pointer.

Answer 2

A non rigorous argument:

Let $X_i$ ($i=1,2 \cdots n$) be the number of times ball $i$ was picked, in $M>0$ tries. Then $\{ X_i \}$ follow a multinomial distribution. Which, for large $n$ and $M$, is asymptotically equivalent (in some sense - see "Poissonization") to $n$ iid Poisson variables with mean $\lambda = M/n$.

Under this approximation, we have $P_1 = P(X_i=1)=\lambda e^{-\lambda}$, and $E[N_1]=n\lambda e^{-\lambda}$ where $N_1=\sum_i [X_i=1]$ is the amount of balls picked once.

Here we are not given a fixed $M$, but instead it's a random variable that corresponds to the first time where the event $E_1\equiv (N_1=0)$ occurs.

A reasonable assumption is that, in average, this should correspond to a value of $M$ for which $E[N_1] \approx 1$. That is

$$ e^\lambda = n \lambda$$

or $$ \lambda = h(n) \approx \log(n) +\log (\log(n) + \log( \log n \cdots ))) $$

where $h(x)$ is the inverse of $e^x/x$; it can be written as $h(n)=-W_{-1}(-1/n)$ where $W_{-1}$ is the secondary branch of the Lambert function [*].

The event that all balls have been picked has probability $P_0 \equiv P(X_i \ne 0 \, \forall i)= (1- e^{-\lambda})^n$ [**].

Plugging the above value of $\lambda$ we get

$$P_0 \approx 1 - n e^{-\lambda} \approx 1 - \frac{1}{\lambda} \approx 1 - \frac{1}{h(n)}\approx 1 - \frac{1}{\log(n) +\log \log(n) }$$

This lacks rigour, of course. But as the graph shows, the agreement is excelent. The lines ap0 and ap1 correspond to the approximations $h(n) \approx \log n $ and $h(n) \approx \log n + \log \log n$ resp.

---

[*] There is also a solution for small time regime, but it's irrelevant for our purpose, because the event of interest does not happen here, it's effect is asymptotically negligible, it would amount to multiply $P_0$ by a factor $(1 - \frac{1}{n} - o(n^{-2}))$; it might be a useful correction for small values of $n$, though.

[**] Actually we should use here the conditional probability $P(X_i \ne 0 | X_i \ne 1)$ . But the correction is asympotically negligible.

Answer 3

This was a very fun question. I have a plan of attack but I don't think I have the time/expertise/time to gain expertise to formally carry it through:

Let's describe the state of the experiment at any point by two numbers: $n$, the number of items that have not been seen yet, and $c$, the ratio (number of items that have been seen only once)/$n$. We fail if at some point after the first step $c$ becomes $0$. We win if $n$ becomes $0$.

Now I want to make the following 4 claims which together imply that as $x$ grows large, the probability of success goes to 1. I believe I can prove 1,2 and 3 but although 4 ought to be true I don't know how to prove it yet.

1. We won't fail at the beginning:
   As the initial $x$ goes to infinity, the probability that we fail before half the items have been seen goes to $0$.

2. A satisfactory ending:
   Let $p>0$ be any small number. If at some point we are left with the state $n,c \geq (1/(p+1))^n$, then the probability that we fail from this point on is at most $p/(1-p)$. The idea (which works assuming some approximations) is that when $n$ changes to $n-1$ the probability that the state is $n-1, c \geq (1/(p+1))^{n-1}$ is at least $1 - p^n$. By a union bound, the probability we reach the state $0, \geq 1$ is at least $1 - p - p^2 - p^3 ... = 1 - p/(1-p)$. By making $p$ as small as we want, the probability of failure here goes to $0$.

Now from here on is the part where I might be way way off. Although I am quite confident about this point #3.

3. The ratio becomes arbitrarily large:
   For any target ratio $d$, there exists an $N$ such that the probability that you reach a state $n\geq N,c\geq d$ tends to $1$ as $x$ tends to infinity. The intuition behind this is based on how the ratio changes as you go from $n$ to $n-1$. The ratio approximately increases by $(c+1)/(n-1) - t/n$ where $t$ is a random variable distributed according to the geometric distribution with parameter $1/(c+1)$ (so it has mean $c$). This isn't exact because the ratio $c$ would change in the middle but if $n$ is large enough, I hope it shouldn't matter. The expected increase in the ratio is $\approx 1/(n-1)$. With many trials and a concentration bound, I think it would follow that the ratio does increase. It's a slow increase, but the sum of $1/n$ diverges, and so does the sum of $1/cn$ (since the increase happens only once every $c$ steps). This is why I think it might be true that the ratio does get arbitrarily large.

It is worth noting that if the above point were false, the probability of finding all balls would not tend to 1. A prerequisite for the probability tending to 1 is that for any $k$, you reach the state $n=1,c \geq k$ since even from this point onwards there is a $1/k$ probability of failure. Point 3 does not tell you that you reach such a state, but combining it with points 4 and 2, it will. Before we move on to point 4, the above reasoning about point 3 also tells us why we would expect the probability of failure starting with $x$ balls to be at least $1/\ln(x)$. The ratio only increases at approximately the rate of a harmonic sum, so it only reaches $\approx \ln(x)$, and so when you have $n=1$, there's still a $1/\ln(x)$ chance of failure.

4. Going from point 3 to a state where you can use point 2:
   Fix a target $n$ and a target ratio $c = ((p+1)/p)^n$. Then there exists $N,c'$ such that for any $M>N$ if you start with the state $M,\geq c'$, then the probability of not making it to the state $n, \geq c$ goes to $0$ as $n$ increases.
   Here the problem is that $n$ is a constant, I can't assume it is large and that the niceness that you have with large numbers holds. However, the hope is that you can prove that the ratio almost surely won't decrease from $c'$ to $c$, especially as $n$ becomes larger. If point 3 were true we already have that for any $c'$, with probability $1$ you will reach a ratio larger than $c'$. So this point takes point 3 to a place where point 2 kicks in.

Sorry I can't give you closure.