First-order linear differential equations whose data are complex random variables: Probabilistic solution and stability analysis via densities

1.   Introduction

The study of differential equations in the complex domain is a consolidated field in the general theory of differential equations ^[1,2,3,4]. The rigorous analysis of this type of differential equations involves profound tools, including Nevanlinna theory ^[5] and the Painlevé transcendents ^[6]. Some recent advances about equations in the complex domain also include difference and integral equations and their strong relationship with differential equations ^[7].

The aim of this paper is to contribute to the analysis of complex differential equations with uncertainties. According to [8, Sec. 4.7], we first point out that in the extant literature one mainly distinguishes two classes of differential equations, namely, stochastic differential equations (SDEs) ^[9,10], and random differential equations (RDEs) ^[11,12]. SDEs are those whose uncertainty is driven by a stochastic process having highly irregular trajectories (e.g., nowhere differentiable). In the case that uncertainty is modelled by the Wiener process, they are termed Itô type SDEs. Solving this class of SDEs requires using special tools as Itô calculus and, particularly, its cornerstone result, the so called Itô lemma. The formulation of these equations can be motivated from its deterministic counterpart by considering Gaussian perturbations of some model parameters, that may be affected by uncertainties. SDEs have demonstrated to be useful tools to model some real-world problems, specifically in Finance and Biology ^[13].

Complementary to SDEs, RDEs are those in which random effects are directly manifested in their coefficients, initial/boundary conditions and/or source terms (hereinafter referred to as input data), which are assumed to have a regular sample behaviour (e.g., differentiable) with respect to the independent variable (usually time and/or space) ^[11,12].

As indicated in [8, p. 96], the techniques and tools required to rigorously analyze SDEs and RDEs are distinctly different. Throughout this paper, we will deal with RDEs only. Despite the considerable influence of RDEs ^[14], to the best of our knowledge, there are not previous studies dealing with these equations in the case that their input data are complex random variables. So, as a first step, we here deal with the simplest case corresponding to the non-homogeneous first-order linear RDE

with initial condition

Here, the coefficient $a(\omega)$, the non-homogeneous term, $b(\omega)$, and the initial condition, $z_0(\omega)$, are assumed to be complex absolutely continuous random variables defined in a common complete probability space $(\Omega, \mathcal{F}, \mathbb{P})$. So, it is important to observe that we are dealing with problem (1.1)–(1.2) where the parameters are complex random variables and the solution is a complex stochastic process parametrized with the real parameter $t$. Recall that a complex random variable, say $x = x(\omega)$, is defined as $x:\Omega \longrightarrow \mathbb{C}$ that can be written as $x(\omega) = x_{R}(\omega)+\mathrm{i}\, x_{C}(\omega)$, where the real and the imaginary (or complex) contributions, $x_{R}(\omega)$ and $x_{C}(\omega)$, are real random variables, and $\mathrm{i} = \sqrt{-1}$ is the imaginary unit. In other words, any complex random variable can always be regarded as a pair whose components are two random variables, the so-called real and the imaginary parts. As a consequence, the distribution of any complex random variable can be considered as the joint distribution of two real random variables ^[15]. For convenience, hereinafter we will utilize this representation. So, the probabilistic information of a complex random variable is given by the joint probability density function (PDF), $f_{RC}(x_{R}, x_{C})$, of the random vector $(x_{R}(\omega), x_{C}(\omega))$, ^{[15,16,17,18]}. Complex random variables are applied in a number of disciplines, including information theory and digital signal processing ^[19]. Notice that, the full probabilistic analysis of the initial value problem (IVP) (1.1)–(1.2), where the parameters of the first order differential equation and the initial condition are assumed to be real absolutely continuous random variables, has been done in a previous contribution ^[20]. Further studies about random ordinary differential equation on the real domain can be found in ^[11,12].

For the sake of generality, in our subsequent analysis for the random IVP formulated in (1.1)–(1.2), we will assume that the random vector

