Skip to main content
added 181 characters in body
Source Link
Iosif Pinelis
  • 127.7k
  • 8
  • 107
  • 229

Let $S_n:=X_1+\dots+X_n$ (with $S_0:=0$), where $X_1,X_2,\dots$ are iid Bernoulli random variables with parameter $1/2$. Let then \begin{equation*} T:=\inf\{n\ge0\colon S_n>f(n)\}. \end{equation*} You are interested in $P(T=n)$.

We have \begin{align*} P(T=n)&=P(S_0\le f(0),\dots,S_{n-1}\le f(n-1),S_n>f(n)) \\ &=\frac1{2^n}\,\sum_{x_1=0}^1\dots\sum_{x_n=0}^1 I_n(x_1,\dots,x_n), \tag{*} \end{align*}\begin{align*} P(T=n)&=P(S_0\le f(0),\dots,S_{n-1}\le f(n-1),S_n>f(n)) \\ &=\frac1{2^n}\,\sum_{x_1=0}^1\dots\sum_{x_n=0}^1 I_n(x_1,\dots,x_n), \tag{$*$} \end{align*} where \begin{equation*} I_n(x_1,\dots,x_n):= 1(s_0\le f(0),\dots,s_{n-1}\le f(n-1),s_n>f(n)) \end{equation*} and $s_k:=x_1+\dots+x_k$ (with $s_0:=0$).

The $n$-fold sum in () provides an analytic expression for $P(T=n)$. However, this expression is rather complicated, and the calculation based on ($*$) provides an analytic expression for $P(T=n)$. However, this expression is rather complicated, and the calculation based on ($*$) requires an exponential in $n$ number of arithmetic operations.

It is much more effective to compute $P(T=n)$ recursively. First here, note that for natural $n$ \begin{equation*} P(T=n)=P_{n-1}-P_n,\tag{0} \end{equation*} where \begin{equation*} P_n:=P(T>n)=P(T\ge n+1). \end{equation*} In turn, \begin{equation*} P_n=\sum_{x=0}^{g(n)}p_{n,x},\tag{1} \end{equation*} where \begin{equation*} g(n):=\min(n,f(n)) \end{equation*} and \begin{equation*} p_{n,x}:=P(T>n,S_n=x). \end{equation*} For $n=1,2,\dots$ \begin{align*} p_{n,x}&=\sum_{y=0}^\infty P(T>n-1,S_{n-1}=y,S_n=x)1(x\le f(n)) \\ &=P(T>n-1,S_{n-1}=x,S_n=x)1(x\le f(n)) \\ &+P(T>n-1,S_{n-1}=x-1,S_n=x)1(x\le f(n)) \\ &=P(T>n-1,S_{n-1}=x,X_n=0)1(x\le f(n)) \\ &+P(T>n-1,S_{n-1}=x-1,X_n=1)1(x\le f(n)) \\ &=P(T>n-1,S_{n-1}=x)P(X_n=0)1(x\le f(n)) \\ &+P(T>n-1,S_{n-1}=x-1)P(X_n=1)1(x\le f(n)) \\ &=\frac{p_{n-1,x}+p_{n-1,x-1}}2\,1(x\le f(n)). \end{align*} Thus, we have a recursive difference scheme to determine $p_{n,x}$: for all $x=0,1,\dots$, \begin{equation*} p_{0,x}=1(x\le f(0),x=0) \end{equation*} and, for all $n=1,2,\dots$, \begin{align*} p_{n,x}=\frac{p_{n-1,x}+p_{n-1,x-1}}2\,1(x\le f(n)). \end{align*} So, to compute all nonzero values of $p_{n,x}$ for all $x=0,1,\dots$ by this scheme, we only need $O(\sum_{k=0}^n g(n))=O(n^2)$ arithmetic operations. Having computed the values of $p_{n,x}$, we use (0) and (1) to quickly finish the calculation of the probability of interest, $P_T(n):=P(T=n)$. This way, Mathematica computes $P_T(1),\dots,P_T(100)$ (for $f(n)\equiv n/2+\sqrt n$) in about 0.11 sec, and then $P_T(1),\dots,P_T(200)$ in about 0.37 sec:

enter image description here

In contrast, using formula (*$*$), Mathematica takes about 0.56 sec to compute just $P_T(1),\dots,P_T(10)$ (for the same $f$) and about 3.2 sec to compute $P_T(1),\dots,P_T(12)$.

