Someone on 知乎 invited me to answer this question. 中文版在此: https://www.zhihu.com/question/420632736/answer/1469911740

Problem Description

A point started at the origin, each step it can either move to left or right by one unit with equal probability.

Q: what is the probability for after m steps, this point never touch any negative point?

Solution

Classes

We are studying string like “+-+-” in this problem. Define the set of all possible strings that satisfying the condition of this problem to be \mathcal{A}, so \mathcal{A} = \{\epsilon, +, ++, +-, \dots \}. Introduce the size function for the strings to be size = s \mapsto len(s) which means we use the string length to be size function. With size function attached, we have the math structure class. Using the same size function, we define more classes:

• \mathcal{A}: all strings (including empty string) that match the problems condition (any partial leading sum is not negative)
• \mathcal{B}: all strings (including empty string) that match the problems condition and also end at the original point
• \mathcal{C}: all strings (including empty string) that match the problems condition and also end at some positive point

Construction Equations

It is straightforward to get \mathcal{A} = \mathcal{B} + \mathcal{C}. And B is known as the Catalan class. So now we only need one more equation to find all the secrete of A.

Let’s study the recurrent structure of \mathcal{A}. The element in \mathcal{A} is:

• empty string,
• or a “+” append at the end of some element from \mathcal{B}
• otherwise, is a “+” or “-” append at the end of some element from \mathcal{C}

Generating function

We can write the above as \mathcal{A} = \epsilon + \mathcal{B} * \{-\} + \mathcal{C} * \{+,-\}. Then turn all the equations in to generating function equations:

• A(z) = B(z) + C(z)
• A(z) = 1+ zB(z) + 2zC(z)

After some computation, we get the generating function for class \mathcal{A}:

Catalan class’s generating function is Catalan(z) = \frac{1-\sqrt{1-4z}}{2z}, this is a little different from the class \mathcal{B}, because the size of \mathcal{B} is doubled. Thus the length of all elements in \mathcal{B} is even. It is just insert 0 into Catalan class: 1, 0, 1, 0, 2,0, 5, 0,14,0, 42,\dots. With the above analysis, we get:

And finally the generating function of \mathcal{A}:

Compute [z^N]A(z)

During the analysis, I find the problem’s even term and odd term can be taken separately. We can adopt this idea to compute [z^N]A(Z).

It is easy to show that the first part contributes only even terms in the final result, and the second part contributes only odd terms in the final results.

For the odd terms, let’s first try to resolve the function f(z) = \frac{1-\sqrt{1-z}}{\sqrt{1-z}}

Thus for the odd part’s generating function:

The computation with {-\frac{1}{2} \choose N} is using the following formula:

We reach the final answer: [z^N]A(z) = {N \choose \lfloor N/2 \rfloor} .

A: the probability is \frac{{m \choose \lfloor m/2 \rfloor}}{2^m}

Asymptotical Results

Let’s do some asymptotical analysis of the previous result. The math tool needed here is Gamma function and Stirling formula:

• {x \choose y} = \frac{x!}{y!(x-y)!} = \frac{\Gamma(x+1)}{\Gamma(y+1)\Gamma(x-y+1)}
• \Gamma(x+1) \sim \sqrt{2\pi x}(\frac{x}{e})^x

With the above tools, we have

And the probability asymptotic result is \frac{2}{\sqrt{2\pi m}}

Python Program Simulation Results

Print some values in the sequence and search in OEIS, I find it is https://oeis.org/A001405.

Let’s write a program to print the exact values and the asymptotic results:

Python program:

from math import sqrt, pi

def asym(m):
    return 2**(m+1)/sqrt(2*pi*m)

exact_ans = [1, 2, 3, 6, 10, 20, 35, 70,
             126, 252, 462, 924, 1716, 3432,
             6435, 12870, 24310, 48620, 92378,
             184756, 352716, 705432, 1352078,
             2704156, 5200300, 10400600, 20058300,
             40116600]

asym_ans = [asym(i) for i in range(1,len(exact_ans)+1)]

m = 1
for e,a in zip(exact_ans, asym_ans):
    diff = abs(e-a)/e*100
    print("m=%d: Exact=%d, Asym=%lf, Diff=%lf%%, AsymProbability=%lf" %
          (m, e, a, diff, a/(2**(m+1))))
    m += 1