has a joint PDF, $f_{0}(z_{0R}, z_{0C}, b_R, b_C, a_R, a_C)$, i.e. independence between input data, $z_0(\omega)$, $a(\omega)$ and $b(\omega)$ is not assumed, including the real and complex parts of everyone of the foregoing random variables. Obviously, here $z_0(\omega) = z_{0R}(\omega)+\mathrm{i} z_{0C}(\omega)$, $b(\omega) = b_R(\omega)+\mathrm{i}b_C(\omega)$ and $a(\omega) = a_R(\omega)+\mathrm{i} a_C(\omega)$.

As in the deterministic theory, a main goal when dealing with its random counterpart is to obtain the solution of the random IVP (1.1)–(1.2). However, as in our context the solution, say $z(t; \omega)$, is a stochastic process (or random function), obtaining $z(t; \omega)$ to the IVP (1.1)–(1.2) is not enough, so one further needs to determine the main probabilistic properties of $z(t; \omega)$, mainly the mean and the variance. In the complex case, the mean of $z(t; \omega)$ is defined as

where $\mathbb{E}[\cdot]$ denotes the expectation operator, and, $z_R(t; \omega)$ and $z_C(t; \omega)$ stands for the real and imaginary parts of $z(t; \omega)$, respectively. So, notice that the expectation of $z(t; \omega)$ does not exist if the expectation of its real or complex part do not. The mean function satisfies the following property

where, for $w\in \mathbb{C}$, $\overline{w}$ stands for the complex conjugate. Indeed,

Using property (1.4) and that $|z(t; \omega)|^2 = z(t; \omega) \overline{z(t; \omega)}$, it is easy to check that the variance function can be calculated by

From this property, one can obtain the following representation of the variance function of a complex stochastic process in terms of the variances of its real and complex parts (which are real stochastic processes),

As we shall see later, this formula will be more fruitful throughout our subsequent development To prove (1.6), let us observe that

In this paper, we will go a step further by determining an explicit expression of the first probability density function (1-PDF) of the solution, $f_1(z; t)$. For each $t$, $f_1(z; t)$ is defined as the joint PDF of the random vector $(z_R(t; \omega), z_C(t; \omega))$, i.e. $f_1(z; t): = f_1(z_R, z_C; t)$. Notice that the 1-PDF of the real and complex contributions, denoted by $f_{1}^R(z_R; t)$ and $f_{1}^C(z_C; t)$, respectively, can be calculated by marginalizing $f_1(z_R, z_C; t)$ as follows

Then, the expectations of the real stochastic processes $z_R(\omega; t)$ and $z_C(\omega; t)$ are

Furthermore, thanks to the representation (1.6) of the variance $\sigma_z^2(t)$ of $z(t; \omega)$, in terms of the variances, $\sigma_{z_{R}}^2(t)$ and $\sigma_{z_{C}}^2(t)$, of the real random processes, $z_R(t; \omega)$ and $z_C(t; \omega)$, respectively, we can see that $\sigma_z^2(t)$ can also be calculated as

The previous development clearly shows the importance of determining $f_1(z_R, z_C; t)$, that will be our main goal henceforth. Apart from the mean and the variance, and similarly as it happens in the real setting, the 1-PDF, $f_1(z_R, z_C; t)$, also permits calculating any one-dimensional moment of the random stochastic process, $z(t; \omega)$ at every time instant $t$, provided they exist. For $t = \hat{t}$ fixed, the PDF, $f_1(z_R, z_C; \hat{t})$, allows us to calculate the probability that the solution at $t = \hat{t}$ lies in any specific region, say $\mathcal{S}\subset \mathbb{R}^2$, of interest

Finally, we also point out that $f_1(z_R, z_C; t)$ can also be applied to construct confidence regions at any confidence level $1-\alpha$, with $\alpha \in (0, 1)$ and $t = \hat{t}$ fixed. To this end, we look for $\zeta > 0$ such that