Let $S_n:=X_1+\dots+X_n$ (with $S_0:=0$), where $X_1,X_2,\dots$ are iid Bernoulli random variables with parameter $1/2$. Let then \begin{equation*} T:=\inf\{n\ge0\colon S_n>f(n)\}. \end{equation*} You are interested in $P(T=n)$.

We have \begin{align*} P(T=n)&=P(S_0\le f(0),\dots,S_{n-1}\le f(n-1),S_n>f(n)) \\ &=\frac1{2^n}\,\sum_{x_1=0}^1\dots\sum_{x_n=0}^1 I_n(x_1,\dots,x_n), \tag{*} \end{align*} where \begin{equation*} I_n(x_1,\dots,x_n):= 1(s_0\le f(0),\dots,s_{n-1}\le f(n-1),s_n>f(n)) \end{equation*} and $s_k:=x_1+\dots+x_k$ (with $s_0:=0$).

The $n$-fold sum in () provides an analytic expression for $P(T=n)$. However, this expression is rather complicated, and the calculation based on () requires an exponential in $n$ number of arithmetic operations.

It is much more effective to compute $P(T=n)$ recursively. First here, note that for natural $n$ \begin{equation*} P(T=n)=P_{n-1}-P_n,\tag{0} \end{equation*} where \begin{equation*} P_n:=P(T>n)=P(T\ge n+1). \end{equation*} In turn, \begin{equation*} P_n=\sum_{x=0}^{g(n)}p_{n,x},\tag{1} \end{equation*} where \begin{equation*} g(n):=\min(n,f(n)) \end{equation*} and \begin{equation*} p_{n,x}:=P(T>n,S_n=x). \end{equation*} For $n=1,2,\dots$ \begin{align*} p_{n,x}&=\sum_{y=0}^\infty P(T>n-1,S_{n-1}=y,S_n=x)1(x\le f(n)) \\ &=P(T>n-1,S_{n-1}=x,S_n=x)1(x\le f(n)) \\ &+P(T>n-1,S_{n-1}=x-1,S_n=x)1(x\le f(n)) \\ &=P(T>n-1,S_{n-1}=x,X_n=0)1(x\le f(n)) \\ &+P(T>n-1,S_{n-1}=x-1,X_n=1)1(x\le f(n)) \\ &=P(T>n-1,S_{n-1}=x)P(X_n=0)1(x\le f(n)) \\ &+P(T>n-1,S_{n-1}=x-1)P(X_n=1)1(x\le f(n)) \\ &=\frac{p_{n-1,x}+p_{n-1,x-1}}2\,1(x\le f(n)). \end{align*} Thus, we have a recursive difference scheme to determine $p_{n,x}$: for all $x=0,1,\dots$, \begin{equation*} p_{0,x}=1(x\le f(0),x=0) \end{equation*} and, for all $n=1,2,\dots$, \begin{align*} p_{n,x}=\frac{p_{n-1,x}+p_{n-1,x-1}}2\,1(x\le f(n)). \end{align*} So, to compute all nonzero values of $p_{n,x}$ for all $x=0,1,\dots$ by this scheme, we only need $O(\sum_{k=0}^n g(n))=O(n^2)$ arithmetic operations. Having computed the values of $p_{n,x}$, we use (0) and (1) to quickly finish the calculation of the probability of interest, $P_T(n):=P(T=n)$. This way, Mathematica computes $P_T(1),\dots,P_T(100)$ (for $f(n)\equiv n/2+\sqrt n$) in about 0.11 sec, and then $P_T(1),\dots,P_T(200)$ in about 0.37 sec:

enter image description here

In contrast, using formula (*), Mathematica takes about 0.56 sec to compute just $P_T(1),\dots,P_T(10)$ (for the same $f$) and about 3.2 sec to compute $P_T(1),\dots,P_T(12)$.

Let $S_n:=X_1+\dots+X_n$ (with $S_0:=0$), where $X_1,X_2,\dots$ are iid Bernoulli random variables with parameter $1/2$. Let then \begin{equation*} T:=\inf\{n\ge0\colon S_n>f(n)\}. \end{equation*} You are interested in $P(T=n)$.

