Fractional Immigration-Death Processes

TL;DR

This paper investigates explicit solutions and stochastic representations for fractional immigration-death processes, analyzing their properties and limit distributions using spectral and spectral methods.

Contribution

It introduces a novel spectral approach to solve fractional difference-differential equations related to immigration-death processes and provides a stochastic representation for these solutions.

Findings

Explicit strong solutions derived using spectral methods

Stochastic representation via stable time-changed processes

Limit distribution analysis of the time-changed process

Abstract

In this paper we study explicit strong solutions for two difference-differential fractional equations, defined via the generator of an immigration-death process, by using spectral methods. Moreover, we give a stochastic representation of the solutions of such difference-differential equations by means of a stable time-changed immigration-death process and we use this stochastic representation to show boundedness and then uniqueness of these strong solutions. Finally, we study the limit distribution of the time-changed process.

Equations

$$C_n(x,\alpha )=C_x(n,\alpha ),n,x\in {\mathbb{N}}_0,$$

$$p(t,x;y)=P(N(t+s)=x|N(s)=y),x,y=0,1,2,\dots ;t,s\ge 0,$$

$$P^{\prime }(t)=AP(t),P^{\prime }(t)=P(t)A,$$

$$A(x,x+1)=B(x)x\ge 0,$$

$$A(x,x-1)=D(x)x\ge 1,$$

$$P_{i,j}(t)\ge 0,\sum _{j=0}^{+\infty }{P}_{i,j}(t)\le 1.$$

$$\begin{array}{rl}Gf(x)& =(B(x)-D(x)){\nabla }^{+}f(x)+D(x)\Delta f(x)\\ & =(B(x)-D(x)){\nabla }^{-}f(x)+B(x)\Delta f(x),\end{array}$$

$${\nabla }^{+}f(x)$$

$${\nabla }^{-}f(x)$$

$$\Delta f(x)$$

$${\nabla }^{+}(fg)(x)$$

$${\nabla }^{-}(fg)(x)$$

$$\Delta (fg)(x)$$

