# Compute probability of seeing all the balls

Given $$\x\$$ distinguishable balls (say they have different colors), sample with replacement repeatedly until all the balls that have been seen, have been seen at least twice.

# Challenge

The input is the integer value $$\x \geq 2\$$.

The challenge is to compute 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\$$. You should output the probability as an exact fraction. For languages that don't have easy fraction arithmetic, your code can output a sum of fractions instead (e.g. $$\1/3 - 23/83 + 14/17\$$).

# Examples

If $$\x = 2\$$ then the only way to stop before seeing both balls is if you sample the same ball twice in a row. This happens with probability $$\1/2\$$.

The probability for $$\x = 3\$$ is $$\4/9\$$.

Why is this? Let us compute the probability that we don't see all the balls. There is a $$\1/3\$$ probability that the same ball is chosen in the first two steps and then we stop. If this doesn't happen, there is then $$\2/3\$$ probability that we don't select the third ball straight away. In this case we will have selected one ball twice and another one a single time. So now we want to compute the probability that we stop before seeing the third ball. We can write this as a recursive equation. Let $$\p\$$ be this probability.

$$\p = 1/3 + p/3\$$. Solving this gives us $$\p = 1/2\$$. Putting it all together we get $$\1/3 + 2/3(2/3 \cdot 1/2) = 5/9\$$. To get the desired probability that we do see all the balls we need $$\1-5/9\$$.

The probability for $$\x = 4\$$ is $$\43/96\$$.

The probability for $$\x = 5\$$ is $$\3517/7500\$$.

The probability for $$\x = 6\$$ is $$\17851/36000\$$.

• By "exact fraction," do you mean we need to output a pair of numbers: numerator and denominator? Is outputting a floating point number disallowed? Commented Dec 7, 2023 at 9:57
• May i output fractions that are not reduced? For example, 8/18 instead of 4/9?
– tsh
Commented Dec 7, 2023 at 11:59
• @KevinCruijssen worked example added
– Simd
Commented Dec 7, 2023 at 12:04
• @tsh Yes, that is fine.
– Simd
Commented Dec 7, 2023 at 12:04
• @Evargalo That seems a step too far.
– Simd
Commented Dec 8, 2023 at 7:11

# Python, 100 94 bytes

from fractions import*
P=lambda x,d=0:Fraction(d*P((d>1)*x,d-1)+x*P(x-1,d+1),d+x)if x else d>0


Attempt This Online!

x is the number of balls not sampled yet, d the number of balls sampled exactly once. Once a ball has been sampled twice it can be ignored in the further analysis, as drawing it more times effects neither the outcome nor the sampling process.

• x and ... or d>0 saves a byte Commented Dec 8, 2023 at 14:14
• -4 bytes by switching from number not sampled to total number Commented Dec 8, 2023 at 21:02
• x<1or Fraction((d>1and d*P(x,d-1))+x*P(x-1,d+1),d+x) also saves 4 bytes but in a different way.
– Neil
Commented Dec 9, 2023 at 16:58

# R, 110100959391 88 bytes

-3 bytes thanks to @Giuseppe and @pajonk