We have \begin{align*} P(T=n)&=P(S_0\le f(0),\dots,S_{n-1}\le f(n-1),S_n>f(n)) \\ &=\frac1{2^n}\,\sum_{x_1=0}^1\dots\sum_{x_n=0}^1 I_n(x_1,\dots,x_n), \tag{$*$} \end{align*} where \begin{equation*} I_n(x_1,\dots,x_n):= 1(s_0\le f(0),\dots,s_{n-1}\le f(n-1),s_n>f(n)) \end{equation*} and $s_k:=x_1+\dots+x_k$ (with $s_0:=0$).

The $n$-fold sum in ($*$) provides an analytic expression for $P(T=n)$. However, this expression is rather complicated, and the calculation based on ($*$) requires an exponential in $n$ number of arithmetic operations.

It is much more effective to compute $P(T=n)$ recursively. First here, note that for natural $n$ \begin{equation*} P(T=n)=P_{n-1}-P_n,\tag{0} \end{equation*} where \begin{equation*} P_n:=P(T>n)=P(T\ge n+1). \end{equation*} In turn, \begin{equation*} P_n=\sum_{x=0}^{g(n)}p_{n,x},\tag{1} \end{equation*} where \begin{equation*} g(n):=\min(n,f(n)) \end{equation*} and \begin{equation*} p_{n,x}:=P(T>n,S_n=x). \end{equation*} For $n=1,2,\dots$ \begin{align*} p_{n,x}&=\sum_{y=0}^\infty P(T>n-1,S_{n-1}=y,S_n=x)1(x\le f(n)) \\ &=P(T>n-1,S_{n-1}=x,S_n=x)1(x\le f(n)) \\ &+P(T>n-1,S_{n-1}=x-1,S_n=x)1(x\le f(n)) \\ &=P(T>n-1,S_{n-1}=x,X_n=0)1(x\le f(n)) \\ &+P(T>n-1,S_{n-1}=x-1,X_n=1)1(x\le f(n)) \\ &=P(T>n-1,S_{n-1}=x)P(X_n=0)1(x\le f(n)) \\ &+P(T>n-1,S_{n-1}=x-1)P(X_n=1)1(x\le f(n)) \\ &=\frac{p_{n-1,x}+p_{n-1,x-1}}2\,1(x\le f(n)). \end{align*} Thus, we have a recursive difference scheme to determine $p_{n,x}$: for all $x=0,1,\dots$, \begin{equation*} p_{0,x}=1(x\le f(0),x=0) \end{equation*} and, for all $n=1,2,\dots$, \begin{align*} p_{n,x}=\frac{p_{n-1,x}+p_{n-1,x-1}}2\,1(x\le f(n)). \end{align*} So, to compute all nonzero values of $p_{n,x}$ for all $x=0,1,\dots$ by this scheme, we only need $O(\sum_{k=0}^n g(n))=O(n^2)$ arithmetic operations. Having computed the values of $p_{n,x}$, we use (0) and (1) to quickly finish the calculation of the probability of interest, $P_T(n):=P(T=n)$. This way, Mathematica computes $P_T(1),\dots,P_T(100)$ (for $f(n)\equiv n/2+\sqrt n$) in about 0.11 sec, and then $P_T(1),\dots,P_T(200)$ in about 0.37 sec:

enter image description here

In contrast, using formula ($*$), Mathematica takes about 0.56 sec to compute just $P_T(1),\dots,P_T(10)$ (for the same $f$) and about 3.2 sec to compute $P_T(1),\dots,P_T(12)$.

added 181 characters in body
Source Link
Iosif Pinelis
  • 127.7k
  • 8
  • 107
  • 229

Let $S_n:=X_1+\dots+X_n$ (with $S_0:=0$), where $X_1,X_2,\dots$ are iid Bernoulli random variables with parameter $1/2$. Let then \begin{equation*} T:=\inf\{n\ge0\colon S_n>f(n)\}. \end{equation*} You are interested in $P(T=n)$.

We have \begin{align*} P(T=n)&=P(S_0\le f(0),\dots,S_{n-1}\le f(n-1),S_n>f(n)) \\ &=\frac1{2^n}\,\sum_{x_1=0}^1\dots\sum_{x_n=0}^1 I_n(x_1,\dots,x_n), \tag{*} \end{align*} where \begin{equation*} I_n(x_1,\dots,x_n):= 1(s_0\le f(0),\dots,s_{n-1}\le f(n-1),s_n>f(n)) \end{equation*} and $s_k:=x_1+\dots+x_k$ (with $s_0:=0$).

