The distribution function for the maximal height of $N$ non-intersecting Bessel paths

TL;DR

This paper derives explicit formulas for the distribution of the maximum height of N non-intersecting Bessel paths starting at a point and ending at zero, using orthogonal polynomial techniques.

Contribution

It provides new explicit formulas for the maximum height distribution of non-intersecting Bessel paths, linking it to Hankel determinants and orthogonal polynomials.

Findings

Distribution functions expressed via Hankel determinants.

Formulas depend on initial starting point a.

Results applicable to non-intersecting Bessel path models.

Abstract

In this paper, we consider $N$ non-intersecting Bessel paths starting at $x=a\geq 0$, and conditioned to end at the origin $x=0$. We derive the explicit formula of the distribution function for the maximum height. Depending on the starting point $a>0$ or $a=0$, the distribution functions are also given in terms of the Hankel determinants associated with the multiple discrete orthogonal polynomials or discrete orthogonal polynomials, respectively.

Equations

$$p_{abs}(t,y|x)=\frac{1}{\sqrt{2\pi t}}\left(e^{-\frac{(x-y{)}^2}{2t}}-e^{-\frac{(x+y{)}^2}{2t}}\right).$$

$$\begin{split}&\mathbb{P}(\max_{0<t<1}b_{N}^{(BE)}(t)<M)\\ &\qquad=\frac{2^{-N/2}\pi^{2N^{2}+N/2}}{M^{N(2N+1)}N!\prod_{k=0}^{N-1}(2k+1)!}\sum_{\mathbf{x}\in\mathbb{Z}^{N}}(\Delta(\mathbf{x^{2}}))^{2}\left(\prod_{j=1}^{N}x_{j}^{2}\right)\exp\biggl{\{}\frac{-\pi^{2}}{2M^{2}}\sum_{j=1}^{N}x_{j}^{2}\biggr{\}},\end{split}$$

$$\Delta (\mathbf{x})=\prod _{1\le j<k\le N}(x_k-x_j)\text{with}\mathbf{x}=(x_1,x_2,\cdots ,x_N).$$

$$p_{ref}(t,y|x)=\frac{1}{\sqrt{2\pi t}}\left(e^{-\frac{(x+y{)}^2}{2t}}+e^{-\frac{(x-y{)}^2}{2t}}\right),$$

$$\begin{split}&\mathbb{P}(\max_{0<t<1}b_{N}^{(R)}(t)<M)\\ &\qquad=\frac{2^{-N/2}\pi^{2N^{2}-3N/2}}{M^{N(2N-1)}N!\prod_{k=0}^{N-1}(2k)!}\sum_{\mathbf{x}\in\{\mathbb{Z}-1/2\}^{N}}(\Delta(\mathbf{x^{2}}))^{2}\exp\biggl{\{}\frac{-\pi^{2}}{2M^{2}}\sum_{j=1}^{N}x_{j}^{2}\biggr{\}};\end{split}$$

$$R(t)=||\mathbf{X}(t)|{|}_2=\sqrt{X_1(t{)}^2+\cdots +X_d(t{)}^2},t\ge 0,$$

$$p(t,y|x)=\frac{y}{x}\frac{1}{\sqrt{2\pi t}}\left(e^{-\frac{(x-y{)}^2}{2t}}-e^{-\frac{(x+y{)}^2}{2t}}\right),$$

$$\begin{array}{rl}& b_1(0)=b_2(0)=\cdots =b_N(0)=a,\\ & b_1(1)=b_2(1)=\cdots =b_N(1)=0,\\ & 0\le b_1(t)<b_2(t)<\cdots <b_N(t)\text{for}0<t<1.\end{array}$$

$$p(t,y|x)$$

$$p(t,y|0)$$

$$I_{\alpha }(z)=\sum _{k=0}^{\infty }\frac{(\frac{z}{2}{)}^{2k+\alpha }}{k!\Gamma (k+\alpha +1)}.$$

$$x_{1,\alpha }<x_{2,\alpha }<\cdots <x_{n,\alpha }<\cdots$$

$$\begin{array}{rl}& \mathbb{P}(\underset{0<t<1}{\max }{b}_N(t)<M)\\ & =c_N(\alpha )M^{-\frac{N(3N+2\alpha +1)}{2}}\underset{1\le i,j\le N}{\det }\left[\sum _{n=1}^{\infty }\frac{(-1{)}^{i-1}{x}_{n,\alpha }^{i+2j+\alpha -3}{J}_{\alpha }^{(i-1)}(\frac{a}{M}{x}_{n,\alpha })e^{-\frac{x_{n,\alpha }^2}{2M^2}}}{J_{\alpha +1}^2(x_{n,\alpha })}\right],\end{array}$$

$$c_N(\alpha )=\frac{2^{\frac{N(3-N)}{2}}{e}^{\frac{a^{2}N}{2}}}{a^{\frac{N(N+2\alpha -1)}{2}}{\prod }_{j=1}^N\Gamma (j)}.$$