$$\{\begin{array}{l}{p}^{\prime }(t,x;y)=Gp(t,x;y)\\ p(0,x;y)={\delta }_{x,y},\end{array}$$

$${\delta }_{x,y}=\{\begin{array}{ll}1& x=y\\ 0& \text{mbox}otherwise.\end{array}$$

$$\begin{array}{rl}Lf(x)& =-{\nabla }^{-}((B(\cdot )-D(\cdot ))f)(x)+\Delta (D(\cdot )f)(x)\\ & =-{\nabla }^{+}((B(\cdot )-D(\cdot ))f)(x)+\Delta (B(\cdot )f)(x),\end{array}$$

$$\{\begin{array}{l}{p}^{\prime }(t,x;y)=Lp(t,x;y)\\ p(0,x;y)={\delta }_{x,y},\end{array}$$

$$G=p_1(x){\nabla }^{+}+p_2(x)\Delta ,$$

$$Gf(x)=-\lambda f(x),$$

$${\nabla }^{+}(p_2(x)m(x))=p_1(x)m(x).$$

$$G=(a-bx){\nabla }^{-}+a\Delta ;$$

$$\frac{du}{dt}(t,x)=Gu(t,x).$$

$$B(x)=a,$$

$$Lf(x)=-{\nabla }^{+}((a-bz)f(z))(x)+a\Delta f(x),$$

$$\begin{array}{c}G(x,x-1)=bxG(x,x)=-(a+bx)G(x,x+1)=a\\ G(0,0)=-aG(0,1)=a\end{array}$$

$$\begin{array}{c}L(x,x-1)=aL(x,x)=-(a+bx)L(x,x+1)=b(x+1)\\ L(0,0)=-aL(0,1)=b.\end{array}$$

$$m(\{x\})=e^{-\alpha }\frac{{\alpha }^x}{x!},x=0,1,2,\dots .$$

$${\Vert f\Vert }_{{\ell }^{\infty }}=\underset{x\in {\mathbb{N}}_0}{\max }|f(x)|;$$

$${\Vert f\Vert }_{{\ell }^1}=\sum _{x=0}^{+\infty }|f(x)|<+\infty .$$

$${\Vert f\Vert }_{{\ell }^2}^2:=\sum _{x=0}^{+\infty }{f}^2(x)<+\infty$$

$$⟨f,g{⟩}_{{\ell }^2}=\sum _{x=0}^{+\infty }f(x)g(x)$$

$${\Vert f\Vert }_{{\ell }^2(m)}^2:=\sum _{x=0}^{+\infty }m(\{x\})f^2(x)<+\infty$$

$$⟨f,g{⟩}_{{\ell }^2(m)}=\sum _{x=0}^{+\infty }m(\{x\})f(x)g(x).$$

Full text

N. Leonenko was supported in particular by Australian Research Council’s Discovery Projects funding scheme (project DP160101366), and by project MTM2015-71839-P of MINECO, Spain (co-funded with FEDER funds).

1. Introduction

Birth-death processes constitute an important class of continuous time Markov chain (CTMC). They are widely used, for instance, in population and evolutionary dynamics (see [34, 35]), queueing theory (see [42]) and in epidemiology (see [3]). A complete classification and characterization of birth-death processes is due to Karlin and McGregor, whose papers [16, 17] are the starting point of the study of family of classical orthogonal polynomials linked to such processes.

Classical orthogonal polynomials are widely used to study the solutions of Kolmogorov equations as in the case in which the state space of the process is continuous, as well as in the discrete one. In the continuous case, the families of classical orthogonal polynomials are used to give a spectral decomposition of Kolmogorov equations induced by the generators of Pearson diffusions [12]. In the discrete case, the discrete analogue of Pearson diffusions is given by a certain class of solvable birth-death processes [21]. Moreover one can associate to any family of classical orthogonal polynomials of discrete variable another particular family, called the dual family [33]. In some cases, a family of classical orthogonal polynomials of discrete variable could be in duality with itself: in this case it is called self-dual family [41]. Among self-dual families, the simplest one is the family of Charlier polynomials, whose self-duality is induced by the following formula

[TABLE]

called duality formula for Charlier polynomials (see Section 3 for the definition of Charlier polynomials).

Charler polynomials are really useful in the study of immigration-death processes (or $M/M/\infty$ queues) [41] and in their general version on $1$-dimensional lattice, called Charlier processes [1, 21]. Indeed, one can give a spectral decomposition of the strong solutions of Kolmogorov equations induced by the generator of the immigration-death processes in terms of such polynomials.

For Pearson diffusions, the classical orthogonal polynomials are powerful tools to study strong solutions of fractional Kolmogorov equations and characterize a stochastic representation of such solutions via time-changed (through the inverse of a Lèvy subordinator) Markov processes [13, 23, 24, 25]. In the discrete case, fractional (time-changed) processes have been widely considered via different approaches. First of all, a fractional version of the Poisson process has been introduced using Mittag-Leffler distributed inter-jump times instead of exponential ones [2, 22, 26, 28, 29] (this approach has been also applied to general counting processes [10]). Such process can be also obtained using a fractional differential-difference equations approach [6, 7] and by means of a time-change [30].

With the same approach, some classes of fractional birth-death processes have been introduced and studied [36, 37, 38]: in these papers, properties of these processes are deduced from a fractional version of their Kolmogorov forward equation.

Here, following the approach of [23], we show the existence of strong solutions for the time-fractional counterpart of the Kolmogorov backward and forward equations of immigration-death processes with the aid of Charlier polynomials and link them to a time-changed immigration-death process.

In particular:

• •

in Section 2 we give some basics on birth-death processes;

• •

in Section 3 we give some notions on the classical immigration-death process, defining its generator and its forward operator;

• •

in Section 4 we show the existence of strong solutions of the time-fractional Kolmogorov backward and forward equations under suitable assumption on the initial data;

• •

in Section 5 we introduce a fractional immigration-death process and show how the strong solutions of the time-fractional Kolmogorov backward and forward equations can be interpreted by using such process;

• •

in Section 6, we show the uniqueness of such strong solutions by using the aforementioned stochastic representation and a uniqueness criterion for uniformly bounded solutions [4], under suitable assumptions on the initial data;

• •

finally, in Section 7 we give the limit distribution of the constructed fractional immigration-death process and we discuss its autocovariance function.

2. Birth-death processes

Let us give some information about general birth-death processes, following the lines of [16, 17]. We say that a time-homogeneous continuous time Markov chain $N(t)$ defined on $\mathbb{N}_{0}=\{0,1,2,\dots\}$ is a birth and death process if and only if, denoting with

[TABLE]

the transition probability functions and $P(t)=(p(t,x;y))_{x,y\geq 0}$ the transition probability matrix, it is solution of the following two differential equations

[TABLE]

with initial condition $P(0)=I$ and the infinite matrix $\operatorname{\mathcal{A}}=(A(x,y))_{x,y\geq 0}$ is such that:

[TABLE]

where $B(x)>0$ for any $x\geq 0$, $D(x)>0$ for any $x\geq 1$ and $D(0)\geq 0$. Equations (2.1) are called respectively backward and forward Kolmogorov equation. In order to obtain $P(t)$ we need to impose other two properties:

[TABLE]

In particular it is possible to show that $N(t)$ is a birth-death process if and only if its generator is given by:

[TABLE]

for $x=0,1,2,\dots$ and $f(-1)=0$, where the difference-type operators $\nabla^{\pm}$ and $\Delta$ are defined as

[TABLE]

The following discrete versions of the Leibnitz rule will be useful

[TABLE]

The backward Kolmogorov equation becomes, for fixed $x\in\mathbb{N}_{0}$

[TABLE]

where $\operatorname{\mathcal{G}}$ works on $y$ and

[TABLE]

is Kronecker symbol.

Moreover we can find a forward operator

[TABLE]

so that for fixed $y\in\mathbb{N}_{0}$ the forward Kolmogorov equation becomes

[TABLE]

where $\operatorname{\mathcal{L}}$ works on $x$.

We will focus on the case in which the generator is in the form:

[TABLE]

where $p_{1}(x)$ and $p_{2}(x)$ are polynomials such that $\deg p_{1}(x)\leq 1$ and $\deg p_{2}(x)\leq 2$. Then we can find the classical orthogonal polynomials of discrete variable as solution of the equation

[TABLE]

for some $\lambda$, which is an hypergeometric type difference equation. The values that these polynomials assume on a lattice $\{D_{1},D_{1}+1,\dots,D_{2}\}$ for some $D_{1},D_{2}$ fully characterize the transition probability and the solutions of the backward and forward Kolmogorov equations. Moreover, these polynomials respect an orthogonality relation in $\ell^{2}(\operatorname{\mathbf{m}})$ for some measure $\operatorname{\mathbf{m}}$ called the spectral measure, which is an atomic measure on the lattice. In this case, the spectral measure coincides with the invariant measure of the process $N(t)$ and its mass function $m(x)=\operatorname{\mathbf{m}}(\{x\})$ is solution of a discrete analogue of the Pearson equation

[TABLE]

Following the lines of [21], for $p_{1}(x)=a-bx$, we can recognize the following three class of solvable birth-death processes:

• •

For $p_{2}(x)=bx$ we have the Immigration-Death process;

• •

For $p_{2}(x)=\frac{1}{2}\sigma^{2}x$ where $\frac{1}{2}\sigma^{2}\not=b$ we have a negative binomial process;

• •

For $p_{2}(x)=\frac{1}{2}\sigma^{2}x(A-x)$ we have a hypergeometric process.

However, we will focus only on the first case for the choice of the polynomials $p_{1}$ and $p_{2}$.

3. Immigration-death processes

Fix $a,b>0$ the operator

[TABLE]

which is a discrete version of the Ornstein-Uhlenbeck generator on $\mathbb{N}_{0}$.

A continuous time Markov chain $N(t)$ defined on $\mathbb{N}_{0}$ that admits $\operatorname{\mathcal{G}}$ as generator will be called immigration-death process (or also $M/M/\infty$ queue: see, for instance, [41]). This process can be generalized to a particular birth-death process with values on a $1$-dimensional lattice called Charlier process (see [1]), but we will focus on the $\mathbb{N}_{0}$-valued one. For such process, the backward Kolmogorov equations are in the form

[TABLE]

Moreover, from $\operatorname{\mathcal{G}}$ we can recognize the birth and death parameters as

[TABLE]

and thus the forward operator as

[TABLE]

where with $\nabla^{+}((a-bz)f(z))(x)$ we intend the operator $\nabla^{+}$ applied to the function $z\mapsto(a-bz)f(z)$ and then evaluated in $x$.

The operators $\operatorname{\mathcal{G}}$ and $\operatorname{\mathcal{L}}$ can be represented as infinite matrices. In particular we have $\operatorname{\mathcal{G}}=(G(x,y))_{x,y\geq 0}$ where, for $x>0$

[TABLE]

and $\operatorname{\mathcal{L}}=(L(x,y))_{x,y\geq 0}$ where, for $x>0$

[TABLE]

The stationary measure of the process $N(t)$ is the Poisson distribution of parameter $\alpha$, given by:

[TABLE]

Now let us introduce the main Banach sequence spaces we will use through this paper:

• •

Let us denote with $\ell^{\infty}$ the Banach space of bounded functions $f:\mathbb{N}_{0}\to\operatorname{\mathbb{R}}$ equipped with the norm

[TABLE]

• •

Let us denote with $c_{0}$ the subspace of $\ell^{\infty}$ of bounded functions $f:\mathbb{N}_{0}\to\operatorname{\mathbb{R}}$ such that $\lim_{x}f(x)=0$;

• •

Let us denote with $\ell^{1}$ the Banach space of the functions $f:\mathbb{N}_{0}\to\operatorname{\mathbb{R}}$ such that

[TABLE]

• •

Let us denote with $\ell^{2}$ the Hilbert space of functions $f:\mathbb{N}_{0}\to\operatorname{\mathbb{R}}$ such that

[TABLE]

equipped with the scalar product

[TABLE]

• •

Let us denote with $\ell^{2}(\operatorname{\mathbf{m}})$ the Hilbert space of functions $f:\mathbb{N}_{0}\to\operatorname{\mathbb{R}}$ such that

[TABLE]

equipped with the scalar product

[TABLE]

Remark 3.1*.*

Let us observe that $\ell^{2}$ is continuously included in $\ell^{2}(\operatorname{\mathbf{m}})$. Consider a function $f\in\ell^{2}$. Then

[TABLE]

Now, let us observe that the sequence $x\mapsto\frac{\alpha^{x}}{x!}$ converges to [math] as $x\to+\infty$, hence there exists a constant $C(\alpha)$ such that $\frac{\alpha^{x}}{x!}\leq C(\alpha)$. Thus

[TABLE]

Moreover, since $\ell^{1}$ is continuously included in $\ell^{2}$ (see [44]), we have that $\ell^{1}$ is also continuously included in $\ell^{2}(\operatorname{\mathbf{m}})$.

From the matrix representation of the generator $\operatorname{\mathcal{G}}$ and the forward operator $\operatorname{\mathcal{L}}$ one can prove the following Lemma.

Lemma 3.2**.**

The operators $\operatorname{\mathcal{G}}:\ell^{2}(\operatorname{\mathbf{m}})\mapsto\ell^{2}(\operatorname{\mathbf{m}})$ and $\operatorname{\mathcal{L}}:\ell^{2}(\operatorname{\mathbf{m}})\mapsto\ell^{2}(\operatorname{\mathbf{m}})$ are continuous.

Proof.

The proof is a straightforward consequence of Schur’s test (see [15]) ∎

Moreover, another interesting property that follows from the matrix representation of $\operatorname{\mathcal{G}}$ is given by the following Lemma.

Lemma 3.3**.**

The process $N(t)$ is a Feller process.

Proof.

The proof is a straightforward consequence of [11, Corollary $3.2$]. ∎

Let us also observe that the spectrum of $\operatorname{\mathcal{G}}$ is given by the sequence $\lambda_{n}=-bn$ ([21]), while the eigenfunctions are defined as $x\mapsto C_{n}(x,\alpha)$ where $\alpha=\frac{a}{b}$ and $C_{n}$ are the Charlier polynomials (see [33, 41]), which are defined by the generating function

[TABLE]

or via the three terms recurrence relations:

[TABLE]

where $C_{0}(x,\alpha)\equiv 1$ and $C_{-1}(x,\alpha)\equiv 0$, or

[TABLE]

The first few Charlier polynomials are

[TABLE]

The orthogonality relation between the polynomials $C_{n}$ is given by

[TABLE]

where $\delta_{n,m}$ is the Kronecker delta symbol. Thus, posing $d_{n}^{2}=\frac{n!}{\alpha^{n}}$, we have that

[TABLE]

Let us then define an orthonormal system of polynomials given by

[TABLE]

Let us also recall that we can exploit the decomposition of a function $g\in\ell^{2}(\operatorname{\mathbf{m}})$ by means of the orthonormal basis $\{Q_{n}\}_{n\in\mathbb{N}_{0}}$. Indeed for any $g\in\ell^{2}(\operatorname{\mathbf{m}})$, given the decomposition $g(x)=\sum_{n=0}^{+\infty}g_{n}Q_{n}(x)$ where $g_{n}=\langle g,Q_{n}\rangle_{\ell^{2}(\operatorname{\mathbf{m}})}$, the sequence $\{g_{n}\}_{n\in\mathbb{N}_{0}}\in\ell^{2}$.

By using such orthonormal system of polynomials, it is well known (see [16, 17] but also [21] for a review) that the transition probability function of the immigration-death process is given by

[TABLE]

where $m(x)=\operatorname{\mathbf{m}}(\{x\})$ and is the fundamental solution of the backward Kolmogorov equation, that is to say that the Cauchy problems

[TABLE]

and

[TABLE]

with $g,f/m\in\ell^{2}(\operatorname{\mathbf{m}})$ admit strong solutions $v$ given by

[TABLE]

and

[TABLE]

where $g(x)=\sum_{n=0}^{+\infty}g_{n}Q_{n}(x)$ and $f(x)/m(x)=\sum_{n=0}^{+\infty}f_{n}Q_{n}(x)$ and the convergence is uniform. In particular from (3.2) one easily obtains that

[TABLE]

4. Strong solutions in the fractional case

Let us introduce the fractional derivative operator (see [27]). Fix $\nu\in(0,1)$ and consider the Caputo fractional derivative given by

[TABLE]

that, if $u$ is differentiable in $t$, can be written also as

[TABLE]

and pose for $\nu=1$ $\frac{\partial^{\nu}u}{\partial t^{\nu}}=\frac{\partial u}{\partial t}$. Note that the classes of functions for which the Caputo fractional derivative is well defined are discussed in [31, Section $2.2$ and $2.3$] (in particular one can use the class of absolutely continuous functions).

Denote with

[TABLE]

the one-sided Laplace transform of $u$ with respect to $t$. Thus we have that the Laplace transform of $\frac{\partial^{\nu}u}{\partial t^{\nu}}$ is given by

[TABLE]

We want to find strong solutions for fractional Cauchy problems in the form:

[TABLE]

for $g\in\ell^{2}(\operatorname{\mathbf{m}})$ with the decomposition $g(x)=\sum_{n=0}^{+\infty}g_{n}Q_{n}(x)$. The main idea is to find a solution via separation of variables. Indeed, we can suppose that $u(t,x)=T(t)\varphi(x)$ and then observing that, if $u$ is solution of the first equation of 4.2, then

[TABLE]

that leads, if $\varphi$ and $T$ do not vanish, to the two coupled equations:

[TABLE]

which are two eigenvalue problems. In particular we have observed that the first one admits a non zero solution if and only if $\lambda=-bn$ for some $n\in\mathbb{N}_{0}$ and in that case we can consider $\varphi(x)=Q_{n}(x)$. Moreover, the second problem admits a solution in the form

[TABLE]

where $E_{\nu}$ is the Mittag-Leffler function defined as

[TABLE]

(see, for instance, [19]). Thus the idea is to find a solution in the form

[TABLE]

Moreover, the initial condition suggests that

[TABLE]

so we have $u_{n}=g_{n}$ and then we expect the solution to be

[TABLE]

These heuristic arguments have shown us how should the solution look like, hence we have to prove that such function $u$ is the solution we are searching for.

With the following Lemma, we will first exhibit the fundamental solution of the fractional Cauchy problem in Eq. (4.2).

Lemma 4.1**.**

Consider the series

[TABLE]

where $Q_{n}$ and $E_{\nu}$ are the functions defined in Equations (3.1) and (4.3) and $m(x)=\operatorname{\mathbf{m}}(\{x\})$. Then such series converges for fixed $t>0$ and $x,y\in\mathbb{N}_{0}$.

Proof.

To show the convergence of $p_{\nu}(x,t;y)$, we need the following self-duality property of the Charlier polynomials (see [33, Equation $2.7.10a$]):

[TABLE]

From this relation we have

[TABLE]

hence we need to show the convergence of the series

[TABLE]

Now, let us observe that equation (4.6) made us fix the degrees of the polynomials involved in the series. Thus, let us denote with $z_{x}$ and $z_{y}$ the last real zeroes of $C_{x}(\cdot,\alpha)$ and $C_{y}(\cdot,\alpha)$ and then let us consider $n_{0}>\max\{z_{x},z_{y}\}$.

We will equivalently prove that the series

[TABLE]

converges. To do this, we need to recall another property of the Charlier polynomials. In particular it is known (see [33, Table $2.3$]) that the director coefficient of $C_{n}(\cdot,\alpha)$ is given by

[TABLE]

In particular, recalling that $\alpha=\frac{a}{b}$, $\alpha>0$ since $a,b>0$ and then $c_{n}>0$ if $n$ is even and $c_{n}<0$ if $n$ is odd.

By using this observation, we can distinguish two cases:

• i

If $x+y$ is even, then, since $c_{x}c_{y}>0$, for any $n\geq n_{0}$ $C_{x}(n,\alpha)C_{y}(n,\alpha)>0$ and then the series (4.7) admits only positive summands. Recalling that $E_{\nu}(-bnt^{\nu})\leq 1$ we obtain

[TABLE]

where the RHS series converges since

[TABLE]

• ii

If $x+y$ is odd, then, since $c_{x}c_{y}<0$, for any $n\geq n_{0}$ $C_{x}(n,\alpha)C_{y}(n,\alpha)<0$ and then the series (4.7) admits only negative summands. As before, we obtain

[TABLE]

where the RHS series converges for equation (4.8).

∎

With Lemma 4.1, we have exploited the fundamental solution of the equation in (4.2). Now we have to show that a function in the form (4.4) is a solution for such fractional Cauchy problem. To do this, let us first show a technical lemma.

Lemma 4.2**.**

For any $t_{0}>0$, there exists a constant $K(t_{0},\nu)$ such that

[TABLE]

Proof.

Let us use the uniform estimate for the Mittag-Leffler function given in [43, Theorem $4$]:

[TABLE]

Consider the function

[TABLE]

Thus we have

[TABLE]

hence the function $f$ is strictly increasing. So we have

[TABLE]

and then

[TABLE]

∎

Now let us exhibit a strong solution for our fractional Cauchy problem.

Theorem 4.3**.**

Let $g\in\ell^{2}(\operatorname{\mathbf{m}})$ with decomposition $g(x)=\sum_{n=0}^{+\infty}g_{n}Q_{n}(x)$. Then the fractional difference-differential Cauchy problem

[TABLE]

admits a strong solution $u$ in the form

[TABLE]

Proof.

First let us observe that obviously if $u$ is in the form (4.10), then $u(0,x)=g(x)$.

Now, let us notice that

[TABLE]

Hence we need to show that the series in (4.10) is convergent at least uniformly in $t$ and that we can change the series with the operators.

Starting from the convergence of the series, by using Cauchy-Schwartz inequality we have

[TABLE]

hence the series in (4.10) totally converges.

Now we need to show that one can exchange the operators with the series. To do that, let us first observe that

[TABLE]

and since $(t-\tau)^{1-\nu}$ is strictly decreasing in $[0,t]$ we can use [39, Theorem $7.16$] with the total convergence of the series (4.10) to obtain

[TABLE]

Now we want to use the following relation:

[TABLE]

but to do this, by using [39, Theorem $7.17$], we need to show the uniform convergence of

[TABLE]

in any compact interval included in $(0,+\infty)$. Hence, by definition of Caputo fractional derivative, as given in (4.1), we really need to show the uniform convergence of

[TABLE]

in any interval of the form $[t_{0},+\infty)$. To do this, let us recall that $\frac{d^{\nu}}{dt^{\nu}}E_{\nu}(-bnt^{\nu})=-bnE_{\nu}(-bnt^{\nu})$ and thus we need to show the uniform convergence of

[TABLE]

Thus, fix $t_{0}>0$ and observe that

[TABLE]

where the first inequality follows from Lemma 4.2 and the second inequality from Cauchy-Schwartz inequality as done before in (4.11). Hence we have shown the total convergence of (4.12) in any interval of the form $[t_{0},+\infty)$.

We have already shown that $\sum_{n=0}^{+\infty}E_{\nu}(-bnt^{\nu})Q_{n}(x)g_{n}$ totally converges with respect to $t$: in the same way we have that also $\sum_{n=0}^{+\infty}E_{\nu}(-bnt^{\nu})Q_{n}(x-1)g_{n}$ and $\sum_{n=0}^{+\infty}E_{\nu}(-bnt^{\nu})Q_{n}(x+1)g_{n}$ totally converge with respect to $t$.

Now, observe that

[TABLE]

and in the same way one can show that

[TABLE]

By using these last two relations, it is easy to show that

[TABLE]

Finally we have that

[TABLE]

and we have concluded the proof. ∎

The same strategy can be used to exhibit a strong solution to the fractional forward Kolmogorov equation.

Theorem 4.4**.**

Let $f$ be a function such that $f/m\in\ell^{2}(\operatorname{\mathbf{m}})$ with decomposition $f(x)/m(x)=\sum_{n=0}^{+\infty}f_{n}Q_{n}(x)$. Then the fractional difference-differential Cauchy problem

[TABLE]

admits a strong solution $u=u(t,x)$ given by

[TABLE]

Proof.

Since $\{f_{n}\}_{n\in\mathbb{N}}\in\ell^{2}$, then, from the previous theorem, we already know that we can exchange operators and series. We only need to prove that the single summand of the series is a solution of the equation and that $u(0,x)=f(x)$.

Let us first notice that

[TABLE]

thus the function $u$ satisfies the given initial condition.

To show that the single summand is solution of the equation, let us write $\operatorname{\mathcal{L}}$ as

[TABLE]

for a generic function $h$.

Moreover, let us observe that

[TABLE]

Let us also recall that $m$ solves a discrete Pearson equation:

[TABLE]

Now, let us observe that

[TABLE]

hence we will only study $\operatorname{\mathcal{L}}(m(\cdot)Q_{n}(\cdot))$. In particular we have

[TABLE]

hence, by using the Discrete Leibnitz Rule ((2.3) and (2.4)), we obtain

[TABLE]

First let us observe that $\Delta=\nabla^{-}\nabla^{+}$, then

[TABLE]

since $m$ satisfies equation (4.15). Moreover

[TABLE]

while

[TABLE]

thus

[TABLE]

Finally, we obtain:

[TABLE]

∎

Remark 4.5*.*

It is easy to see that $p_{\nu}(t,x;y)$ is strong solution of the fractional backward equation

[TABLE]

where $\operatorname{\mathcal{G}}$ operates on $y$, and is also strong solution of the fractional forward equation

[TABLE]

where $\operatorname{\mathcal{L}}$ operates on $x$. In particular, as shown by Theorems 4.3 and 4.4, it is the fundamental solution of such equations.

5. Stochastic representation of the solutions

Now we want to exhibit a process whose “transition probability” is the fundamental solution $p_{\nu}(t,x;y)$ we have described previously.

To do this, let us consider a classical immigration-death process $N_{1}(t)$ (as defined before). Let us also consider a $\nu$-stable subordinator $\sigma_{\nu}(t)$ with Laplace transform

[TABLE]

and its inverse process (or first passage time process) $L_{\nu}(t)$ defined as

[TABLE]

The latter admits density (see [2, 32])

[TABLE]

where $g_{\nu}(x)$ is the density of $\sigma_{\nu}(1)$ given by

[TABLE]

Alternatives for $f_{\nu}(y,t)$ are given in [18, 20].

Thus, let us define the fractional immigration-death process as $N_{\nu}(t):=N_{1}(L_{\nu}(t))$. This is a semi-Markov process as defined in [14]. However, we say that such process admits a transition probability mass $p_{\nu}(t,x;y)$ if for any $B\subseteq\mathbb{N}_{0}$:

[TABLE]

Hence, we can use such process to characterize the fundamental solution we found in the previous section.

Theorem 5.1**.**

The process $N_{\nu}(t)$ admits a transition probability mass $p_{\nu}(t,x;y)$ in the form (4.5).

Proof.

Let us first recall that (see, for instance, [32]) the process $L_{\nu}(t)$ admits a density $f_{t}(\tau)=\operatorname{\mathbb{P}}(L_{\nu}(t)\in d\tau)$. Moreover, let us recall (see [8]) that

[TABLE]

Now, observe that for any $B\subseteq\mathbb{N}_{0}$, since $N_{1}(t)$ admits a transition probability mass, we have

[TABLE]

Now, if $B$ is a finite set, we have

[TABLE]

If $B$ is infinite, let us consider the sets $I_{m}:=\{x\in\mathbb{N}_{0}:\ x\leq m\}$ and $B_{m}:=B\cap I_{m}$. Thus, $B_{m}$ is finite and then

[TABLE]

Since $p_{1}(\tau,x;y)f_{t}(\tau)$ is non-negative, we can use the monotone convergence theorem to obtain, taking the limit as $m\to+\infty$

[TABLE]

Now let us only consider

[TABLE]

and recall that (see, for instance, [21] for the specific case of the immigration-death process, but in general such decomposition can be found in [16, 17])

[TABLE]

Hence we have

[TABLE]

Now we have to show that we can exchange integral and series. To do this, let us first observe that

[TABLE]

Let us consider $z_{x}$ and $z_{y}$ the last real zeros of $C_{x}(n,\alpha)$ and $C_{y}(n,\alpha)$ and consider a $n_{0}\in\mathbb{N}$ such that $n_{0}>\max\{z_{x},z_{y}\}$. Thus we have

[TABLE]

and

[TABLE]

Now, fix $\tau_{0}>0$ and observe that for $\tau>\tau_{0}$ and $n\geq n_{0}+1$ the function

[TABLE]

does not change sign, by Fubini’s theorem (see [40, Theorem $8.8$]) we have that

[TABLE]

Now we have to pass to the limit as $\tau_{0}\to 0$. To do this, let us observe that

[TABLE]

and let us distinguish two cases.

• $i)$

If $x+y$ is even,

[TABLE]

and in particular we have

[TABLE]

where

[TABLE]

as we observed before. Then we can use dominated convergence theorem to take the limit as $\tau_{0}\to 0$ and obtain

[TABLE]

• $ii)$

If $x+y$ is odd, then

[TABLE]

and in particular we have

[TABLE]

where

[TABLE]

as we observed before. Then we can use dominated convergence theorem to take the limit as $\tau_{0}\to 0$ and obtain

[TABLE]

Hence in general we have for any $x,y\in\mathbb{N}_{0}$

[TABLE]

and then

[TABLE]

Finally we have

[TABLE]

Thus $p_{\nu}(t,x;y)$ exists and

[TABLE]

∎

Now we are ready to prove the following Theorem.

Theorem 5.2**.**

Let $g\in c_{0}$ such that $g(x)=\sum_{n=0}^{+\infty}g_{n}Q_{n}(x)$. Then the function

[TABLE]

is a solution of (4.9).

Proof.

First of all observe that

[TABLE]

hence $g\in\ell^{2}(\operatorname{\mathbf{m}})$ and we are under the hypotheses of Theorem 4.3.

Consider a generic $f\in c_{0}$ and define the family of operators

[TABLE]

In particular, by Lemma 3.3, we know that $N_{1}(t)$ is a Feller process, hence $(T_{t})_{t\geq 0}$ is a Feller semigroup. Moreover, strong continuity of $(T_{t})_{t\geq 0}$ follows from [9, Lemma $1.4$]. Then, by using [5, Theorem $3.1$], we know that, since $\operatorname{\mathcal{G}}$ is the generator of $T_{t}$, the function

[TABLE]

where $g_{\nu}(s)$ is the density of $\sigma_{\nu}(1)$, is a solution of (4.9). But if we use the change of variables $\tau=\left(\frac{t}{s}\right)^{\nu}$, and the fact that $f_{t}(\tau)=\frac{t}{\nu}\tau^{-1-\frac{1}{\nu}}g_{\nu}(t\tau^{-\frac{1}{\nu}})$ for $\tau\geq 0$ (see, for instance, [32]), we obtain

[TABLE]

∎

Finally, we can provide the stochastic representation of solutions of (4.14).

Corollary 5.3**.**

Let $p_{\nu}(t,x;y)$ be the transition density of $N_{\nu}(t)$. Then, for any $f$ such that $\frac{f}{m}\in\ell^{2}(\operatorname{\mathbf{m}})$ with decomposition $f(x)/m(x)=\sum_{n=0}^{+\infty}f_{n}Q_{n}(x)$. Thus

[TABLE]

is a solution of (4.14).

Proof.

This easily follows from Theorems 5.1 and 4.4. ∎

We can use the last Corollary to exploit the asymptotic behaviour of the density of the process $N_{\nu}(t)$.

6. Uniqueness of strong solutions

In this section, we aim to show that the strong solutions of (4.9) and (4.14) are unique under some hypotheses.

Proposition 6.1**.**

The function $p_{\nu}(t,x;y)$ given in (4.5) is the unique global solution of (4.16) (for fixed $x$) and (4.17) (for fixed $y$).

Proof.

Let us first notice that, by Theorem 5.1 we know that $0\leq p_{\nu}(t,x;y)\leq 1$ and

[TABLE]

hence $\left\|p_{\nu}(t,\cdot;y)\right\|_{\ell^{1}}=\left\|p_{\nu}(t,x;\cdot)\right\|_{\ell^{1}}=1$. Thus $p_{\nu}(t,x;y)$ is bounded in $\ell^{1}$ uniformly with respect to $t\geq 0$ for fixed $x$ or $y$. Since $\ell^{1}$ is continuously embedded in $\ell^{2}(\operatorname{\mathbf{m}})$ (see Remark 3.1), then $p_{\nu}(t,x;y)$ is uniformly bounded also in $\ell^{2}(\operatorname{\mathbf{m}})$. Moreover, we have shown in Lemma 3.2 that the operators $\operatorname{\mathcal{G}}$ and $\operatorname{\mathcal{L}}$ are continuous. Hence by [4, Corollary $2$] we can conclude that $p_{\nu}(t,x;y)$ is the unique global solution of (4.16) and (4.17). ∎

Now let us show the uniqueness of the solutions of the backward equation (4.9).

Proposition 6.2**.**

Let $g\in\ell^{\infty}$ such that $g(x)=\sum_{n=0}^{+\infty}g_{n}Q_{n}(x)$. Then the strong solution $u(t,x)$ of (4.9) is in $\ell^{\infty}$, hence also in $\ell^{2}(\operatorname{\mathbf{m}})$, for any $t\geq 0$ and it is the unique global solution in $\ell^{2}(\operatorname{\mathbf{m}})$.

Proof.

First of all, let us observe that if $g\in\ell^{\infty}$, then $g\in\ell^{2}(\operatorname{\mathbf{m}})$ too, so we are under the hypotheses of Theorem 4.3. Moreover we have, by using Theorem 5.1 and Jensen inequality:

[TABLE]

and then

[TABLE]

obtaining the uniform bound for $x\mapsto u(t,x)$. Hence $u(t,\cdot)\in\ell^{2}(\operatorname{\mathbf{m}})$ for any $t\geq 0$. Moreover, since $\operatorname{\mathcal{G}}$ is a continuous operator, by [4, Corollary $2$], it is the unique global solution of (4.9). ∎

We can also obtain the uniqueness of solutions of (4.14).

Proposition 6.3**.**

Let $f:\operatorname{\mathbb{N}}_{0}\to\operatorname{\mathbb{R}}$ be a function such that $\frac{f}{m}\in\ell^{2}(\operatorname{\mathbf{m}})$ and $\frac{f(x)}{m(x)}=\sum_{n=0}^{+\infty}f_{n}Q_{n}(x)$. Then the strong solution $u(t,x)$ of (4.14) is in $\ell^{\infty}$, hence in $\ell^{2}(\operatorname{\mathbf{m}})$, and it is the unique global solution in $\ell^{2}(\operatorname{\mathbf{m}})$.

Proof.

Let us observe that

[TABLE]

Thus, by Jensen inequality, we have

[TABLE]

Now, from Theorem 5.1, we know that $p_{\nu}(t,x;y)\leq 1$, then

[TABLE]

Finally, we have

[TABLE]

thus, since $\operatorname{\mathcal{L}}$ is a continuous operator, from [4, Corollary $2$] we have that $u(t,x)$ is the unique global solution of (4.14). ∎

Remark 6.4*.*

The condition $f/m\in\ell^{2}(\operatorname{\mathbf{m}})$ is stronger than $f\in\ell^{2}$ for any probability measure $\operatorname{\mathbf{m}}$ on $\mathbb{N}_{0}$. Indeed we can show that the following two properties

• $a)$