The $n$-fold sum in () provides an analytic expression for $P(T=n)$. However, this expression is rather complicated, and the calculation based on () requires an exponential in $n$ number of arithmetic operations.

It is much more effective to compute $P(T=n)$ recursively. First here, note that for natural $n$ \begin{equation*} P(T=n)=P_{n-1}-P_n,\tag{0} \end{equation*} where \begin{equation*} P_n:=P(T>n)=P(T\ge n+1). \end{equation*} In turn, \begin{equation*} P_n=\sum_{x=0}^{g(n)}p_{n,x},\tag{1} \end{equation*} where \begin{equation*} g(n):=\min(n,f(n)) \end{equation*} and \begin{equation*} p_{n,x}:=P(T>n,S_n=x). \end{equation*} For $n=1,2,\dots$ \begin{align*} p_{n,x}&=\sum_{y=0}^\infty P(T>n-1,S_{n-1}=y,S_n=x)1(x\le f(n)) \\ &=P(T>n-1,S_{n-1}=x,S_n=x)1(x\le f(n)) \\ &+P(T>n-1,S_{n-1}=x-1,S_n=x)1(x\le f(n)) \\ &=P(T>n-1,S_{n-1}=x,X_n=0)1(x\le f(n)) \\ &+P(T>n-1,S_{n-1}=x-1,X_n=1)1(x\le f(n)) \\ &=P(T>n-1,S_{n-1}=x)P(X_n=0)1(x\le f(n)) \\ &+P(T>n-1,S_{n-1}=x-1)P(X_n=1)1(x\le f(n)) \\ &=\frac{p_{n-1,x}+p_{n-1,x-1}}2\,1(x\le f(n)). \end{align*} Thus, we have a recursive difference scheme to determine $p_{n,x}$: for all $x=0,1,\dots$, \begin{equation*} p_{0,x}=1(x\le f(0),x=0) \end{equation*} and, for all $n=1,2,\dots$, \begin{align*} p_{n,x}=\frac{p_{n-1,x}+p_{n-1,x-1}}2\,1(x\le f(n)). \end{align*} So, to compute all nonzero values of $p_{n,x}$ for all $x=0,1,\dots$ by this scheme, we only need $O(\sum_{k=0}^n g(n))=O(n^2)$ arithmetic operations. Having computed the values of $p_{n,x}$, we use (0) and (1) to quickly finish the calculation of the probability of interest, $P_T(n):=P(T=n)$. This way, Mathematica computes $P_T(1),\dots,P_T(100)$ (for $f(n)\equiv n/2+\sqrt n$) in about 0.11 sec, and then $P_T(1),\dots,P_T(200)$ in about 0.37 sec:

enter image description here

In contrast, using formula (*), Mathematica takes about 0.56 sec to compute just $P_T(1),\dots,P_T(10)$ (for the same $f$) and about 3.2 sec to compute $P_T(1),\dots,P_T(12)$.

Let $S_n:=X_1+\dots+X_n$ (with $S_0:=0$), where $X_1,X_2,\dots$ are iid Bernoulli random variables with parameter $1/2$. Let then \begin{equation*} T:=\inf\{n\ge0\colon S_n>f(n)\}. \end{equation*} You are interested in $P(T=n)$.

We have \begin{align*} P(T=n)&=P(S_0\le f(0),\dots,S_{n-1}\le f(n-1),S_n>f(n)) \\ &=\frac1{2^n}\,\sum_{x_1=0}^1\dots\sum_{x_n=0}^1 I_n(x_1,\dots,x_n), \tag{*} \end{align*} where \begin{equation*} I_n(x_1,\dots,x_n):= 1(s_0\le f(0),\dots,s_{n-1}\le f(n-1),s_n>f(n)) \end{equation*} and $s_k:=x_1+\dots+x_k$ (with $s_0:=0$).

The $n$-fold sum in () provides an analytic expression for $P(T=n)$. However, this expression is rather complicated, and the calculation based on () requires an exponential in $n$ number of arithmetic operations.