Then, the confidence region is given by the interior points limited by the plane curve $\{(z_R, z_C)\in \mathbb{R}^2: f_1(z_R, z_C; \hat{t}) = \zeta \}$. One typically takes $\alpha = 0.05$, then $f_1(z_R, z_C; \hat{t}) = \zeta$ determines a set in the plane with $(1-\alpha)\times 100 = 95\%$ of confidence level.

To obtain an expression of the 1-PDF of the solution $z(t; \omega)$, and also of the real and complex contributions, we will extensively apply the Random Variable Transformation (RVT) technique. This method is a very powerful technique which allows us to compute the PDF of a random vector that results after transforming another random vector whose PDF is known ^{[11,21,22,23]}.

Theorem 1 (Random Variable Transformation method, ^[11]). Let $\mathbf{u}(\omega) = (u_1(\omega), \ldots, u_n(\omega))$ and $\mathbf{v}(\omega) = (v_1(\omega), \ldots, v_n(\omega))$ be two $n$-dimensional absolutely continuous random vectors. Let $\mathbf{r}: \mathbb{R}^n \rightarrow \mathbb{R}^n$ be a one-to-one deterministic transformation of $\mathbf{u}$ into $\mathbf{v}$, i.e., $\mathbf{v} = \mathbf{r}(\mathbf{u})$. Assume that $\mathbf{r}$ is continuous in $\mathbf{u}$ and has continuous partial derivatives with respect to $\mathbf{u}$. Then, if $f_{\mathbf{U}}(\mathbf{u})$ denotes the joint PDF of vector $\mathbf{u}(\omega)$, and $\mathbf{s} = \mathbf{r}^{-1} = (s_1(v_1, \ldots, v_n), \ldots, s_n(v_1, \ldots, v_n))$ represents the inverse mapping of $\mathbf{r} = (r_1(u_1, \ldots, u_n), \ldots, r_n(u_1, \ldots, u_n))$, the joint PDF of vector $\mathbf{v}(\omega)$ is given by

where $\left| J_n \right|$ is the absolute value of the Jacobian, which is defined by

The paper is organized as follows. Section 2 is devoted to obtain the 1-PDF of the solution stochastic process, $z(t; \omega)$, to the random IVP (1.1)–(1.2). For sake of completeness, in Section 2, explicit expressions for the 1-PDF for the real and imaginary components of the solution are also determined. The study is divided into two cases for the sake of clarity, namely the so-called exponential solution (Subsection 2.1) and the linear solution (Subsection 2.2). In Subsection 2.3, an analysis of the probabilistic stability of the solution is discussed in terms of the densities. Section 3 is addressed to illustrate the theoretical results obtained in Subsection 2 by means of several numerical examples. Conclusions and future work are presented in Section 4.

2.   Determining an explicit expression for the first probability density function of the solution stochastic process

As it has been indicated in Section 1, a RDE in the complex domain is a natural generalization of its deterministic counterpart. Therefore, for the first-order IVP (1.1)–(1.2), its solution is similarly obtained as in the deterministic case [2, Sec. 5.1]

where, as usual, "w.p.'' stands for "with probability''. Since $a(\omega)$ is assumed to be an absolutely continuous random variable, the first case (exponential solution) will happen w.p. 1, while the second case correspond to the linear solution and it corresponds to the case that the random IVP whose RDE is $z'(t; \omega) = b(\omega)$, $\omega \in \Omega$. For the sake of clarity, henceforth, we separately study each case.

2.1.   Exponential solution

In order to determine the 1-PDF of $z(t; \omega)$, we first rewrite (2.1) in terms of the real and complex parts of each random variable in (1.1), i.e. $z_{0R}(\omega)$, $z_{0C}(\omega)$, $b_R(\omega)$, $b_C(\omega)$, $a_R(\omega)$ and $a_C(\omega)$. From the properties of complex numbers, in particular the Euler's formula, ${\rm{e}}^{x_R+\mathrm {i} x_C} = {\rm{e}}^{x_R} \left(\cos (x_C)+\mathrm {i} \sin (x_C)\right)$, after some algebra calculations, the real and complex components of $z(t; \omega)$ can be expressed as

where