$$\mathbb{P}(\underset{0<t<1}{\max }{b}_N(t)<M)={\tilde{c}}_N(\alpha )M^{-2N(N+\alpha )}\underset{1\le i,j\le N}{\det }\left[\sum _{n=1}^{\infty }\frac{x_{n,\alpha }^{2i+2j+2\alpha -4}{e}^{-\frac{x_{n,\alpha }^2}{2M^2}}}{J_{\alpha +1}^2(x_{n,\alpha })}\right],$$

$${\tilde{c}}_N(\alpha )=\frac{2^{2N-\alpha N-N^2}}{{\prod }_{j=1}^N\Gamma (j)\Gamma (\alpha +j)}.$$

$$(-1{)}^n{z}^n{J}_{\alpha }^{(n)}(z)=\sum _{k=0}^n{c}_{n,k}{z}^k{J}_{\alpha +k}(z)\text{for }n\in \mathbb{N},$$

$$\begin{array}{rl}\mathbb{P}(\underset{0<t<1}{\max }{b}_N(t)<M)& =\frac{2^{\frac{N(3-N)}{2}}{e}^{\frac{a^{2}N}{2}}}{a^{\frac{N(N+2\alpha -1)}{2}}{\prod }_{j=1}^N\Gamma (j)}{M}^{-\frac{N(3N+2\alpha +1)}{2}}\\ & \times \underset{1\le i,j\le N}{\det }\left[\sum _{n=1}^{\infty }{x}_{n,\alpha }^{i+2j-3}\frac{J_{\alpha +i-1}(\frac{a}{M}{x}_{n,\alpha })x_{n,\alpha }^{\alpha }{e}^{-\frac{x_{n,\alpha }^2}{2M^2}}}{J_{\alpha +1}^2(x_{n,\alpha })}\right].\end{array}$$

$$\underset{1\le i,j\le N}{\det }\left[\sum _{n=1}^{\infty }{x}_{n,\alpha }^{i+2j-3}\frac{J_{\alpha +i-1}(\frac{a}{M}{x}_{n,\alpha })x_{n,\alpha }^{\alpha }{e}^{-\frac{x_{n,\alpha }^2}{2M^2}}}{J_{\alpha +1}^2(x_{n,\alpha })}\right]\sim \underset{1\le i,j\le N}{\det }\left[\frac{(\frac{a}{2M}{)}^{\alpha +i-1}}{\Gamma (\alpha +j)}\sum _{n=1}^{\infty }{x}_{n,\alpha }^{2i+2j-4}\frac{x_{n,\alpha }^{2\alpha }{e}^{-\frac{x_{n,\alpha }^2}{2M^2}}}{J_{\alpha +1}^2(x_{n,\alpha })}\right]$$

$$\mathbb{P}(\underset{0<t<1}{\max }b(t)<M)=\frac{2^{1-\alpha }}{\Gamma (\alpha +1)}{M}^{-2(\alpha +1)}\sum _{n=1}^{\infty }\frac{x_{n,\alpha }^{2\alpha }{e}^{-\frac{x_{n,\alpha }^2}{2M^2}}}{J_{\alpha +1}^2(x_{n,\alpha })}.$$

$$w_1(x)=\frac{x^{\frac{\alpha }{2}}{J}_{\alpha }(\frac{a}{M}\sqrt{x})}{J_{\alpha +1}^2(\sqrt{x})}{e}^{-\frac{x}{2M^2}},w_2(x)$$

$$m_k^{(1)}=\sum _{x\in \{x_{n,\alpha }^2,n\in \mathbb{N}\}}{x}^k{w}_1(x),m_k^{(2)}=\sum _{x\in \{x_{n,\alpha }^2,n\in \mathbb{N}\}}{x}^k{w}_2(x).$$

$$H_N:=\det \left(\begin{array}{cccc}{m}_0^{(1)}& m_1^{(1)}& \dots & m_{N-1}^{(1)}\\ m_1^{(1)}& m_2^{(1)}& \dots & m_N^{(1)}\\ ⋮& ⋮& \ddots & ⋮\\ m_{n_1-1}^{(1)}& m_{n_1}^{(1)}& \dots & m_{N+n_1-2}^{(1)}\\ m_0^{(2)}& m_1^{(2)}& \dots & m_{N-1}^{(2)}\\ m_1^{(2)}& m_2^{(2)}& \dots & m_N^{(2)}\\ ⋮& ⋮& \ddots & ⋮\\ m_{n_2-1}^{(2)}& m_{n_2}^{(2)}& \dots & m_{N+n_2-2}^{(2)}\end{array}\right)$$