?=\(y,s=1,[=if)y[s[(y*(y-1?s+1)+s*(y?s-1))/(y+s),0],1]
f=\(x)MASS::fractions(?x-1)


s is the number of balls seen exactly once.

y is the number of balls never seen.

Try it online!

• Using Attempt this online would let you use the shorter function syntax. I think using '['='if' would also save a byte. Commented Dec 7, 2023 at 17:03
• Like this Commented Dec 8, 2023 at 16:14
• Renaming p to ? saves further 2 bytes to @Giuseppe's version Try it online! Commented Dec 8, 2023 at 19:28

# JavaScript (ES6), 81 bytes

This is based on ovs' answer.

Returns [numerator, denominator].

f=(x,d=0,q)=>[P=x?d*f(d>1&&x,d-1)[0]*f(x-1,d+1,q=Q)[1]+x*P*q:d>0,Q=x?(d+x)*q*Q:1]


Try it online!

Or with some post-processing to get the simplified fractions: Try it online!

### Method

Since we don't have native fraction support, we have to compute:

$$\left(\frac{d\cdot p_0}{q_0}+\frac{x\cdot p_1}{q_1}\right)/(d+x)=\frac{d\cdot p_0\cdot q_1+x\cdot p_1\cdot q_0}{q_0\cdot q_1(d+x)}$$

where $$\(p_0,q_0)\$$ are given by a first recursive call and $$\(p_1,q_1)\$$ are given by a second recursive call.

In the JS implementation, we save the numerator and denominator computed by the last call into the global variables P and Q respectively.

This allows us to compute the new numerator as:

//                                 q0
//                                 |
P = d * f(...)[0]   *   f(..., q = Q)[1] + x * P * q
//      \_______/       \______________/       |   |
//   [p0, q0][0] = p0   [p1, q1][1] = q1       p1  q0


and the new denominator as:

Q = (d + x) * q * Q
//            |   |
//            q0  q1


# PARI/GP, 54 bytes

p(x,d)=if(x,(d*p((d>1)*x,d-1)+x*p(x-1,d+1))/(d+x),d,1)


Attempt This Online!

A port of @ovs's Python answer.

• I am surprised pari gives you fractions by default.
– Simd
Commented Dec 8, 2023 at 18:17

# Charcoal, 63 bytes

Ｎθ⊞υ¹Ｆθ«≔⟦⟧ηＦ⁻θι⊞η▷”V5⊘？º≦kgmＥ\0⦃W‹”⟦⁺∧⊖κ×κ↨η⁰×⊕ι§υ⊕κ⁺⊕ικ⟧≔ηυ»υ


Attempt This Online! Link is to verbose version of code. Explanation: Another port of @ovs' Python answer, except using dynamic programming instead of recursion and I also replace d*P((d>1)*x,d-1) with (d-1 and d*P(x,d-1)) as that way I can make P(0,d)=1.

I wanted to remove the dependence on fractions.Fraction as even compressed that takes up 17 bytes but the best I could do was 66 bytes:

Ｎθ⊞υ¹Ｆθ«≔⟦⟧ηＦ⁻θι«Ｆ‹κ²≔Ｅ²μζ≔⟦Σ××⟦κ⊕ι⟧ζ⮌§υ⊕κ×⁺κ⊕ι⊟×ζ§υ⊕κ⟧ζ⊞ηζ»≔ηυ»Ｉυ


Attempt This Online! Link is to verbose version of code.

• Not sure if that's of any help in Charcoal, but I've added some details about what I'm doing in my answer. Commented Dec 8, 2023 at 11:21
• @Arnauld Thanks, that was very helpful, if only to find out that it still takes longer in Charcoal...
– Neil
Commented Dec 9, 2023 at 20:05

# Pyth, 36 bytes

M?G?Hs*VLV*,HG]+GH,ghGtHgtGH] B1Yg1t


Try it online!

Uses @ovs's simplification method. But since no fraction arithmetic, returns a list of fractions to add, each fraction being in the form [numerator, denominator]

### Explanation

M?G?Hs*VLV*,HG]+GH,ghGtHgtGH] B1Yg1tQ    # implicitly add Q
# implicitly assign Q = eval(input())
M                                        # define g(G,H):
?G                                      #   if g == 0:
Y        #     return an empty list
?H                                    #   if h == 0:
] B1         #     return [[1, 1]]
,ghGtHgtGH             #   [g(G+1, h-1), g(G-1, H)]
*,HG]+GH                       #   [[H, G+H], [G, G+H]]
V                               #   vectorize (apply to first of each list, then second)
L                                #   map constant pair over list returned by g
*V                                 #   vectorized multiplication
s                                   #   merge both results into one list
g1tQ    # print g(1, Q-1)


# Wolfram Language(Mathematica), 70 bytes

A port of @ovs' Python answer in Mathematica.

Golfed version. Try it online!

x_~P~d_:0:=If[x<1,Boole[d>0],(d P[Boole[d>1]x,d-1]+P[x-1,d+1]x)/(d+x)]