It is much more effective to compute $P(T=n)$ recursively. First here, note that for natural $n$ \begin{equation*} P(T=n)=P_{n-1}-P_n,\tag{0} \end{equation*} where \begin{equation*} P_n:=P(T>n)=P(T\ge n+1). \end{equation*} In turn, \begin{equation*} P_n=\sum_{x=0}^{g(n)}p_{n,x},\tag{1} \end{equation*} where \begin{equation*} g(n):=\min(n,f(n)) \end{equation*} and \begin{equation*} p_{n,x}:=P(T>n,S_n=x). \end{equation*} For $n=1,2,\dots$ \begin{align*} p_{n,x}&=\sum_{y=0}^\infty P(T>n-1,S_{n-1}=y,S_n=x)1(x\le f(n)) \\ &=P(T>n-1,S_{n-1}=x,S_n=x)1(x\le f(n)) \\ &+P(T>n-1,S_{n-1}=x-1,S_n=x)1(x\le f(n)) \\ &=P(T>n-1,S_{n-1}=x,X_n=0)1(x\le f(n)) \\ &+P(T>n-1,S_{n-1}=x-1,X_n=1)1(x\le f(n)) \\ &=P(T>n-1,S_{n-1}=x)P(X_n=0)1(x\le f(n)) \\ &+P(T>n-1,S_{n-1}=x-1)P(X_n=1)1(x\le f(n)) \\ &=\frac{p_{n-1,x}+p_{n-1,x-1}}2\,1(x\le f(n)). \end{align*} Thus, we have a recursive difference scheme to determine $p_{n,x}$: for all $x=0,1,\dots$, \begin{equation*} p_{0,x}=1(x\le f(0),x=0) \end{equation*} and, for all $n=1,2,\dots$, \begin{align*} p_{n,x}=\frac{p_{n-1,x}+p_{n-1,x-1}}2\,1(x\le f(n)). \end{align*} So, to compute all nonzero values of $p_{n,x}$ for all $x=0,1,\dots$ by this scheme, we only need $O(\sum_{k=0}^n g(n))=O(n^2)$ arithmetic operations. Having computed the values of $p_{n,x}$, we use (0) and (1) to quickly finish the calculation of the probability of interest, $P_T(n):=P(T=n)$. This way, Mathematica computes $P_T(1),\dots,P_T(100)$ (for $f(n)\equiv n/2+\sqrt n$) in about 0.11 sec, and then $P_T(1),\dots,P_T(200)$ in about 0.37 sec:

enter image description here

Let $S_n:=X_1+\dots+X_n$ (with $S_0:=0$), where $X_1,X_2,\dots$ are iid Bernoulli random variables with parameter $1/2$. Let then \begin{equation*} T:=\inf\{n\ge0\colon S_n>f(n)\}. \end{equation*} You are interested in $P(T=n)$.

We have \begin{align*} P(T=n)&=P(S_0\le f(0),\dots,S_{n-1}\le f(n-1),S_n>f(n)) \\ &=\frac1{2^n}\,\sum_{x_1=0}^1\dots\sum_{x_n=0}^1 I_n(x_1,\dots,x_n), \tag{*} \end{align*} where \begin{equation*} I_n(x_1,\dots,x_n):= 1(s_0\le f(0),\dots,s_{n-1}\le f(n-1),s_n>f(n)) \end{equation*} and $s_k:=x_1+\dots+x_k$ (with $s_0:=0$).

The $n$-fold sum in () provides an analytic expression for $P(T=n)$. However, this expression is rather complicated, and the calculation based on () requires an exponential in $n$ number of arithmetic operations.