$$n_1=\{\begin{array}{ll}\frac{N}{2},& \text{mbox}ifN\text{mbox}iseven,\\ \frac{N+1}{2},& \text{mbox}ifN\text{mbox}isodd,\end{array}\text{and}{n}_2=\{\begin{array}{ll}\frac{N}{2},& \text{mbox}ifN\text{mbox}iseven,\\ \frac{N-1}{2},& \text{mbox}ifN\text{mbox}isodd.\end{array}$$

$$\tilde{w}(x)=\frac{x^{\alpha }{e}^{-\frac{x}{2M^2}}}{J_{\alpha +1}^2(\sqrt{x})},M>0,$$

$${\tilde{H}}_N:=\det \left(\begin{array}{cccc}{\tilde{m}}_0& {\tilde{m}}_1& \dots & {\tilde{m}}_{N-1}\\ {\tilde{m}}_1& {\tilde{m}}_2& \dots & {\tilde{m}}_N\\ ⋮& ⋮& \ddots & ⋮\\ {\tilde{m}}_{N-1}& {\tilde{m}}_N& \dots & {\tilde{m}}_{2N-2}\end{array}\right),$$

$${\tilde{m}}_k=\sum _{x\in \{x_{n,\alpha }^2,n\in \mathbb{N}\}}{x}^k\tilde{w}(x).$$

$$\mathbb{P}(\underset{0<t<1}{\max }{b}_N(t)<M)=c_N(\alpha )(-1{)}^{\frac{N(N-3)+2n_1^2}{4}}{M}^{-\frac{N(3N+2\alpha +1)}{2}}{H}_N,$$

$$\mathbb{P}(\underset{0<t<1}{\max }{b}_N(t)<M)={\tilde{c}}_N(\alpha )M^{-2N(N+\alpha )}{\tilde{H}}_N,$$

$$J_{-\frac{1}{2}}(x)=\sqrt{\frac{2}{\pi x}}\cos x,J_{\frac{1}{2}}(x)=\sqrt{\frac{2}{\pi x}}\sin x\text{and}{I}_{-\frac{1}{2}}(x)=\sqrt{\frac{2}{\pi x}}\cosh x;$$

$$\underset{1\le i,j\le N}{\det }\left[\sum _{n=1}^{\infty }{x}_{n,\alpha }^{2i+2j-4}\frac{x_{n,\alpha }^{2\alpha }{e}^{-\frac{x_{n,\alpha }^2}{2M^2}}}{J_{\alpha +1}^2(x_{n,\alpha })}\right]=\underset{1\le i,j\le N}{\det }\left[\sum _{n=1}^{\infty }\frac{{\pi }^{2i+2j-3}}{2}(n-\frac{1}{2}{)}^{2j+2i-4}{e}^{-\frac{(n-\frac{1}{2}{)}^2{\pi }^2}{2M^2}}\right].$$

Full text

2010 Mathematics Subject Classification. Primary 60J65, 33C47.

Keywords and phrases: non-intersecting Bessel paths; maximum distribution; multiple orthogonal polynomials; Hankel determinant.

$\ast$ Department of Mathematics, City University of Hong Kong, Hong Kong.

Email: dandai@cityu.edu.hk

${\dagger}$ Department of Mathematics, City University of Hong Kong, Hong Kong.

Email: lumingyao2-c@my.cityu.edu.hk

1 Introduction

Random walk models play an important role in various areas in physics, chemistry and computer sciences. For example, they are introduced to explain basic concepts in statistical physics [22] and stochastic algorithms [19]. When a proper diffusion scaling limit is taken, the random walk models are reduced to Brownian motion models. The standard one-dimensional Brownian motion $\{W_{t}:t\geq 0\}$ initiated at $x$ is a stochastic process with the following properties (cf. [3]):

• •

$W_{0}=x$ almost surely;

• •

$W_{t}$ is continuous almost surely;

• •

for all $0=t_{0}<t_{1}<\cdots<t_{n}$, the increments $W_{t_{i}}-W_{t_{i-1}}$ are independent and normally distributed, with $\mathbf{E}(W_{t_{i}}-W_{t_{i-1}})=0$ and $\mathbf{E}(W_{t_{i}}-W_{t_{i-1}})^{2}=t_{i}-t_{i-1}$.

There is a vast literature in the study of Brownian motions, among which people are interested in the maximal height of the outermost path in different non-intersecting Brownian motions. For example, let $0\leq b_{1}^{(BE)}(t)<b_{2}^{(BE)}(t)<\cdots<b_{N}^{(BE)}(t)$ be $N$ non-intersecting Brownian excursions, i.e. non-intersecting Brownian motions with an absorbing wall located at $x=0$, conditioned to have the same starting point $x=0$ and return to the origin at the end. The transition probability for a single Brownian motion with an absorbing wall at $x=0$, passing from $x$ to $y$ over the time interval $t$, is given by

[TABLE]

Let $M$ be a positive constant. Then, the probability that the maximal height of the outermost path $\max\limits_{0<t<1}b_{N}^{(BE)}(t)$ is less than $M$ is given explicitly as follows:

[TABLE]

where $\Delta(\mathbf{x})$ is the Vandermonde determinant

[TABLE]

The formula (1.2) was first obtained by Schehr et al. [24] through a path integral approach. Later, it was derived by using the Karlin-McGregor formula in Kobayashi et al. [14], and by lattice paths in Feierl [5]. One may also consider $N$ non-intersecting Brownian motions $0\leq b_{1}^{(R)}(t)<b_{2}^{(R)}(t)<\cdots<b_{N}^{(R)}(t)$ with a reflecting wall located at $x=0$, conditioned to have the same starting point $x=0$ and return to the origin at the end. In this case, the transition probability is

[TABLE]

and the probability for the maximal height is given by

[TABLE]

see Liechty [17, Eq. (1.18)]. It is worth mentioning that the above two probabilities (1.2) and (1.5) can be put in terms of Hankel determinants of discrete Gaussian orthogonal polynomial. Based on this observation, Liechty [17] applied the Riemann-Hilbert approach to rigorously prove that, in a proper scaling limit as $N\to\infty$, the limiting distribution of the maximal height in both cases converges to the Tracy-Widom distribution for the largest eigenvalue of the Gaussian orthogonal ensemble. This justifies the formal derivations first given by Forrester, Majumdar and Schehr [6].

In this paper, we focus on the distribution of the maximal height of the outermost path in $N$ non-intersecting Bessel paths. Recall that if $\{\mathbf{X}(t):t\geq 0\}$ is a $d$-dimensional Brownian motion, then the diffusion process

[TABLE]

is called a $d$-dimensional Bessel process or a Bessel process of order $\alpha$ with parameter $\alpha=\frac{d}{2}-1$; see [3, 10, 23]. Bessel processes have important application in financial mathematics due to their close relation to financial models such as geometric Brownian motion and Cox-Ingersoll-Ross processes; see Göing-Jaeschke and Yor [7] and references therein. When $d=1$, or equivalently $\alpha=-\frac{1}{2}$, the Bessel process $R(t)$ reduces to the Brownian motion with a reflecting wall located at the origin (1.4). While when $d=3$, or equivalently $\alpha=\frac{1}{2}$, the transition probability of $R(t)$ becomes

[TABLE]

which is closely related to the Brownian motion with an absorbing wall located at the origin (1.1); see Katori and Tanemura [12, 13]. In recent years, there has been considerable interest in the study of non-intersecting (squared) Bessel paths; for example, see [4, 11, 15, 16].

Similar to the non-intersecting Brownian motions mentioned above, we consider a model of $N$ non-intersecting Bessel paths: $\{b_{j}(t)\}_{j=1}^{N}$. All paths start at $x=a\geq 0$, remain non-negative and are conditioned to end at time $t=1$ at $x=0$. That is,

[TABLE]

The transition probability for a single particle passes form $x$ to $y$ over the time interval $t$ is given by (cf. [3, 9])

[TABLE]

for $y\geq 0$, $t>0$ and $\alpha>-1$, where $I_{\alpha}$ is the modified Bessel function of the first kind

[TABLE]

One can see from (1.6) that when $d=2(\alpha+1)$ is an integer, the Bessel process is the distance to the origin of a $d$-dimensional standard Brownian motion.

The main results of the paper are stated in the next section.

2 Statement of results

2.1 The maximal distribution

Our first result is the explicit expression of the probability of the maximal height of the outermost path $\max\limits_{0<t<1}b_{N}(t)$.

Theorem 1**.**

For $\alpha>-1$, let

[TABLE]

be the zeros of the Bessel function $J_{\alpha}(x)$. For $N\geq 1$, $M\geq a$ and the starting point $a>0$, we have

[TABLE]

where $c_{N}(\alpha)$ is a constant independent of $M$:

[TABLE]

When the starting point is the origin, i.e., $a=0$, we have

[TABLE]

where the constant $\widetilde{c}_{N}(\alpha)$ is given by

[TABLE]

Remark 1*.*

Although the expressions (2.2) and (2.4) appear different, they indeed agree with each other. To see it, from properties of the Bessel function [20, Eq. (10.6.2)], one gets the following expression for the higher order derivative of $J_{\alpha}(z)$:

[TABLE]

with the coefficient $c_{n,n}=1$. We replace $J_{\alpha}^{(i-1)}$ in (2.2) by the right-hand side of the above formula, then apply row operations to eliminate $J_{\nu}$ for $\nu<\alpha+i-1$. This gives us

[TABLE]

Since $J_{\alpha}(z)\sim\frac{z^{\alpha}}{2^{\alpha}\Gamma(\alpha+1)}$ as $z\to 0$, the determinant in the above formula behaves like

[TABLE]

as $a\to 0$. Then, (2.4) follows from the above two formulae.

Remark 2*.*

When the starting point $a$ is the origin and the number of Bessel paths $N$ is equal to 1, the distribution formula (2.4) reduces to

[TABLE]

This agrees with the results obtained by Pitman and Yor [21, Eq. (5)].

2.2 Relation to multiple orthogonal polynomials

The maximum distributions in (1.2) and (1.5) can be put in terms of Hankel determinants of discrete Gaussian orthogonal polynomial; see Liechty [17]. For the present model, we have similar results. Depending on the starting point $a>0$ or $a=0$, the distribution functions are given in terms of the Hankel determinants associated with multiple discrete orthogonal polynomials of type II or discrete orthogonal polynomials, respectively. For more properties of multiple orthogonal polynomials, one may refer to [2] and [8, Ch. 23].

In the case where the starting point $a$ is positive, the corresponding two discrete weight functions for the multiple orthogonal polynomial are

[TABLE]

which are supported on the nodes $\{x_{n,\alpha}^{2}\}_{n=1}^{\infty}$ with $x_{n,\alpha}$ being zeros of the Bessel function $J_{\alpha}(x)$ given in (2.1). The moments $m_{k}^{(j)},j=1,2$, are given by

[TABLE]

The Hankel determinant is

[TABLE]

with $N=n_{1}+n_{2}$ and

[TABLE]

In the case where the starting point is located at the origin, the discrete weight function for the orthogonal polynomial is

[TABLE]

which is also supported on the nodes $\{x_{n,\alpha}^{2}\}_{n=1}^{\infty}$. The corresponding Hankel determinant is

[TABLE]

where the moments $\widetilde{m}_{k}$ are given by

[TABLE]

Then, we have the following result.

Theorem 2**.**

Let the Hankel determinants $H_{N}$ and $\widetilde{H}_{N}$ be defined in (2.11) and (2.14), respectively. For $a>0$, we have

[TABLE]

where $c_{N}(\alpha)$ is given in (2.3). When $a=0$, we have

[TABLE]

where $\widetilde{c}_{N}(\alpha)$ is given in (2.5).

2.3 Special cases: $\alpha=\pm\frac{1}{2}$

As mentioned in Section 1, the Bessel paths reduce to the Brownian motion with a reflecting wall at the origin when $\alpha=-\frac{1}{2}$. Indeed, recall that the Bessel functions reduce to the elementary functions when $\alpha=\pm\frac{1}{2}$:

[TABLE]

see [20, Sec. 10.16]. Then, the transition probability (1.9) becomes (1.4). When the starting point $a=0$, the model (1.8) exactly reduces to the non-intersecting Brownian motion with a reflecting wall at the origin, conditioned to have the same starting point $x=0$ and return to the origin at the end. Then, the probability for the maximum height given in (2.4) also reduces to that of the Brownian motion with a reflecting wall given in (1.5). To see it, note that the zeroes in (2.1) are $x_{n,\alpha}=(n-\frac{1}{2})\pi$, $n\in\mathbb{N}$ when $\alpha=-\frac{1}{2}$. Then, the determinant in (2.4) becomes

[TABLE]

Moreover, as $\Gamma(j)\Gamma(-\frac{1}{2}+j)=2^{2-2j}\sqrt{\pi}\,\Gamma(2j-1)$, the constant in (2.5) is

[TABLE]

With the above two formulae, the probability in (2.4) reduces to

[TABLE]

which agrees with (1.5).

When $\alpha=\frac{1}{2}$, the Bessel paths are closely related to the Brownian motions with an absorbing wall at the origin; see the relation between their transition probabilities in (1.7). Although there is an additional factor $y/x$ in (1.7), the maximum distributions are exactly the same. Indeed, since the zeroes in (2.1) are $x_{n,\alpha}=n\pi$, $n\in\mathbb{Z}$ when $\alpha=\frac{1}{2}$, we adopt similar computations as above and obtain from (2.4)

[TABLE]

which agrees with (1.2).

The rest of the paper is organized as follows. In Section 3, we apply the Karlin-McGregor formula and express the probability $\mathbb{P}(\max\limits_{0<t<1}b_{N}(t)<M)$ in terms of a ratio of two determinants. Properties for some general determinants are also discussed. In Section 4, we derive the asymptotics for the determinants appearing in $\mathbb{P}(\max\limits_{0<t<1}b_{N}(t)<M)$. Then, our main results are proved with the asymptotic results. Finally, in Section 5, we present some numerical computations and conclude with a summary.

3 Some preliminary work

3.1 The Karlin-McGregor formula

Introduce the notations

[TABLE]

and the constant vectors

[TABLE]

According to the Karlin-McGregor formula in the affine Weyl alcove of height $M$ (see [13]) and the definition of $N$ non-intersecting Bessel processes, we obtain

[TABLE]

with

[TABLE]

Here, $p(t,y|x)$ is the transition probability defined in (1.9) and (1.10), and $p^{M}(t,y|x)$ is the transition probability of the particle passes from $x$ to $y$ over the time interval $t$ with an absorbing wall at the position $M$.

To find $p^{M}(t,y|x)$, let us consider the following diffusion equation

[TABLE]

see similar equations in [1, 3]. Using the method of separation of variables, one gets the unique solution to the above equation

[TABLE]

which is the transition density with respect to the speed measure $m_{\alpha}(dy)=2y^{2\alpha+1}dy$. Therefore, we obtain

[TABLE]

where $J_{\alpha}(z)$ is the Bessel function of the first kind

[TABLE]

and $x_{n,\alpha}$’s are its zeros given in (2.1).

Remark 3*.*

When $\alpha=-\frac{1}{2}$, the zeroes in (2.1) are $x_{n,-\frac{1}{2}}=(n-\frac{1}{2})\pi$, $n\in\mathbb{N}$. Then, the transition probability (3.7) becomes

[TABLE]

By the Poisson summation formula, we have

[TABLE]

It is easily seen from the above formula that $\lim\limits_{M\to\infty}p^{M}(t,y|x)=p(t,y|x)$. This, together with (3.3), gives us the desired result $\lim\limits_{M\to\infty}\mathbb{P}(\max\limits_{0<t<1}b_{N}(t)<M)=1$.

3.2 Some properties for determinants

We need some preliminary results for determinants.

Lemma 1**.**

Let $f$ and $g$ be functions with two independent variables, and $\mathbf{n}=(n_{1},n_{2},\cdots,n_{N})$ be an $N$-dimensional vector. Then, we have

[TABLE]

Proof.

From the Leibniz formula for the determinant, we have

[TABLE]

where $\sigma$ is a permutations of the set $\{1,2,\cdots,N\}$. Substituting the above formula into the left-hand side of (3.11), we get

[TABLE]

Let us change the index $n_{j}$ to $n_{\tau(j)}$, where $\tau$ is an arbitrary permutation of the set $\{1,2,\cdots,N\}$. Since the summation is taken for all $n\in\mathbb{N}^{N}$, we have

[TABLE]

Consider a new permutation $\rho$ defined as $\rho:=\tau\circ\sigma^{-1}$. Denote $l=\sigma(j)$, then $\tau(j)=\tau(\sigma^{-1}(l))=\rho(l)$. We change the index in the first product in the above formula and obtain

[TABLE]

Separating the summations about $\tau$ and $\rho$, we get

[TABLE]

With the Leibniz formula for the determinants again, we obtain (3.11) from the above formula.

This completes the proof of the lemma. ∎

Remark 4*.*

One may take matrix transpose in (3.11) and rewrite the formula as

[TABLE]

Moreover, from the proof, it is easy to verify that Lemma 1 still holds if an additional factor $\prod_{i=1}^{N}h(n_{i})$ appears in the summation. More precisely, we have

[TABLE]

To study the limits of determinants in (3.3), we will prove one more lemma below. First, let us consider the following Schur function

[TABLE]

with $\mathbf{x}$ being the vector defined in (3.1); see [18]. Here $\boldsymbol{\mu}=(\mu_{1},\mu_{2},\cdots,\mu_{N})$ is a sequence of non-negative integers in decreasing order $\mu_{1}\geq\mu_{2}\geq\cdots\geq\mu_{N}$, which is called a partition. Obviously, both the numerator and denominator are Vandermonde determinants. For example, we have

[TABLE]

Note that $s_{\boldsymbol{\mu}}(\mathbf{0})=0$ unless $\boldsymbol{\mu}=\mathbf{0}$. Moreover, we have the following property

[TABLE]

Then, we have the following result.

Lemma 2**.**

Let $f$ be a smooth function, then we have

[TABLE]

as $\mathbf{x}\to\mathbf{a}$.

Proof.

We expand $f(x_{i}y_{j})$ at $x_{i}=a$ into a Taylor series

[TABLE]

With (3.11) and the above formula, we have

[TABLE]

where $\mathbb{N}_{0}=\{0,1,2,\cdots\}$. We order the index $\mathbf{k}=(k_{1},k_{2},\cdots,k_{N})$ such that $\{k_{j}\}_{j=1}^{N}$ is a decreasing sequence. This gives us

[TABLE]

Changing the index from $\mathbf{k}$ to $\boldsymbol{\mu}$ with $\mu_{j}=k_{j}-N+j$, we get from the definition of the Schur function in (3.14) that

[TABLE]

Finally, using the asymptotics of the Schur function in (3.16), we obtain (3.17). ∎

4 Derivation of the distribution function

In this section, we will prove our main theorem by studying some properties of $q(\mathbf{x},\mathbf{y})$ and $q^{M}(\mathbf{x},\mathbf{y})$ defined in (3.4).

4.1 Asymptotics of $q(\mathbf{x},\mathbf{y})$ and $q^{M}(\mathbf{x},\mathbf{y})$

From the transition probabilities $p(t,y|x)$ in (1.9) and $p^{M}(t,y|x)$ in (3.7), we have

[TABLE]

To obtain the probability of the maximal height $\displaystyle\mathbb{P}(\max_{0<t<1}b_{N}(t)<M)$ in (3.3), we need the asymptotics of the above two functions as $\mathbf{x}\to\mathbf{a}$ and $\mathbf{y}\to\mathbf{0}$.

Proposition 1**.**

For $a>0,N\geq 1$, $M\geq a$, as $\mathbf{x}\to\mathbf{a}$ and $\mathbf{y}\to\mathbf{0}$, we have

[TABLE]

and

[TABLE]

where $\Delta(\mathbf{x})$ is defined in (1.3) and

[TABLE]

When $a=0$, we have

[TABLE]

and

[TABLE]

Proof.

By the multilinearity of the determinant, we have from (4.1):

[TABLE]

With the definition of $I_{\alpha}(z)$ in (1.11), we get

[TABLE]

Moreover, from Lemma 1, we obtain the following expression

[TABLE]

Since $\displaystyle\det_{1\leq i,j\leq N}[x_{i}^{2k_{j}}]=0$ and $\displaystyle\det_{1\leq i,j\leq N}[y_{i}^{2k_{j}}]=0$ for any $k_{i}=k_{j}$, we may order the index $\mathbf{k}=(k_{1},k_{2},\cdots,k_{N})$ such that $k_{1}>k_{2}>\cdots>k_{N}$. Moreover, set $\mu_{j}=k_{j}-N+j$. Obviously, we have $\sum\limits_{j=1}^{N}k_{j}=\sum\limits_{j=1}^{N}\mu_{j}+N(N-1)/2$. Replacing the indexes $k_{j}$ by $\mu_{j}$ in the above formula, we obtain from the above formula and the definition of the Schur function in (3.14):

[TABLE]

Recall that $s_{\boldsymbol{\mu}}(\mathbf{0})=0$ unless $\boldsymbol{\mu}=\mathbf{0}$. As $\mathbf{y}\to\mathbf{0}$, the leading term in the above summation comes from the index $\boldsymbol{\mu}=\mathbf{0}$. Moreover, since $s_{\boldsymbol{\mu}}(\mathbf{x}^{2})$ and $s_{\boldsymbol{\mu}}(\mathbf{y}^{2})$ always choose the same index $\boldsymbol{\mu}$, we have

[TABLE]

as $\mathbf{y}\to\mathbf{0}$ uniformly for all $\mathbf{x}$. Since

[TABLE]

we obtain (4.3) for $a>0$ and (4.6) for $a=0$, respectively, by using (4.8) and the above three formulae.

Next, we consider $q^{M}(\mathbf{x},\mathbf{y})$ in (4.2). Similar to the above computations, we obtain from the multilinearity of determinants and Lemma 1,

[TABLE]

Regarding the determinant $\displaystyle\det_{1\leq i,j\leq N}[J_{\alpha}(\frac{x_{n_{j},\alpha}}{M}x_{i})]$, we first expand $J_{\alpha}(\frac{x_{n_{j},\alpha}}{M}x_{i})$ into a Taylor series at $x_{i}=a$. Then, from Lemma 2, we have

[TABLE]

as $\mathbf{x}\to\mathbf{a}$. If the starting point $a=0$, the expansion near $x=a$ is different as $J_{\alpha}(x)$ has an algebraic singularity at the origin; cf. (3.8). In this case, we adopt similar computations in deriving (4.9) to obtain

[TABLE]

Similar to (4.10), the above formula yields

[TABLE]

as $\mathbf{x}\to\mathbf{0}$. Of course, similar formula also holds for $\det\limits_{1\leq i,j\leq N}\left[J_{\alpha}(\frac{x_{n_{j},\alpha}}{M}y_{i})\right]$ when $\mathbf{y}\to\mathbf{0}$.

Finally, with (4.12), (4.13) and (4.15), we obtain the desired asymptotics in (4.4) for $a>0$ and (4.7) for $a=0$, respectively. ∎

4.2 Proof of Theorem 1

The asymptotics for $q(\mathbf{x},\mathbf{y})$ and $q^{M}(\mathbf{x},\mathbf{y})$ in the previous are enough for us to derive the probability for the maximum height.

Recall the expression of $\displaystyle\mathbb{P}(\max_{0<t<1}b_{N}(t)<M)$ in terms of $q(\mathbf{x},\mathbf{y})$ and $q^{M}(\mathbf{x},\mathbf{y})$ in (3.3). When $a>0$, using the asymptotics obtained in (4.3) and (4.4), we have

[TABLE]

where $c_{N}(\alpha)$ is given in (2.3). Now, let us focus on the summation in the above formula and put it into a determinantal form. From the definition of Vandermonde determinants, we have

[TABLE]

Using (3.13), we get

[TABLE]

Moving the preceding factors into the determinant and using the multilinearity of determinant, we obtain (2.2) from the above formula.

When $a=0$, from (4.6) and (4.7), the distribution formula can be written as

[TABLE]

where $\widetilde{c}_{N}(\alpha)$ is given in (2.5). With the aid of (3.13) again, we obtain (2.4) from above formula by similar computations.

This completes the proof of Theorem 1. ∎

4.3 Proof of Theorem 2

Before proving Theorem 2, we need one more property for the derivatives of the Bessel function $J_{\alpha}(z)$.

Lemma 3**.**

For $j\in\mathbb{N}$, we have the following differential relations

[TABLE]

where $P_{j}(z)$, $Q_{j}(z)$, $\widetilde{P}_{j}(z)$ and $\widetilde{Q}_{j}(z)$ are polynomials of degree $j$. Moreover, $P_{2j-1}(z)$ and $\widetilde{P}_{2j-1}(z)$ are odd functions, while $Q_{2j}(z)$ and $\widetilde{Q}_{2j}(z)$ are even functions.

Proof.

We prove this lemma by induction. It is well-known that the Bessel functions satisfy the following recurrence relations

[TABLE]

see [20, Eq. (10.6.2)]. From the above formulae, it is easy to see

[TABLE]

The above formulae show that the lemma holds for $j=1$.

Assume the lemma holds for $j=k$. Then, we have

[TABLE]

Taking derivative on both sides of the above formula and using (4.20) again, we obtain

[TABLE]

with

[TABLE]

So, (4.18) holds for $j=k+1$. Differentiating (4.22) one more time, we also get (4.19) for $j=k+1$ through similar computations.

This finish the proof of the lemma. ∎

With the above preparation, we prove our second theorem below.

Proof of Theorem 2. Let us focus on the determinant in (2.2). The entries in the first row are given by

[TABLE]

For the second row, by using (4.20), we replace $J_{\alpha}^{\prime}(\frac{a}{M}x_{n,\alpha})$ with $J_{\alpha+1}(\frac{a}{M}x_{n,\alpha})$ and $J_{\alpha}(\frac{a}{M}x_{n,\alpha})$. This gives us

[TABLE]

Note that the coefficient of $J_{\alpha}$ term is independent of $n$. Applying a row operation, we eliminate the $J_{\alpha}$ term in the second row and obtain

[TABLE]

Regarding the third row, we use (4.19) to replace $x_{n,\alpha}^{2}J_{\alpha}^{\prime\prime}(\frac{a}{M}x_{n,\alpha})$ and obtain

[TABLE]

Note that $d_{0}$ and $d_{1}$ in the above formulae are two $n$-independent constants. Then, with the first two rows given in (4.23) and (4.24), we again apply row operations to simplify the entries in the third row as

[TABLE]

One may repeat the above procedures by first replacing the $(-1)^{i-1}x_{n,\alpha}^{i-1}J_{\alpha}^{(i-1)}(\frac{a}{M}x_{n,\alpha})$ term with $J_{\alpha+1}$ and $J_{\alpha}$ functions, then apply row operations to simplify the entries. Here, we also need the properties that the polynomials $P_{2k-1}$, $\widetilde{P}_{2k-1}$ and $\widetilde{Q}_{2k}$, $Q_{2k}$ in (4.18) and (4.19) are even or odd functions. Therefore, the determinant in (2.2) can be put into a form as follows:

[TABLE]

where

[TABLE]

Note that the entries in the right-hand side of (4.26) are indeed the moments $m_{k}^{(i)},\ i=1,2,$ defined in (2.10). Finally, we obtain (2.16) by interchanging rows in (4.26), such that all the ones containing $J_{\alpha}(\frac{a}{M}x_{n,\alpha})$ are moved to the top block, while the ones with $J_{\alpha+1}(\frac{a}{M}x_{n,\alpha})$ are moved to the bottom block.

For the case $a=0$, it is straightforward to see that the entries in the determinant of (2.4) are entries associated with the discrete orthogonal polynomials defined in (2.15). Then, the formula (2.17) follows immediately.

This completes the proof of Theorem 2. ∎

5 Numerical simulations and conclusions

5.1 Numerical simulations

Based on the explicit formulae in Theorem 1 and 2, we compute $\displaystyle\mathbb{P}(\max_{0<t<1}b_{N}(t)<M)$ numerically. We set the number of Bessel paths $N$ to be $10$ and the order of Bessel functions $\alpha$ to be 1 as an illustration.

As expected, when the starting point $a$ is fixed, the probability $\displaystyle\mathbb{P}(\max_{0<t<1}b_{N}(t)<M)$ is positive and increase with respect to $M$. When $M$ is large, the probability tends to 1; see Figure 2 and Table 5.1. When the upper constraint $M$ is fixed, the probability decreases with respect to $a$, and tends to 0 when $a$ is large; see Figure 2 and Table 5.1.

[TABLE]

\captionof

tableThe probability with $N=10,\alpha=1$ and $a=1$.

[TABLE]

\captionof

tableThe probability with $N=10,\alpha=1$ and $M=5$.

5.2 Conclusions

In summary, we derive the probability $\displaystyle\mathbb{P}(\max_{0<t<1}b_{N}(t)<M)$ for the maximum height of $N$ non-intersecting Bessel paths exactly. The probability can be put in terms of Hankel determinants associated with the multiple discrete orthogonal polynomials or discrete orthogonal polynomials, respectively, depending on the starting point $a>0$ or $a=0$. It is well-known that the Deift-Zhou nonlinear steepest descent method for Riemann-Hilbert problems is a powerful tool to study asymptotic problems related to (multiple) orthogonal polynomials. Therefore, this relation is important for us to derive the asymptotics of $\displaystyle\mathbb{P}(\max_{0<t<1}b_{N}(t)<M)$ as the number of paths $N\to\infty$. With suitable scaling, we expect the limiting distribution converges to the Tracy-Widom distribution for the largest eigenvalue of the Gaussian orthogonal ensemble due to its close relation to the Airy2 process. We will leave the asymptotic study in a forthcoming publication.

Acknowledgments

This work was partially supported by grants from the Research Grants Council of the Hong Kong Special Administrative Region, China (Project No. CityU 11300115, CityU 11303016), and by grants from City University of Hong Kong (Project No. 7004864, 7005032).