Let $t > 0$ fixed, we apply the RVT method to compute the PDF of the random vector $\mathbf{v}(\omega) = \left(v_1(\omega), v_2(\omega), v_3(\omega), v_4(\omega), v_5(\omega), v_6(\omega)\right) = (z_R(t; \omega), z_C(t; \omega), a_R(\omega), a_C(\omega), b_R(\omega), b_C(\omega))$ in terms of the known PDF of the (random) input data, $f_0(\mathbf{u})$, where $\mathbf{u} = (z_{0R}, z_{0C}, a_R, a_C, b_R, b_C)$. First, for the sake of clarity, we rewrite formulas (2.2) and (2.3) as follows

being

and $h_i$, $i = 1, 2, 3, 4$, given in (2.4). Let us define the deterministic mapping $\mathbf{r}:\mathbb{R}^6 \to \mathbb{R}^6$

where $z_R(t)$ and $z_C(t)$ are given in (2.2) and (2.3), respectively. Isolating $z_{0R}$ and $z_{0C}$ (which appear in the definitions of $z_R(t)$ and $z_C(t)$), the inverse of the mapping $\mathbf{r}$, $\mathbf{s}:\mathbb{R}^6\to \mathbb{R}^6$, is

Denoting by $\mathbf{I}_4$ the identity matrix of size 4, the Jacobian of the inverse mapping $\mathbf{s} = \mathbf{r}^{-1}$ is given by

Notice that we have used that $h_3^2+h_4^2 = {\rm{e}}^{2v_3 t}$. Then, applying the RVT technique (stated in Theorem 1.1), the joint PDF of the random vector $\mathbf{v}(\omega)$, in terms of the PDF, $f_0(\cdot)$, of the random input parameters, $\mathbf{u}(\omega)$, is

As, for each time instant $t$, we are interested in the PDF of the complex random variable $z(t; \omega)$, which is given by the joint PDF of the random vector $(z_R(t; \omega), z_C(t; \omega))$. So, it is obtained by marginalizing the joint PDF, $f_{\mathbf{v}}(\mathbf{v})$, given in formula (2.5), with respect to $(v_3(\omega), v_4(\omega), v_5(\omega), v_6(\omega)) = (a_R(\omega), a_C(\omega), b_R(\omega), b_C(\omega))$,

Then, taking $t > 0$ arbitrary, the 1-PDF of the solution stochastic process, $z(t; \omega)$, of the random IVP (1.1)–(1.2) is

Remark 1. Assuming independence between the complex RV $z_0(\omega)$, with PDF $f_{z_0}(z_{0R}, z_{0C})$, and the random vector $(a_R(\omega), a_C(\omega), b_R(\omega), b_C(\omega))$, the 1-PDF given by the integral expression $(2.6)$ can be rewritten in terms of the expectation operator

Notice that $\mathbb{E}_{(a_R, a_C, b_R, b_C)} \left[\cdot \right]$ stands for the expectation with respect to the random vector $(a_R(\omega), a_C(\omega), b_R(\omega), b_C(\omega))$. Thus, taking samples of the random vector $(a_R(\omega), a_C(\omega), b_R(\omega), b_C(\omega))$, the 1-PDF, $f_1(z_R, z_C; t)$, can be approximated by computing the mean of $f_{z_0}(s_1(z_R, z_C, a_R, a_C, b_R, b_C), s_2(z_R, z_C, a_R, a_C, b_R, b_C)){\rm{e}}^{-2a_C t}$ over all the samples. Recall that the absolute error of the approximation decreases as $1/\sqrt{N}$, being $N$ the number of simulations.

Remark 2. Recall that from expression $(2.6)$, we can also calculate the 1-PDF of the real and complex components, denoted by $f_{1}^R(z_R; t)$ and $f_{1}^C(z_R; t)$, respectively, as it was shown in $(1.7)$.

2.2.   Linear solution

In this case, the RDE writes

The solution stochastic process of the random IVP (2.7) and (1.2) is

