On the Entropy of Fractionally Integrated Gauss–Markov Processes

: This paper is devoted to the estimation of the entropy of the dynamical system { X α ( t ) , t ≥ 0 } , where the stochastic process X α ( t ) consists of the fractional Riemann–Liouville integral of order α ∈ ( 0, 1 ) of a Gauss–Markov process. The study is based on a speciﬁc algorithm suitably devised in order to perform the simulation of sample paths of such processes and to evaluate the numerical approximation of the entropy. We focus on fractionally integrated Brownian motion and Ornstein–Uhlenbeck process due their main rule in the theory and application ﬁelds. Their entropy is speciﬁcally estimated by computing its approximation (ApEn). We investigate the relation between the value of α and the complexity degree; we show that the entropy of X α ( t ) is a decreasing function of α ∈ ( 0, 1 ) .


Introduction
In the study of a biological system, whose time evolution is modeled by a stochastic process that depends on a certain parameter α, often there is a need to find how a change in the value of α affects the qualitative behavior of the system, as well as its complexity degree, or entropy. Another useful information is the knowledge of a stochastic ordering, with respect to expectation of functionals of the process (e.g., its mean and variance), when varying α.
As a case study, we are interested to the qualitative behavior of the fractional integral of a Gauss-Markov (GM) process, when varying the order α of the fractional integration.
Actually, GM processes and their fractional integrals over time are very relevant in various application fields, especially in Biology-e.g., in stochastic models for neuronal activity (see [1]). In particular, the fractional integral of order α ∈ (0, 1) of a GM process, say X α (t), is suitable to describe certain stochastic phenomena with long range memory dynamics, involving correlated input processes (see [2]).
As an example of application, one can consider the model for the neuronal activity, based on the coupled differential equations: Here, D α stands for the Caputo fractional derivative (see [3]); η(t) is in place of the white noise, usually utilized in the stochastic differential equation, which describes a Leaky Integrate-and-Fire (LIF) neuronal model (see, for example, [4]). The colored noise process η(t) is the correlated process obeying the second of Equation (1) and it is the input for the first one; it is indeed a time-non-homogeneous GM process of Ornstein-Uhlenbeck (OU)-type (see Section 2). The stochastic process V(t) represents the voltage of the neuronal membrane, whereas C m is the membrane capacitance, g L the leak conductance, V L the resting (equilibrium) level potential, I(t) the synaptic current (deterministic function), τ is the correlation time of η(t) and B(t) the noise (standard BM). As we can see, the process V(t), which is the solution of (1) belongs to the class of fractional integrals of GM processes. Indeed, it is a specific example of X α (t) process, V(t) being the Caputo fractional integral of the η(t) process [5]. The biophysical motivation in the above model is to describe a neuronal activity as a perfect integrator (without leakage), from an initial time until to the current time, of the process η(t), representing the time dependent input. The use of fractional operators allows us to regulate the time scale by choosing the fractional order of integration suitably adherent to the neuro-physiological evidences. Indeed, such a model can be useful, for instance, in the investigation and simulation of synchronous/asynchronous communications in networks of neurons [6].
To introduce the terms of our investigation, we recall some definitions. A continuous GM process Y(t) is a stochastic process of the form: where is a monotone increasing, differentiable and non-negative function.
For a continuous function f (t), its Riemann-Liouville (RL) fractional integral of order α ∈ (0, 1) is defined as (see [7]): where Γ is the Gamma Euler function-i.e., Γ(z) = +∞ 0 t z−1 e −t dt , z > 0. We recall also that the Caputo fractional derivative of order α of a function f (t) is defined by (see [3]): where f denotes the ordinary derivative of f . Notice that, taking the limit for α → 1, one gets Referring to the neuronal model (1), assuming that V(0) = 0 (and, in some cases, also η(0) = 0), the RL fractional integral I α is used as the left-inverse of the Caputo derivative D α (see [8,9]). In this way, we find that the solution V(t) of (1) involves the RL fractional integral process of the GM process η(t), specifically: Thus, V(t) turns out to be written in terms of the fractional integral of η(t).
From this consideration, in the framework of general stochastic models involving correlated processes, it appears useful to investigate the properties of X α (t) = I α (Y)(t),-i.e., the fractional integral of a GM process Y(t), as varying α ∈ (0, 1). Although X α (t) is not Markov, we have showed in [2] that it is still a Gaussian process with mean µ α (t) and variance σ 2 α (t); for instance, the fractional integral of BM has mean µ α (t) = 0 and variance (for closed formulae of the mean µ α (t) and variance σ 2 α (t) of the fractional integral of a general GM process, see [2]). For fixed α, σ 2 α (t) turned out to be increasing, as a function of t. Moreover, in [2] we found that for small values of time t the variances of the considered fractionally integrated GM processes become ever lower, as α increases (i.e., the variance decreases as a function of α); for large values of t this behavior is overturned, and the variance increases with α (see [2]).
In this paper, we aim to characterize the qualitative behavior of the dynamical system X α (t), α ∈ (0, 1) by means of its entropy. Indeed, the entropy is widely used for this purpose in many fields (see [10][11][12][13][14]). In Biology, entropy is useful to characterize the behavior of, for example, Leaky Integrate-and-Fire (LIF) neuronal models (see [4]). In finance, Kelly in [15] introduced entropy for gambling on horse races, and Breiman in [16] for investments in general markets. Finally, the admissible self-financing strategy achieving the maximum entropy results in a growth optimal strategy (see [17]).
In order to specify the entropy for the processes considered in this paper, we first note that, for a fixed time s the r.v. X α (s) is normally distributed with mean µ α (s) and variance σ 2 α (s), so recalling that the entropy of a r.v. X with density f (x) is given by where, by calculation, it easily follows that the entropy of the normal r.v. X = X α (s) with fixed s, depends only on σ 2 α (s) and it is given by (see [18], p. 181): Thus, the larger the variance σ 2 α (s), the larger the entropy of X α (s) for a fixed time s. In this paper we are interested in studying a different quantity: for a certain value of α ∈ (0, 1), and T > 0, our aim is to find the entropy of trajectories of X α (t), t ∈ [0, T], which involves all the points of the trajectories up to time T, and to show that the entropy is a decreasing function of α.
We do not actually compute the entropy of X α (t), but its approximate entropy ApEn (see [19]), obtained by using several long enough simulated trajectories (they were previously obtained in [2], for the fractional integral of some noteworthy GM processes Y(t), namely, BM and Ornstein-Uhlenbeck (OU)). In fact, Pincus [19] has showed that ApEn is suitable to quantify the concept of changing complexity, being able to distinguish a wide variety of system behaviors. Indeed, for general time series, it can potentially separate those coming from deterministic systems and stochastic ones, and those coming from periodic systems and chaotic ones; moreover, for a homogeneous, ergodic Markov chain, ApEn coincides with Kolmogorov-Sinai entropy. Thus, though X α (t) is not a Markov process, its approximate entropy ApEn is able to characterize the complexity degree of the system, when varying α.
As we said, we previously found that, in all the considered cases of GM processes, for large t the variance σ 2 α (t) of their fractional integral X α (t) is an increasing function of α, while for small t it decreases with α; instead, the covariance function has more diversified behaviors (see [2]).
In the present article, we show that, for small values of α ∈ (0, 1), X α (t) exhibits a large value of the complexity degree; a possible explanation is that, for small α the trajectories of the process X α (t) become more jagged, giving rise to a greater value of the complexity degree. In fact, our estimates of ApEn show that it is a decreasing function of α ∈ (0, 1). This behavior appears for the fractional integral of BM (FIBM), as well as for the fractional integral of the OU process (FIOU).