$f\in\ell^{2}$;

• $b)$

$f/\sqrt{m}\in\ell^{2}(\operatorname{\mathbf{m}})$;

are equivalent: this can be done simply observing that

[TABLE]

Moreover, if we consider the property

• $c)$

$f/m\in\ell^{2}(\operatorname{\mathbf{m}})$;

we can see that $c)\Rightarrow a)$. Indeed we have, since $m(x)\leq 1$

[TABLE]

However, if we consider $f(x)=\sqrt{m(x)}$, it is easy to verify that $f\in\ell^{2}$ but $f/m\not\in\ell^{2}(\operatorname{\mathbf{m}})$.

7. Limit distribution of $N_{\nu}(t)$

In this section we want to give some results on the limit distribution of $N_{\nu}(t)$. In particular we have

Theorem 7.1**.**

Let $p_{\nu}(t,x;y)$ be the transition probability mass of $N_{\nu}(t)$. Then, given any initial probability mass $f$ such that $f/m\in\ell^{2}(\operatorname{\mathbf{m}})$ and $f(x)/m(x)=\sum_{n\in\mathbb{N}_{0}}f_{n}Q_{n}(x)$, the probability mass of $N_{\nu}(t)$ asymptotically converges towards a Poisson measure, that is to say if $p_{\nu}(t,x)=\sum_{y\in\mathbb{N}_{0}}p_{\nu}(t,x;y)f(y)$, then