To compute the 1-PDF, $\hat{f}_1(z; t)$, of $z(t; \omega)$, we first fix $t > 0$, and then we will determine the PDF of the random vector of the auxiliary random vector $\mathbf{v}(\omega) = (b_R(\omega) t+z_{0R}(\omega), b_C(\omega)t+z_{0C}(\omega), b_R(\omega), b_C(\omega))$ in terms of the PDF of the random vector $\mathbf{u}(\omega) = (z_{0R}(\omega), z_{0C}(\omega), b_R(\omega), b_C(\omega))$, denoted by $\hat{f}_0(\mathbf{u}) = \hat{f}_0(z_{0R}, z_{0C}, b_R, b_C)$. To this end, we define the following mapping

Denoting $\mathbf{v} = (v_1, v_2, v_3, v_4)$, then the inverse mapping $\mathbf{s}:\mathbb{R}^4\to \mathbb{R}^4$ of $\mathbf{r}$ and its Jacobian, $J_4$, are

Therefore, according to Theorem 1.1, the PDF of the random vector $\mathbf{v}(\omega)$ is

Since this argument is valid for any $t$, to obtain the 1-PDF of $z(t; \omega)$, we marginalize with respect to the random variables $b_R(\omega)$ and $b_C(\omega)$

Remark 3. Assuming independence between the complex RVs $z_0(\omega)$ and $b(\omega)$, the 1-PDF given by the integral expression $(2.9)$ can be rewritten in terms of the expectation operator with respect to the random vector $(b_R(\omega), b_C(\omega))$

where $\hat{f}_{z_0}(z_{0R}, z_{0C}) = \int_{\mathbb{R}^2} \hat{f}_0(z_{0R}, z_{0C}, b_R, b_C)\, \mathrm{d} b_R \, \mathrm{d} b_C$. Expression $(2.10)$ permits approximating the 1-PDF via the Monte Carlo method by sampling the vector $(b_R(\omega), b_C(\omega))$, $\omega \in \Omega$.

2.3.   Analysis of probabilistic stability

The study of stability of the solution stochastic process is presented for the exponential solution, since in the linear case the solution tends to infinity when $t\to \infty$. Otherwise, the solution is stable only when $\mathbb{P}[\{\omega \in \Omega: b(\omega) = 0 \}] = 1$. In such case, the equilibrium state is the initial condition, $z^{\ast}(\omega) = z_0(\omega)$, and the PDF of the equilibrium point is directly obtained via $f_{0}(z_{0R}, z_{0C})$.

The stability of the exponential solution can be deduced from the expression of the solution stochastic process $z(t; \omega) = z_R(t; \omega)+\mathrm{i} z_C(t; \omega)$, where the real and complex contributions are given in (2.2)–(2.4). We observe that the stability depends on the behaviour of the stochastic processes $h_3(a_R(\omega), a_C(\omega); t)$ and $h_4(a_R(\omega), a_C(\omega); t)$. Let $\omega \in \Omega$ be fixed, then

since it is the limit of an exponential times a bounded (trigonometric) function, and the exponential tends to zero provided the exponent is negative. Then, as we are dealing with random variables, we will have stability with a certain probability

If the solution does not converge, it diverges oscillating, since the divergence depends on the real sine and cosine trigonometric functions.

In case of stability, we observe that the distribution of the equilibrium point, which is a complex random variable, $z^*(\omega) = z_R^*(\omega)+\text{i}z_C^*(\omega)$, is determined by the joint PDF of the random vector $(z_R^*(\omega), z_C^*(\omega))$ defined by

Now, we take advantage of RVT method (Theorem 1.1) to determine an exact expression for the joint PDF of the equilibrium point given in (2.12). To this end, let us fix $\omega \in \Omega$ and define

Let $\mathbf{r}:\mathbb{R}^4 \to \mathbb{R}^4$ be the deterministic mapping

The inverse mapping $\mathbf{s}:\mathbb{R}^4 \to \mathbb{R}^4$ is

The Jacobian of $\mathbf{s}$ is $J_4 = v_1^2+v_2^2\neq 0$. Therefore, the joint PDF of the random vector $\mathbf{v}(\omega)$ is