It is much more effective to compute $P(T=n)$ recursively. First here, note that for natural $n$ \begin{equation*} P(T=n)=P_{n-1}-P_n,\tag{0} \end{equation*} where \begin{equation*} P_n:=P(T>n)=P(T\ge n+1). \end{equation*} In turn, \begin{equation*} P_n=\sum_{x=0}^{g(n)}p_{n,x},\tag{1} \end{equation*} where \begin{equation*} g(n):=\min(n,f(n)) \end{equation*} and \begin{equation*} p_{n,x}:=P(T>n,S_n=x). \end{equation*} For $n=1,2,\dots$ \begin{align*} p_{n,x}&=\sum_{y=0}^\infty P(T>n-1,S_{n-1}=y,S_n=x)1(x\le f(n)) \\ &=P(T>n-1,S_{n-1}=x,S_n=x)1(x\le f(n)) \\ &+P(T>n-1,S_{n-1}=x-1,S_n=x)1(x\le f(n)) \\ &=P(T>n-1,S_{n-1}=x,X_n=0)1(x\le f(n)) \\ &+P(T>n-1,S_{n-1}=x-1,X_n=1)1(x\le f(n)) \\ &=P(T>n-1,S_{n-1}=x)P(X_n=0)1(x\le f(n)) \\ &+P(T>n-1,S_{n-1}=x-1)P(X_n=1)1(x\le f(n)) \\ &=\frac{p_{n-1,x}+p_{n-1,x-1}}2\,1(x\le f(n)). \end{align*} Thus, we have a recursive difference scheme to determine $p_{n,x}$: for all $x=0,1,\dots$, \begin{equation*} p_{0,x}=1(x\le f(0),x=0) \end{equation*} and, for all $n=1,2,\dots$, \begin{align*} p_{n,x}=\frac{p_{n-1,x}+p_{n-1,x-1}}2\,1(x\le f(n)). \end{align*} So, to compute all nonzero values of $p_{n,x}$ for all $x=0,1,\dots$ by this scheme, we only need $O(\sum_{k=0}^n g(n))=O(n^2)$ arithmetic operations. Having computed the values of $p_{n,x}$, we use (0) and (1) to quickly finish the calculation of the probability of interest, $P_T(n):=P(T=n)$. This way, Mathematica computes $P_T(1),\dots,P_T(100)$ (for $f(n)\equiv n/2+\sqrt n$) in about 0.11 sec, and then $P_T(1),\dots,P_T(200)$ in about 0.37 sec:

enter image description here

In contrast, using formula (*), Mathematica takes about 0.56 sec to compute just $P_T(1),\dots,P_T(10)$ (for the same $f$) and about 3.2 sec to compute $P_T(1),\dots,P_T(12)$.

added 16 characters in body
Source Link
Iosif Pinelis
  • 127.7k
  • 8
  • 107
  • 229

Let $S_n:=X_1+\dots+X_n$ (with $S_0:=0$), where $X_1,X_2,\dots$ are iid Bernoulli random variables with parameter $1/2$. Let then \begin{equation*} T:=\inf\{n\ge0\colon S_n>f(n)\}. \end{equation*} You are interested in $P(T=n)$.

We have \begin{align*} P(T=n)&=P(S_0\le f(0),\dots,S_{n-1}\le f(n-1),S_n>f(n)) \\ &=\frac1{2^n}\,\sum_{x_1=0}^1\dots\sum_{x_n=0}^1 I_n(x_1,\dots,x_n), \tag{*} \end{align*} where \begin{equation*} I_n(x_1,\dots,x_n):= 1(s_0\le f(0),\dots,s_{n-1}\le f(n-1),s_n>f(n)) \end{equation*} and $s_k:=x_1+\dots+x_k$ (with $s_0:=0$).

The $n$-fold sum in (*) provides an analytic expression for $P(T=n)$. However, this expression is rather complicated, and the calculation based on (1) provides an analytic expression for $P(T=n)$. However, this expression is rather complicated, and the calculation based on () requires an exponential in $n$ number of arithmetic operations.