The Entropy of the Trajectories of X α (t)
In this section, we study the complexity degree of the trajectories of the process X α (t), in two noteworthy cases of GM processes Y(t), precisely: is the Ornstein-Uhlenbeck (OU) process, driven by the SDE.
which can be expressed as (see [20]): where the equality is meant in distribution, and OU process Y(t) is a GM process of the form (2), with: and covariance Then, X α (t) = I α (Y)(t) is called the fractionally integrated OU (FIOU) process. Both FIBM and FIOU are Gaussian processes whose variance and covariance functions were explicitly obtained in [2] and studied, as functions of α ∈ (0, 1).
To study the complexity degree of the trajectories of the process X α (t), in cases (i) and (ii), we make use of several simulated trajectories of length N, previously obtained in [2], for N large. The sample paths have been obtained by using the R software, with time discretization step h = 0.01 and by means of the same sequence of pseudo-random Gaussian numbers. The simulation algorithm has been realized as an R script. More specifically, we specialize the algorithm to simulate an array of (x 1 , x 2 , . . . , x N ) Gaussian numbers with a specified covariance matrix. Indeed, we first set the time instants t 1 < t 2 < . . . < t N (with t 0 = 0 and t i = t i−1 + h, i = 1, . . . , N) and we evaluate the elements of the covariance matrix C i,j = cov(X α (t i ), X α (t j )). Note that, for each fractionally integrated Gauss-Markov process here considered, we implemented a specific algorithm to be evaluated by numerical procedures the mathematical expression of the covariance according to Equation (3.5) of [2]. Then, we apply the Cholesky decomposition to matrix C in order to determine the lower triangular matrix G, such that C = GG T , where G T is the transposition of G. Finally, we generate N pseudo-Gaussian standard numbers (z 1 , z 2 , . . . , z N ) ≡ z T and we set x i = G i z (for i = 1, . . . , N, with G i the i−th row of matrix G) such that the obtained array (x 1 , x 2 , . . . , x N ) is a simulation of a centered Gaussian distributed N-dimensional r.v. with covariances cov(x i , x j ) = C i,j for i, j = 1, . . . , N.
In particular, referring to algorithms for the generation of pseudo-random numbers (see [21]), the main steps of implementation were the following (for more, see [2]): of an equi-spaced temporal grid. STEP 2 The Cholesky decomposition algorithm is applied to the covariance matrix C in order to obtain a lower triangular matrix G(i, j), such that C = GG T . STEP 3 The N-dimensional array z of standard pseudo-Gaussian numbers is generated. STEP 4 The sequence of simulated values of the correlated fractionally integrated process is constructed as the array x = Gz.