As we are interested in the joint PDF of the random vector $(z_R^*(\omega), z_C^*(\omega))$, we marginalize the last expression with respect to $a_R(\omega)$ and $a_C(\omega)$, then obtaining the desired PDF

Notice that from the very beginning, we are assuming that we know the joint PDF, $f_0(\cdot)$, of the random vector $\left(z_{0R}(\omega), z_{0C}(\omega), a_R(\omega), a_C(\omega), b_R(\omega), b_C(\omega)\right)$. Then, the PDF of the random vector $\left(a_R(\omega), a_C(\omega), b_R(\omega), b_C(\omega)\right)$ is the marginal of $f_0(\cdot)$ with respect the random vector $(z_{0R}(\omega), z_{0C}(\omega))$. As a consequence, the PDF of the equilibrium point finally writes

Remark 4. To calculate the PDF of the equilibrium point, given by $(2.13)$, the Monte Carlo method can be applied under certain hypotheses. For example, assuming independence between the random vectors $(b_R(\omega), b_C(\omega))$ and $\left(z_{0R}(\omega), z_{0C}(\omega), a_R(\omega), a_C(\omega)\right)$, expression $(2.13)$ can be rewritten as

3.   Numerical examples

In this section, the theoretical results established in Section 2 are applied in three numerical examples. In the first example, statistical dependence between some input data is considered, while in the second and in the third examples, we assume independence between all the input data. In all the examples, after fixing of the probabilistic distributions of input data, we compute the 1-PDF, $f_1(z; t)$, of the solution stochastic process, $z(t; \omega)$; the 1-PDFs, $f^{R}_1(z_R; t)$ and $f_1^{C}(z_C; t)$, of the real and complex contributions, $z_R(t; \omega)$ and $z_C(t; \omega)$, respectively, as well as the mean and the variance of these stochastic processes. In the examples, statistical quantities of interest such as the mean, variance and confidence regions are also calculated in order to account for the probabilistic evolution of the solution over the time. Finally, the probability of stability, $p ^*$, is also determined.

In the first numerical example, we have chosen an appropriate distribution of random variable $a_R(\omega)$ to assure the solution is stable, w.p. 1, i.e. $p^* = 1$. Therefore, in this example we show that the mean and the variance of the solution tend to the one corresponding to the equilibrium state, which is a random variable. In addition, for the sake of clarity, we show that the distributions of the real and complex stochastic processes tend to the distribution of the real and complex components, respectively, of the equilibrium as time increases.

In the second example, we have chosen the distribution of the random input $a(\omega)$ so that $0 < p ^* <$ is small. Therefore, the computations show that stability will not likely take place.

In the third example, we consider the linear solution, i.e. $a(\omega) = 0$ w.p. 1. The probabilistic analysis corresponding to this case has been shown in Subsection 2.2. In this case, the solution stochastic process is unstable. We will show this behaviour by means of the PDF of the solution and its main statistics (mean, variance and confidence regions).

Finally, we point out that a wide range of probability distributions has been selected to carry out the numerical experiments to better illustrate the applicability of the theoretical results.

3.1.   Example 3.1

In this example, we consider that some input data are dependent random variables. In particular, we have chosen the following probability distributions:

● $(a_R(\omega), a_C(\omega))$ follows an Uniform distribution on the rectangle $[-10.1, -10]\times[-0.5, 0.5]$, i.e., $(a_R(\omega), a_C(\omega))\sim \text{U}\left([-10.1, -10]\times[-0.5, 0.5]\right)$.

● $(z_{0R}(\omega), z_{0C}(\omega))$ has a multivariate Gaussian distribution with mean $\mu_{z_0} = (0, 0)$ and variance-covariance matrix

Then, $(z_{0R}(\omega), z_{0C}(\omega))\sim \text{N}\left(\mu_{z_0}, \Sigma_{z_0}\right)$.

● $b_{R}(\omega)$ has a Triangular distribution on the interval $[1.7, 2.5]$ with mode at $2$, i.e., $b_R(\omega)\sim \text{Tr}([1.7, 2.5], 2)$.