It is much more effective to compute $P(T=n)$ recursively. First here, note that for natural $n$ \begin{equation*} P(T=n)=P_{n-1}-P_n,\tag{0} \end{equation*} where \begin{equation*} P_n:=P(T>n)=P(T\ge n+1). \end{equation*} In turn, \begin{equation*} P_n=\sum_{x=0}^{g(n)}p_{n,x},\tag{1} \end{equation*} where \begin{equation*} g(n):=\min(n,f(n)) \end{equation*} and \begin{equation*} p_{n,x}:=P(T>n,S_n=x). \end{equation*} For $n=1,2,\dots$ \begin{align*} p_{n,x}&=\sum_{y=0}^\infty P(T>n-1,S_{n-1}=y,S_n=x)1(x\le f(n)) \\ &=P(T>n-1,S_{n-1}=x,S_n=x)1(x\le f(n)) \\ &+P(T>n-1,S_{n-1}=x-1,S_n=x)1(x\le f(n)) \\ &=P(T>n-1,S_{n-1}=x,X_n=0)1(x\le f(n)) \\ &+P(T>n-1,S_{n-1}=x-1,X_n=1)1(x\le f(n)) \\ &=P(T>n-1,S_{n-1}=x)P(X_n=0)1(x\le f(n)) \\ &+P(T>n-1,S_{n-1}=x-1)P(X_n=1)1(x\le f(n)) \\ &=\frac{p_{n-1,x}+p_{n-1,x-1}}2\,1(x\le f(n)). \end{align*} Thus, we have a recursive difference scheme to determine $p_{n,x}$: for all $x=0,1,\dots$, \begin{equation*} p_{0,x}=1(x\le f(0),x=0) \end{equation*} and, for all $n=1,2,\dots$, \begin{align*} p_{n,x}=\frac{p_{n-1,x}+p_{n-1,x-1}}2\,1(x\le f(n)). \end{align*} So, to compute all nonzero values of $p_{n,x}$ for all $x=0,1,\dots$ by this scheme, we only need $O(\sum_{k=0}^n g(n))=O(n^2)$ arithmetic operations. Having computed the values of $p_{n,x}$, we use (0) and (1) to quickly finish the calculation of the probability of interest, $P(T=n)$$P_T(n):=P(T=n)$. This way, Mathematica computes $P_0,\dots,P_{100}$$P_T(1),\dots,P_T(100)$ (for $f(n)\equiv n/2+\sqrt n$) in about 0.11 sec, and then $P_0,\dots,P_{200}$$P_T(1),\dots,P_T(200)$ in about 0.3437 sec:

enter image description hereenter image description here

Let $S_n:=X_1+\dots+X_n$ (with $S_0:=0$), where $X_1,X_2,\dots$ are iid Bernoulli random variables with parameter $1/2$. Let then \begin{equation*} T:=\inf\{n\ge0\colon S_n>f(n)\}. \end{equation*} You are interested in $P(T=n)$.

We have \begin{align*} P(T=n)&=P(S_0\le f(0),\dots,S_{n-1}\le f(n-1),S_n>f(n)) \\ &=\frac1{2^n}\,\sum_{x_1=0}^1\dots\sum_{x_n=0}^1 I_n(x_1,\dots,x_n), \tag{*} \end{align*} where \begin{equation*} I_n(x_1,\dots,x_n):= 1(s_0\le f(0),\dots,s_{n-1}\le f(n-1),s_n>f(n)) \end{equation*} and $s_k:=x_1+\dots+x_k$ (with $s_0:=0$).

The $n$-fold sum in (*) provides an analytic expression for $P(T=n)$. However, this expression is rather complicated, and the calculation based on (1) requires an exponential in $n$ number of arithmetic operations.

It is much more effective to compute $P(T=n)$ recursively. First here, note that for natural $n$ \begin{equation*} P(T=n)=P_{n-1}-P_n,\tag{0} \end{equation*} where \begin{equation*} P_n:=P(T>n)=P(T\ge n+1). \end{equation*} In turn, \begin{equation*} P_n=\sum_{x=0}^{g(n)}p_{n,x},\tag{1} \end{equation*} where \begin{equation*} g(n):=\min(n,f(n)) \end{equation*} and \begin{equation*} p_{n,x}:=P(T>n,S_n=x). \end{equation*} For $n=1,2,\dots$ \begin{align*} p_{n,x}&=\sum_{y=0}^\infty P(T>n-1,S_{n-1}=y,S_n=x)1(x\le f(n)) \\ &=P(T>n-1,S_{n-1}=x,S_n=x)1(x\le f(n)) \\ &+P(T>n-1,S_{n-1}=x-1,S_n=x)1(x\le f(n)) \\ &=P(T>n-1,S_{n-1}=x,X_n=0)1(x\le f(n)) \\ &+P(T>n-1,S_{n-1}=x-1,X_n=1)1(x\le f(n)) \\ &=P(T>n-1,S_{n-1}=x)P(X_n=0)1(x\le f(n)) \\ &+P(T>n-1,S_{n-1}=x-1)P(X_n=1)1(x\le f(n)) \\ &=\frac{p_{n-1,x}+p_{n-1,x-1}}2\,1(x\le f(n)). \end{align*} Thus, we have a recursive difference scheme to determine $p_{n,x}$: for all $x=0,1,\dots$, \begin{equation*} p_{0,x}=1(x\le f(0),x=0) \end{equation*} and, for all $n=1,2,\dots$, \begin{align*} p_{n,x}=\frac{p_{n-1,x}+p_{n-1,x-1}}2\,1(x\le f(n)). \end{align*} So, to compute all nonzero values of $p_{n,x}$ for all $x=0,1,\dots$ by this scheme, we only need $O(\sum_{k=0}^n g(n))=O(n^2)$ arithmetic operations. Having computed the values of $p_{n,x}$, we use (0) and (1) to quickly finish the calculation of the probability of interest, $P(T=n)$. This way, Mathematica computes $P_0,\dots,P_{100}$ (for $f(n)\equiv n/2+\sqrt n$) in about 0.11 sec, and then $P_0,\dots,P_{200}$ in about 0.34 sec:

enter image description here

Let $S_n:=X_1+\dots+X_n$ (with $S_0:=0$), where $X_1,X_2,\dots$ are iid Bernoulli random variables with parameter $1/2$. Let then \begin{equation*} T:=\inf\{n\ge0\colon S_n>f(n)\}. \end{equation*} You are interested in $P(T=n)$.

We have \begin{align*} P(T=n)&=P(S_0\le f(0),\dots,S_{n-1}\le f(n-1),S_n>f(n)) \\ &=\frac1{2^n}\,\sum_{x_1=0}^1\dots\sum_{x_n=0}^1 I_n(x_1,\dots,x_n), \tag{*} \end{align*} where \begin{equation*} I_n(x_1,\dots,x_n):= 1(s_0\le f(0),\dots,s_{n-1}\le f(n-1),s_n>f(n)) \end{equation*} and $s_k:=x_1+\dots+x_k$ (with $s_0:=0$).

The $n$-fold sum in () provides an analytic expression for $P(T=n)$. However, this expression is rather complicated, and the calculation based on () requires an exponential in $n$ number of arithmetic operations.

It is much more effective to compute $P(T=n)$ recursively. First here, note that for natural $n$ \begin{equation*} P(T=n)=P_{n-1}-P_n,\tag{0} \end{equation*} where \begin{equation*} P_n:=P(T>n)=P(T\ge n+1). \end{equation*} In turn, \begin{equation*} P_n=\sum_{x=0}^{g(n)}p_{n,x},\tag{1} \end{equation*} where \begin{equation*} g(n):=\min(n,f(n)) \end{equation*} and \begin{equation*} p_{n,x}:=P(T>n,S_n=x). \end{equation*} For $n=1,2,\dots$ \begin{align*} p_{n,x}&=\sum_{y=0}^\infty P(T>n-1,S_{n-1}=y,S_n=x)1(x\le f(n)) \\ &=P(T>n-1,S_{n-1}=x,S_n=x)1(x\le f(n)) \\ &+P(T>n-1,S_{n-1}=x-1,S_n=x)1(x\le f(n)) \\ &=P(T>n-1,S_{n-1}=x,X_n=0)1(x\le f(n)) \\ &+P(T>n-1,S_{n-1}=x-1,X_n=1)1(x\le f(n)) \\ &=P(T>n-1,S_{n-1}=x)P(X_n=0)1(x\le f(n)) \\ &+P(T>n-1,S_{n-1}=x-1)P(X_n=1)1(x\le f(n)) \\ &=\frac{p_{n-1,x}+p_{n-1,x-1}}2\,1(x\le f(n)). \end{align*} Thus, we have a recursive difference scheme to determine $p_{n,x}$: for all $x=0,1,\dots$, \begin{equation*} p_{0,x}=1(x\le f(0),x=0) \end{equation*} and, for all $n=1,2,\dots$, \begin{align*} p_{n,x}=\frac{p_{n-1,x}+p_{n-1,x-1}}2\,1(x\le f(n)). \end{align*} So, to compute all nonzero values of $p_{n,x}$ for all $x=0,1,\dots$ by this scheme, we only need $O(\sum_{k=0}^n g(n))=O(n^2)$ arithmetic operations. Having computed the values of $p_{n,x}$, we use (0) and (1) to quickly finish the calculation of the probability of interest, $P_T(n):=P(T=n)$. This way, Mathematica computes $P_T(1),\dots,P_T(100)$ (for $f(n)\equiv n/2+\sqrt n$) in about 0.11 sec, and then $P_T(1),\dots,P_T(200)$ in about 0.37 sec:

enter image description here

Source Link
Iosif Pinelis
  • 127.7k
  • 8
  • 107
  • 229
Loading