[TABLE]

Proof.

By Corollary 5.3 we know that $p_{\nu}(t,x)$ is solution of (4.14). Moreover, since $f/m\in\ell^{2}(\operatorname{\mathbf{m}})$, we know that this solution is unique from Proposition 6.3. Finally, from Theorem 4.4, we have that

[TABLE]

In particular we have

[TABLE]

But $Q_{0}(x)=1$ and $f_{0}=\sum_{x\in\mathbb{N}_{0}}f(x)=1$ since $f$ is a probability mass. Thus we have

[TABLE]

Now, by Cauchy-Schwartz inequality, the fact that $E_{\nu}(-bnt^{\nu})\leq 1$ and the duality formula for $C_{n}(x,\alpha)$, we obtain

[TABLE]

hence the second series totally converges. Thus we can take the limit inside the series and, since $\lim_{t\to+\infty}E_{\nu}(-bnt^{\nu})=0$, we have

[TABLE]

∎

From Theorem 7.1 we know that whatever is the distribution of $N_{\nu}(0)$, the limit distribution of $N_{\nu}(t)$ is always $\operatorname{\mathbf{m}}$. Moreover, we can show that $\operatorname{\mathbf{m}}$ is an invariant one-dimensional distribution for $N_{\nu}(t)$, that is to say that if $N_{\nu}(0)$ has distribution $\operatorname{\mathbf{m}}$, then $N_{\nu}(t)$ admits $\operatorname{\mathbf{m}}$ as distribution for any $t>0$.