● $b_{C}(\omega)$ has a Gamma distribution with shape and scale parameters $3$ and $2$, respectively, i.e., $b_C(\omega)\sim \text{Ga}(3;2)$.

In Figure 1, we show 3D graphical representations for the 1-PDF, $f_1(z_R, z_C; t)$, at the time instants $t = 0.05$ and $t = 0.2$. Calculations have been carried out from expression (2.6). Since $a_R$ lies in negative interval w.p. 1, then $z(t; \omega)\longrightarrow(z_R^{\ast}(\omega), z_C^{\ast}(\omega))$ as $t\to \infty$, where, the equilibrium point $(z_R^{\ast}(\omega), z_C^{\ast}(\omega))$, given by (2.12), is asymptotically stable. In Figure 2, we have plotted the PDF of $(z_R^{\ast}(\omega), z_C^{\ast}(\omega))$ using (2.13). Comparing the plots shown in Figures 1 and 2, we can graphically observe that convergence of PDF $f_1(z_R, z_C; t)$ to $f_1(z_R^{\ast}, z_C^{\ast})$ is fast. Finally, in Figure 3, we show the evolution of the mean of the solution $\mathbb{E}[z(t; \omega)]$ (blue points) at different time instants, $t\in \{0, 0.05, 0.1, 0.2 \}$, together with $50\%$ and $90\%$ confidence regions. They have been computed by expressions (1.3) and (1.8) for the mean, and (1.10) with $\alpha = 0.5$ and $\alpha = 0.1$ for the confidence regions (red and black solid curves, respectively). The same statistics are also shown for the equilibrium point $(z_R^{\ast}(\omega), z_C^{\ast}(\omega))$ (black point). We can observe that both the mean and confidence regions of $z(t; \omega)$ rapidly approximate the corresponding ones (plotted by dashed curves) associated to $(z_R^{\ast}(\omega), z_C^{\ast}(\omega))$ as $t$ increases, thus illustrating stability. It is interesting to observe the full agreement between the shapes and scales of confidence regions observed in Figure 3 at the time instants: $t = 0.05$ and Figure 1 (left panel); $t = 0.2$ and Figure 1 (right panel), and $t = \infty$, corresponding to the equilibrium (Figure 2). We complete the graphical analysis by plotting in Figure 4 the evolution of the PDF of the real component, $f_1^R(z_R; t)$ (left panel), and of the imaginary component, $f_1^C(z_R; t)$ (right panel), at the time instants $t\in \{0, 0.05, 0.1, 0.15, 0.2 \}$, together with the corresponding PDFs of the components of the equilibrium. From these plots we can graphically confirm the stability, as expected. Again, it must be pointed out that the results shown in Figures 1–3 and in Figure 4 are in full agreement. In Figure 5 we show the expectation and the variance of the real and complex parts of the solution stochastic process together with the corresponding statistics of the equilibrium state. Observe that these graphical results completely agree with the ones shown in Figure 4.

3.2.   Example 3.2

In this example we assume independence between all the input data. Specifically, we consider the following probability distributions for each random variable:

● $a_R(\omega)$ has an Uniform distribution, defined on the interval $[-0.02, 0.4]$, i.e., $a_R(\omega)\sim \text{U}([-0.02, 0.4])$.

● $a_C(\omega)$ has a Beta distribution with shape parameters $2$ and 5, i.e., $\text{Be}(2;5)$.

● The real component of the initial condition has a triangular distribution with mode 0 and defined on the interval $[-0.2, 0.1]$, i.e., $z_{0R}(\omega)\sim \text{T}([-0.2, 0.1];0)$.

● The complex component of the initial condition has a Gamma distribution with parameters 2 and $0.2$, i.e., $z_{0C}(\omega)\sim \text{Ga}(2;0.2)$.

● The parameter $b_R(\omega)$ has a Gaussian distribution with mean $-1$ and standard deviation $0.01$ truncated on the interval $T = [-2, 0]$, i.e., $b_R(\omega) \sim \text{N}_{T}(-1;0.01^2)$.