The Approximate Entropy
In [19] Pincus defined the concept of approximate entropy (ApEn) to measure the complexity of a system, proving also that, for a Markov chain, ApEn equals the entropy rate of the chain. In fact, to measure chaos concerning a given set of data, we have at our disposal Hausdorff and correlation dimension, K-S entropy, and the Lyapunov spectrum (see [19]); indeed, to calculate one of the above parameters, one needs an impractically large amount of data. Instead, calculation of ApEn(m, r) (see below for the definition) only requires relatively few points. Actually, as shown in [19], if one uses only 1000 points, and m is taken as being equal to 2, ApEn(m, r) can characterize a large variety of system behaviors, since it is able to distinguish deterministic systems from stochastic ones, and periodic systems from chaotic ones.
For instance, Abundo et al. [10] used ApEn to obtain numerical approximations of the entropy rate, with the final purpose to investigate the degree of cooperativity of proteins in a Markov model with binomial transition distributions. They showed that the corresponding ApEn is a decreasing function of the degree of cooperativity (for more about approximation of entropy by numerical algorithms, see [12] and references therein). Now, we recall from [19] the definition of ApEn. Let {x 1 , x 2 , ..., x N } be given a time-series of data, equally spaced in time, and fix an integer m > 0 and a positive number r. Next, let us consider a sequence of vectors {v 1 in which the distance d(·, ·) between two vectors is defined by We observe that the C i (m, r) quantities measure up to a tolerance r the frequency of patterns which are similar to a certain pattern whose window length is m. Now, define and Given N data points, formula (16) can be implemented by defining the statistics Heuristically, we can say that ApEn is a measure of the logarithmic likelihood that runs of patterns that are close for m observations, remain close on the next incremental comparison. A greater likelihood of remaining close (i.e., regularity) produces smaller ApEn values, and viceversa. On the basis of simulated data, Pincus showed that, for N = 1000 and m = 2, for values of r, between 0.1 and 0.2 times the standard deviation of the x i data produce reasonable statistical validity of ApEn(m, r, N). Moreover, he showed that, for a homogeneous, ergodic Markov chain, ApEn coincides with the Kolmogorov-Sinai entropy (see [14]), that is where p ij denotes the transition probability of the Markov chain from the state i to the state j, and π j = lim n→∞ p (n) ij is the j−th component of the vector π = (π 1 , π 2 , . . . ) of the stationary probabilities, being p (n) ij the n−step transition probability of the Markov chain from the state i to the state j.

Calculation of the Entropy of Simulated Trajectories of the Process X α (t)
In the case of FIBM and FIOU, for a set of values α ∈ (0, 1), we have performed L (discretized) trajectories (x 1 , x 2 , . . . , x N ) of length N of the process X α (t), by means of the simulation algorithm previously described in STEPS 1-4. In particular, for each simulated path, we follow the remaining steps: STEP 5 Construction of the array {v 1,α , v 2,α , ..., v N−m+1,α } in R m (for a fixed m) by extracting from a given sample path (x 1 , x 2 , . . . , .., X α (t i+m−1 )). STEP 6 Construction of the distance matrix D α i,j whose elements are d α i,j are defined as the follows distance between vectors v i,α and v j,α -i.e., STEP 7 After setting r = 0.1 * S, with S sample deviation of simulated paths, evaluation of array C α whose components are provided as and We have taken the number of paths L large enough and N from 100 to 300, and for each of these L trajectories of length N, corresponding to a value of α, we have estimated ApEn α (i), i = 1, . . . , L by means of the approximation ApEn(2, r, N), where r = 0.1× (the standard deviation of trajectory points); then, the approximate entropy of X α (t) has been obtained by ApEn L α = 1 L ∑ L i=1 ApEn α (i). This allowed us to study the dependence of the entropy of the sample paths of X α (t) = I α (Y)(t) on the parameter α, showing that the entropy-namely a measure of the complexity of the dynamical system X α (t)-is a decreasing function of α ∈ (0, 1).
Since the fractional integral of order zero of Y(t) is nothing but the process Y(t) itself, and the fractional integral of order 1 is the ordinary Riemann integral of Y(t), our result means that fractional integration introduces a greater degree of complexity than that corresponding to ordinary integration; moreover, the maximum degree of complexity is obtained for the original process Y(t) (that is, without integration).
In Figures 1 and 2 we plot the numerical results for ApEn, as a function of α, for FIBM and FIOU, respectively. When the estimates of ApEn have been obtained for N = 100, it appears clear that ApEn is a decreasing function of α.
Moreover, our calculation highlights that, for small values of α, the trajectories of FIBM and FIOU become more jagged, giving rise to a greater value of the complexity degree (see Figure 3).  We also show that the results of ApEn as N increases in Figures 4 and 5. Our investigations show that the estimated values of ApEn for FIOU, for a given α and a given trajectory length, are considerably larger than those for FIBM (compare Figures 4 and 5). This possibly depends on the fact that the trajectories of FIOU are more complicated than those of FIBM, giving rise to a greater complexity degree. Moreover, contrary to the case of FIBM, where for all α the estimated value of ApEn is a decreasing function of the length N of simulated trajectories, in the case of FIOU, for α ≤ 0.5, the estimated value of ApEn appears to be an increasing function of N. Perhaps if one used far longer trajectories (N ≥ 1000) to estimate ApEn, the values obtained in both cases would be comparable and they would exhibit the same behavior as a function of N. Notice, however, that to simulate very long trajectories is impractical from the computational point of view (even for N = 300, the CPU time to evaluate ApEn in the case of FIOU was of order of almost one hour).

Conclusions and Final Remarks
In this paper, we further investigated the qualitative behavior of the fractional integral of order α ∈ (0, 1) of a Gauss-Markov process, that we already studied in [2].
Actually, Gauss-Markov processes and their fractional integrals over time are very relevant in various application fields, especially in Biology-e.g., in stochastic models for neuronal activity (see [1]). In fact, the fractional integral of order α ∈ (0, 1) of a Gauss-Markov process Y(t), say X α (t), is suitable to describe stochastic phenomena with long range memory dynamics, involving correlated input processes, which are very relevant in Biology (see [2]).
While in [2] we have showed that X α (t) is itself a Gaussian process, and we have found its variance and covariance, obtaining that the variance σ 2 α (t) of X α (t) is an increasing function of α, in this paper we have characterized the qualitative behavior of the dynamical system X α (t), α ∈ (0, 1) by means of its complexity degree, or entropy. Actually, for several values of α we have estimated its approximate entropy ApEn, obtained by long enough trajectories of the process X α (t). Specifically, we investigate the problem by means of the implementation of an algorithm based on STEPS 1-8 detailed described in the paper. We have found that ApEn is a decreasing function of α; this behavior appeared for the fractional integral of the Brownian motion, as well as for the fractional integral of Ornstein-Uhlenbeck process. Since the fractional integral of Y(t) of order zero is nothing but the process Y(t) itself, and the fractional integral of order 1 is the Riemann integral of Y(t), our result means that fractional integration introduces a greater degree of complexity than in the case of ordinary integration; moreover, the maximum degree of complexity is obtained for the original Gauss-Markov process Y(t) (that is, without integration).
Furthermore, we remark that the algorithm for computing ApEn uses numerical data, which can be used independently of knowing the process where they come from. However, in our case, we study the process X α (t), when varying the parameter α, so we need to simulate its trajectories, and to make use of the obtained numerical values to estimate ApEn. As regards the possibility of finding out, by using ApEn, if certain data come from a particular class of possible systems, we have not investigated this. Our aim was only to characterize the behavior of fractionally integrated Gauss-Markov process X α (t), as varying the parameter α, by means of the corresponding value of ApEn.
As a future work, we aim to estimate the entropy for other cases of fractionally integrated Gauss-Markov processes X α (t), such as the fractional integral of stationary Ornstein-Uhlenbeck process. Moreover, in order to further characterize the qualitative behavior of X α (t) in terms of α, our investigation will be addressed to estimate the fractal dimension of its trajectories, as a function of α.