Proposition 7.2**.**

Suppose $N_{\nu}(0)$ has distribution $\operatorname{\mathbf{m}}$. Then for any $t>0$, $N_{\nu}(t)$ has distribution $\operatorname{\mathbf{m}}$.

Proof.

Let us observe that the density of $N_{\nu}(t)$ is given by

[TABLE]

From Theorem 5.1 we have that

[TABLE]

then we have

[TABLE]

Recalling that $Q_{0}(y)\equiv 1$ we have that

[TABLE]

hence

[TABLE]

∎

However, since $N_{\nu}(t)$ is not Markovian, this Proposition does not guarantee the stationarity of the process when $N_{\nu}(0)$ admits $\operatorname{\mathbf{m}}$ as distribution. However, it is still possible to compute the autocovariance function of the process $N_{\nu}(t)$.

Proposition 7.3**.**

Suppose $N_{\nu}(t)$ admits $\operatorname{\mathbf{m}}$ as initial distribution. Then

[TABLE]

We omit the proof of this Proposition since it is identical to the one in [24], after observing that if $N_{1}(t)$ admits $\operatorname{\mathbf{m}}$ as initial distribution, then $N_{1}(t)$ is stationary and, from (3.3),

[TABLE]

Remark $3.2$ and $3.3$ of [24] easily apply also to our process $N_{\nu}(t)$. Indeed, since in this case $N_{\nu}(t)$ is distributed as $N_{\nu}(0)$, then the variance $\operatorname{\mathbb{D}}[N_{\nu}(t)]=\operatorname{\mathbb{D}}[N_{\nu}(0)]=\alpha$, which can be obtained from (7.1) when $t=s$ with the same calculations as in [24, Remark $3.2$]. Moreover, $N_{1}(t)$ exhibits long-range dependence, while, with the same calculations of [24, Remark $3.3$], one can show that $\operatorname{Cov}(N_{\nu}(t),N_{\nu}(s))$ decays as a power of $t$, hence it exhibits short-range dependence.