● $b_C(\omega)$ is uniformly distributed on the interval $[-2, -1.8]$, i.e., $b_C(\omega)\sim \text{U}([-2, -1.8])$.

Because of independence of input data, the joint PDF of them can be expressed as the product of the corresponding marginals PDFs, i.e.

In Figure 6, we have plotted the 1-PDF of the solution stochastic process, $f_1(z_R, z_C; t)$, at the time instants $t = 0.5$ and $= 1.5$. In Figure 7, we have represented the mean (blue points) and confidence regions at $50\%$ (red curves) and at $90\%$ (blue curves) at the time instants $t\in \{0, 0.5, 1, 1.5, 2 \}$. From this plot we can see that variability increases over the time. The PDFs of the real and imaginary components, corresponding to the above time instants, are shown in Figure 8. Notice that the aforementioned increasing of the variability as time goes on, is in full agreement with these graphical representations. This feature indicates that in the long run the solution tends to become unstable. Indeed, the probability of reaching stability can be exactly computed via (2.11)

which is small. This situation is illustrated through the mean and the variance of the real and complex components, which are shown in Figure 9. In Figure 10, the same statistics are shown over the larger domain $[0, 120]$.

3.3.   Example 3.3

As previously indicated, in this example we shall assume that $a(\omega) = 0$ w.p. 1, so the linear solution, given in (2.8), is considered. We will choose the same probability distributions as in Example 3.2 for the input data, which are also assumed to be independent. In Figure 11, we show 3D graphical representation of the 1-PDF of the solution process, $f_1(z_R, z_C; t)$ at the time instants $t = 2$ (left) and $t = 10$ (right). We can see that the PDF moves the coordinates of its real and imaginary part towards greater (in absolute value) numbers as times increases. This fact is consistently also observed in Figure 12, where the mean (blue points) and $50\%$ (red curves) and $90\%$ (black curves) confidence regions have been plotted at the time instants $t \in \{0, 1, \ldots, 10 \}$. From this plot we also observe that the variability increases as $t$ does. In Figure 13, the evolution of the PDFs of the real (left) and imaginary (right) components of the solution over the aforementioned time instants are shown. Notice that the plots are in full agreement with the ones shown in Figures 11 and 12. We complete the graphical analysis showing, in Figure 14, the plots of the mean (left panels) and of the variance (right panels) of the real and imaginary parts of the solution stochastic processes. It is clearly observed from these statistics that the solution is unstable. {Observe that it is in full agreement with the fact that this example corresponds to the case termed "linear solution'' (see Section 2.2) whose solution is a linear stochastic process (see expression (2.8)) such that $\mathbb{P}[\{\omega \in \Omega: b(\omega) = 0 \}]\neq 1$ (see also the discussion at the beginning of Section 2.3), so the solution diverges as $t$ tends to infinite, i.e., is unstable.}

4.   Conclusion

In this paper, we have performed a probabilistic study of first order linear differential equations whose coefficients and the initial condition are complex random variables. The solution is a complex valued stochastic process parametrized with respect to a real parameter (identified with the time as usual). We have provided an exact expression for the probability density function of the solution stochastic process under very general hypotheses on the model parameters. Furthermore, our study has included a probabilistic analysis of the stability. Despite the simplicity of the differential equation, to the best of our knowledge, this is the first contribution in this direction and the ideas exhibited throughout the paper can open new avenues for studying other families of differential equations with uncertainties.

Acknowledgments

This work has been supported by the Spanish Ministerio de Economía, Industria y Competitividad (MINECO), the Agencia Estatal de Investigación (AEI) and Fondo Europeo de Desarrollo Regional (FEDER UE) grant PID2020-115270GB-I00 and Generalitat Valenciana (Grant AICO/2021/302).

We give our deepest thanks to the four anonymous referees that have improved the final version of the paper with their valuable comments.

Conflict of interest

The authors declare that there are not conflict of interest.