Thermodynamics of Superdiffusion Generated by Lévy–Wiener Fluctuating Forces

Scale free Lévy motion is a generalized analogue of the Wiener process. Its time derivative extends the notion of “white noise” to non-Gaussian noise sources, and as such, it has been widely used to model natural signal variations described by an overdamped Langevin stochastic differential equation. Here, we consider the dynamics of an archetypal model: a Brownian-like particle is driven by external forces, and noise is represented by uncorrelated Lévy fluctuations. An unperturbed system of that form eventually attains a steady state which is uniquely determined by the set of parameter values. We show that the analyzed Markov process with the stability index α<2 violates the detailed balance, i.e., its stationary state is quantified by a stationary probability density and nonvanishing current. We discuss consequences of the non-Gibbsian character of the stationary state of the system and its impact on the general form of the fluctuation–dissipation theorem derived for weak external forcing.

The normal distribution plays a special role in statistics and physics. This is due to the abundance of the observations displaying Gaussian fluctuations which is explained by the central limit theorem. Additionally, the Gaussian assumption in tandem with the Markov property, in many cases, makes calculations more manageable. Therefore, Gaussian white noise is used as an archetypal process modeling complex interactions of a test particle with its environment. This approximation works perfectly well when interactions are independent and bounded. Both these assumptions can be violated resulting in more general non-Markovian and non-Gaussian processes. Also, vast experimental evidence indicates that fluctuations can be of a more general, non-Gaussian type.
Many collective, complex physical systems are characterized by fluctuations strongly deviating from the "canonical" Gaussian description, often categorized by diverging mean and variance. In terms of the probability calculus, distributions governing this type of behavior are frequently identified with so called Lévy α-stable laws [1], where "stability" signifies that the product of characteristic functions of two independent such laws results in a characteristic function of another variable of the same type. Moreover, according to the generalized central limit theorem [2], these distributions appear as limits of sums of independent and identically distributed random variables with divergent variance, i.e., random variables with probability density functions behaving as the power law 1/|x| α+1 for large |x| values.
Lévy statistics are omnipresent and have been used as valuable models of various data sets-they have been detected in the critical state and in self-organized criticality phenomena. They correctly fit empirical distributions of extreme events, like avalanches and earthquakes, and describe the broadening of spectral lines in plasma, the distribution of gravitational force between randomly localized stars, and relaxation phenomena in disordered systems well [3]. A generalization of Gaussian white noise to its non-Gaussian, Lévy stable counterpart serves as a model of impulsive, large scale variations observed, e.g., in turbulent heat flow [4] and solar flare fluctuations [5,6], hole transport in semiconductors [7], transmission of light in inhomogeneous materials [8], and anomalous diffusive transport [9,10]. Apart from financial mathematics where the Lévy fluctuations became an attractive model for price variation [11], there is also accumulating evidence from biological experiments showing that production of mRNA and proteins occurs in a pulsatile manner and creates non-Gaussian noise in individual cells [12][13][14]. Those bursty events result in high transcriptional activity followed by long periods of inactivity and are characterized by heavy tailed distributions and Lévy-like statistics. Moreover, the cytoplasmatic mechanical activity has been documented to be far from equilibrium [15], and the total intensity of cytoskeletal noise has been estimated to exceed the level of thermal noise (k B T). Since the intrinsic stochastic excitations may play a crucial role in transcriptional regulatory systems [14], it could well be that non-Gaussian Lévy noise should be a proper model of choice for underlying fluctuations in biological systems.
Over the past years, stochastic differential equations with non-Gaussian Lévy noises have gained a lot of interest in regard to modeling the DNA-target search for binding sites [16], active transport within cells [17], and search strategies [18][19][20][21][22][23][24]. In all of those problems, the primary focus is on proposing a stochastic model to address the issues of relaxation and kinetics of the system under investigation. The paradigmatic choice of preference is usually the Langevin equation, for which stochastic energetics have been defined under the action of Gaussian white noise [25,26]. It thus seems plausible to deeply explore the possibilities and limitations of the Langevin approach to study systems with a more general form of fluctuating forces. Before addressing this point, we briefly review the properties of α-stable laws [1,[27][28][29][30][31][32] and discuss the profoundly nonequilibrium properties of Lévy-stable noises. Definition 1. Let X 1 and X 2 be independent identically distributed random variables. The distribution they share is said to be stable if where F X stands for the distribution function of X. If d = 0 for every pair (a, b), then the distribution is said to be strictly stable.
The family of all distributions which satisfies the requirement of stability was constructed by Paul Lévy [27]. These distributions are referred to as α-stable distributions or Lévy α-stable distributions.
In accordance with Samorodnitsky [1], we recall here a representation of standard α-stable distributions with the following parameters: α (identified as a characteristic exponent or a stability index), β (responsible for the skewness of the probability law), and σ-called the scale parameter that control thes overall distribution width. Finally, µ is the location parameter. The stability index α describes the asymptotic (i.e., large x) behavior of the density function (p(x) = d dx F X (x)), which for the symmetric case with α < 2 assumes the form Note that for α < 2, any α-stable density is characterized by the diverging standard deviation. Moreover, for α < 1, the mean value also diverges.
Following this definition, the formal time derivative of the symmetric α-stable motion defines a generalization of the Gaussian white noise, i.e., a symmetric Markov α-stable noise which turns into a standard Gaussian form for α = 2. Lévy flights thus extend the Brownian motion paradigm to self-similar motions ({X(ct), t 0} and {c 1/α X(t), t 0} have the same distributions) with uncorrelated random steps. However, unlike standard Brownian motions for which the mean-squared displacement (MSD) grows in time ( ∆X 2 (t) ∝ t), the dispersion of the position in the Lévy motion diverges, and the width of the resulting asymptotic Lévy (super)-diffusion must be characterized by some fractional moments [30] or the interquantile distance (see Figure 1). This fact indicates that the Lévy flight, i.e., the stochastic process defined by Equation (5), displays some non-physical properties, which are especially visible for α < 2. Due to the heavy tails of α-stable densities, with α < 2, there is a significant probability of extremely long jumps leading to the divergence of variance of Lévy flights. The resulting infinite propagation velocity can be eliminated through the spatio-temporal coupling of jump lengths and associated waiting or travel times resulting in so-called Lévy walks [33,34]. The non-Markovian character of Lévy walks, however, impedes the stochastic thermodynamics analysis of these systems. Here, we focus on the study of fluctuation dissipation relations for Lévy flights and direct the interested reader to the further discussion of differences and analogies between Lévy walks and Lévy flights in Ref. [35].
Lévy flights and Lévy fluctuations are ubiquitous in nature. They occur in turbulent flows [36], incoherent radiation trapping [37] diffusion of particles in random media [38], the spreading of epidemics and human travel behavior [39], in economic and paleoclimatic time series [40,41] and in random movements of the cell cytoskeleton generated by motor proteins [42]. The evidence for bursty, Lévy flight-like behavior in human cognition retrieval from semantic memory and mental searches have been recently documented [18,19,43]. Despite their omnipresence, the thermodynamics of non-Gaussian stable fluctuations poses several objectives [44,45] beyond the standard description typical for Gaussian cases. In the forthcoming brief resume, we comment on some unusual properties of Lévy fluctuations and their thermodynamic consequences. Along with ten exemplary sample trajectories, the median (red) and quantile areas (gray and violet) are plotted. The quantiles were obtained from 10 5 sample trajectories integrated with the Euler method (∆t = 0.01).

Lévy Flights and Detailed Balance
Let (X t ) t 0 be a stochastic process that is solved by the following Langevin equation: where dS (α) t stands for increments of the Lévy-Wiener process (Lévy motion) defined above. Equation (5) describes the overdamped motion of the random walker in the static harmonic potential subjected to the action of the α-stable (Lévy) noise. Since the system is overdamped, its state is fully determined by position x(t) which is a random variable. Alternatively, the probability density function (PDF) of the process p(x, t) obeys [30,[45][46][47] which is conditioned on the value of the stochastic process at time t 0 = 0 (the process is time homogeneous). Parameter a sets the time constant of the deterministic part of the system. For the sake of brevity, in the following, we put a = 1. The fractional operator ∂ α ∂|x| α with 0 < α 2 stands for the Riesz fractional derivative, defined by the Fourier transform [28,29] Accordingly, in the (fractional) diffusion equation (Equation (6)), the term σ α 0 can be interpreted as the generalized diffusion constant [3,48], where σ 0 is the scale parameter that characterizes the underlying Lévy noise in the overdamped Langevin Equation (5). The generalized diffusion coefficient is related to the distribution asymptotics (see Equation (4) and Ref. [48]). Please note, that in situations when finite propagation velocity is assumed, one can calculate the constrained or "moving" mean squared displacement [49] which differs from the generalized diffusion coefficient (compare with Refs. [48,49]). The solution of this SFPE is given by where S α (·) is the standard symmetric α-stable distribution whose Fourier transform is The time-dependent location (µ(t)) and scale (σ(t)) parameters read and The unique stationary state is denoted p s (x) = lim t→∞ p s (x, t|x 0 ). The stationary state is given by the α-stable density with a different scale parameter than in the noise term, i.e., σ ∞ = σ 0 /α 1/α , which is the t → ∞ limit of Equation (10).
A dynamic system driven by fluctuating forces (Equation (5)) is said to be microscopically reversible if the probability of its trajectory and its time reversal are identical, so that the work imposed by external forcing equals the change in free energy [50,51]. Here, we examine the flow of probability between the states and show that the relative probabilities of trajectories in the Lévy-Wiener process do not fulfill the detailed balance condition.

Definition 4.
We say that there is a detailed balance if the following condition is met: ∀ x,y∈R,t>0 : p(x, t|y)p s (y) = p(y, t|x)p s (x). (11) Theorem 1. The detailed balance holds for the solution of (5) if α = 2.
Proof. The following calculations show that the detailed balance condition (11) is trivially fulfilled for α = 2. (For the sake of clarity we have assumed σ 0 = 1/ √ 2 here which corresponds to the standard normal density N(0, 1), see Equation (2). It is trivial to prove that a change in the parameter does not change the result (Theorem 1)).
Now, let us turn to the case with α < 2. From (11) it follows that, in particular, and thus Conditions (12) and (13) are necessary, but not sufficient, for the detailed balance to hold. For α < 2 and |x| > 0, we have (note that lim t→0 Intuitively, this means that for very short times, the propagator does not depend on the deterministic force felt by the particle (see Equation (5)). By plugging this into (12) and (13), we arrive at (x > 0): which shows that the detailed balance does not hold for α < 2 (note that for a symmetric stable variable, h(x) takes its maximum value at x = 0 and decays beyond this value). This observation can be easily checked for the Cauchy-Wiener process (α = 1) whose propagator (p(x, t|x 0 )) satisfying Equation (6) reads Figure 2 displays the sample ratios (g(x, t)) for this case and various values of x. For small t and moderate values of x, clear deviations from g(x, t) = 1 are visible. A similar result was obtained a particle subjected to Cauchy (α = 1) white noise in the parabolic potential well in Ref. [52], where the nonzero flow of probability between the segments of the well was demonstrated.
Note that for a system governed by Gaussian fluctuations (dS (2) t ), the corresponding SFPE takes on the form of the conservation law of the probability (∂ t p(x, t|x 0 ) + ∇J(x, t) = 0) with the current (J(x, t) = −xp(x, t|x 0 ) − σ 2 0 ∂ x p(x, t|x 0 )) and σ 0 , representing the intensity of the noise at the level of the Langevin equation [53]. At steady state, the divergence of the flow has to vanish [54]. For a trivial case (J ≡ 0), this implies that −xp s (x) − σ 2 0 ∂ x p s (x) = 0, and the driving force can be written as with a potential U(x) In other words, under Gaussian fluctuations, the detailed balance condition (J ≡ 0) is equivalent to the requirement that at equilibrium, the Boltzmann relationship between the weight of the state (equilibrium probability) and underlying potential surface exists (clearly, this condition does not need to be satisfied for general nonequilibrium systems [54]).
Condition J ≡ 0 implies a further thermodynamic relationship [25,26] between the probability density function and the relative or distance entropy that measures deviation from the stationary state (p s (x)): where Z ≡ exp(−U(x)/σ 2 0 )dx is the partition function, F ≡ −σ 2 ln Z stands for the free energy function, and the second term on the RHS denotes the Shannon entropy. If the strength of fluctuations (σ 2 ) can be related to the ambient temperature T, the left-hand-side of the expression above achieves a simple interpretation of the instantaneous free energy of the stochastic system [25,50] at hand. Notably, the relative entropy (or the instantaneous free energy of the system) defined above plays a role in the Lyapunov function (H-functional) [44] for a general class of stochastic Markovian systems described by Equations (5) and (6). Through an analogy to the Gibbs equilibrium state (cf. Equations (18) and (19)), it is thus tempting to interpret the functional Φ ≡ − ln p s (x) as a nonequilibrium pseudo-potential [25,50,52,55,56] from which the equivalents of "internal" and "free" energies of the system can be identified as in Equation (19) above. This approach yields Φ = F (x)m and is incompatible with the requirement that the reference (stationary) state should correspond to the Gibbs equilibrium with p s (x) fulfilling Equation (18). Indeed, the issue of thermalization or equilibration of Lévy flights can be only addressed by predefined confinement of the motion in properly "tailored" potentials [31,32,44,[57][58][59], or otherwise described by nonextensive Tsallis' thermodynamics [60]. There are also known "non-Langevin" scenarios [61,62] based on a master equation that result in Levy flights fulfilling a detailed balance with the stationary state of the Boltzmann-Gibbs type.
The scaling properties (Equation (7)) of the time dependent density p(x, t) result for the general case in dynamic (Shannon) entropy (cf. Equation (17)) taking on the form thus indicating the entropy production rate: In the case of the free Lévy flights described by the stability index α, the above rate formula results in d dt S[p(x, t)] ∝ (αt) −1 , i.e., the entropy production rate is positive, decreases for increasing α, and attains a minimal value for α = 2, corresponding to generic Gaussian fluctuations, typical for states close to equilibrium [63,64]. Similarly, for Lévy flights in the quadratic potential, the evaluation of Equation (21) yields d dt S[p(x, t)] = 1 e αt −1 which coincides with the rate of a free Lévy process for short times (t 1). These significant differences in the statistical properties of systems driven by non-Gaussian Lévy fluctuations, and, in particular, divergence of the second moment, imply the lack of a simple Einstein's fluctuation-dissipation relationship between fluctuation strength and the magnitude of dissipation. In a forthcoming section, we review the linear relaxation theory based on the identity derived by Hatano and Sasa [65] for dynamic Markov systems and present an extension of the generalized fluctuation-dissipation theorem to the nonequilibrium system subjected to thermal (Gaussian noise) and nonthermal (external Cauchy noise) fluctuations.

Linear Response and Fluctuation-Dissipation Relationship under Lévy Noise
Relation between externally induced and spontaneous fluctuations in systems close to equilibrium is described by the concept of the response function [66,67]. Linear response of a generic observable X(t) measured at time t to small perturbations f (t) and the correlation between unperturbed observable is characterized by the susceptibility measured in the reference, unperturbed state: where As a result, the celebrated Green-Kubo fluctuation-dissipation theorem Equation (22) links relaxation properties of the system to correlations of spontaneous fluctuations around the equilibrium state.
For systems described by the Hamiltonian H(x) and manipulated by weak forces f (t), the generic observable X(t) is identified as a "conjugate variable" in the Onsager's theory and can be interpreted as the fluctuation of the quantity Here, in line with Equations (18) and (19), the free energy is Onsager's theory and Green-Kubo relationships are grounded on the expansion around equilibrium states, which restricts their applicability to systems close to equilibrium. However, the generalization of the fluctuation-dissipation theorem to out of equilibrium systems has been established in a number of cases [65,68,69]. In their seminal paper from 2001, Hatano and Sasa [65] proposed an identity for Markov dynamic systems evolving between two steady states. The identity has been further explored [55,68] to derive the general form of the fluctuation-dissipation theorem for non-energy conserving dynamics. In what follows, we discuss an extension of this formalism to the generalized fluctuation-dissipation theorem for systems away from equilibrium and under the action of non-equilibrated environment (heat bath). We model interactions of the particle with the reservoir by assuming that there are two statistically independent sources of noise affecting the particle: a white Gaussian noise and a white Cauchy noise. The latter introduces large impulses, leading, effectively, to infinite moments of the position of the particle and also, infinite average energy. Accordingly, the statistical temperature of the system cannot be properly defined, i.e., the system is out of equilibrium even in its stationary state. Yet, despite those peculiarities, by properly identifying variables conjugated with external perturbations, the generalized form of the fluctuation-dissipation theorem can be re-established and the validity of linear response to perturbations can be verified.
We focus our attention on the system described the following Langevin equation where (S (2) t ) t 0 , and (S (1) t ) t 0 denote Wiener and symmetric Cauchy processes, respectively, with the scale parameters σ (2) = σ 0 and σ (1) = γ 0 standing for noise magnitudes and f (t) defining the external deterministic force field. Since the Langevin Equation (25) is linear, its solution depends linearly on two statistically independent noises and the corresponding probability distribution function (p(x, t)) of the dynamic variable (x(t)) attains the convoluted form of two Lévy PDFs with the stability indices α = 1 and α = 2. The corresponding characteristic function is expressed by a product p(k, t) = e ikµ(t)−σ 2 (t)|k| 2 −γ(t)|k| , (26) and fulfills [56,70,71] the generalized Smoluchowski-Fokker-Planck equation The Lévy jump statistics entering dynamics through the noise term S t lead to an inherently nonequilibrium situation (note that due to the presence of the Cauchy noise, the expected value (E[X t ]) does not exist for γ 0 = 0), and the overall probability distribution of the process deviates from the Gibbs-Boltzmann form, exhibiting bulk Gaussian part and asymptotic heavy tails. Inserting Equation (26) where we assume µ(0) = 0 at time t = 0. For a constant force ( f (t) = f ), the long time limit of these equations gives the stationary parameters lim t→∞ µ(t) = f /a, lim t→∞ γ(t) = γ 0 /a, lim t→∞ σ 2 (t) = σ 2 0 /(2a) characterizing a non-equilibrium steady state p s (x, f ) = lim t→∞ p X t (x) of the system. We further assume that the system is initially prepared in this stationary state. The question we are going to address is whether the generalized form of the fluctuation-dissipation theorem can be used here to test the linear response of the system to external drivings. This problem is addressed below in the framework of stochastic thermodynamics [25,26,50,[52][53][54].
In line with Onsager's theory, we first introduce a non-equilibrium pseudo-potential [55,56] Φ(x, f ) = − ln p s (x, f ) and use it to define a variable X GC that is conjugated to the perturbation f : where w(x) = e −x 2 erfc(−ix) is the complex error function [1,29], and the second equality in Equation (29) follows from the linearity of Equation (25). Thus, the conjugate variable X GC defines a new stochastic process related to the original {X t , t 0} by a nontrivial transformation: In Figure 3, the functional dependence (X GC (x)) is displayed for σ 0 = 1, a = 1 and several different values of the Cauchy noise strength (γ 0 ). In contrast to the pure Gaussian case (cf. Equations (26) and (29) with γ(t) ≡ 0), when the conjugate variable reads X G = −2x/σ 2 0 , the variable X GC is a nonmonotonic function of x. This is related to the heavy tails of the noise distribution and not to the mixing of two types of noise, as can be seen from the form of the conjugate variable in the case of σ 0 = 0, a = 1 [55,72] X C (x) = − 2x The fluctuation-dissipation theorem (or its generalized form) is expected to link the autocorrelation function of the conjugate variable X GC (t)X GC (0) to a generalized susceptibility where the averaging is performed over the stationary (unperturbed f = 0) state. The expected response of the conjugate variable to the perturbation f (t) can be then either determined exactly from where p X t (x, t) is calculated including f (t) and under stationary initial conditions or otherwise, estimated from the linear response function The susceptibility and the corresponding linear response (LR) of the conjugate variable (35) can be easily derived analytically [55] if either γ(t) = 0 or σ(t) = 0. Remarkably, although the linear system with pure Cauchy noise (σ(t) = 0, γ(t) = 0) has a nonlinear conjugate variable (32), its susceptibility is proportional to e −t , as in the case of pure Gaussian noise ((σ(t) = 0, γ(t) = 0)) with a linear conjugate variable. We confirmed that prediction with simulations (see Figure 4). In other cases χ(t) can be obtained by averaging over a large number of simulations of the unperturbed stochastic differential equation (SDE). Importantly, the simulations show that the combination of Gaussian and Cauchy noise (25) preserves the exponential time dependence of the susceptibility. This seems to be a general property of the susceptibility of linear systems driven by additive white noise and can be investigated by spectral decomposition of the fractional Fokker-Planck equation [73]. Similarly, exact responses according to Equation (34) can be obtained by means of stochastic simulations of the perturbed SDE. We performed tests of the generalized fluctuation-dissipation theorem by analyzing the response of the system to the sum of a small periodic and a linearly increasing force ( f (t)). The integrals in Equation (35) were evaluated numerically and compared with the exact results of Equation (34) (see Figure 5). A comparison between the exact Equation (34) response and the one evaluated by means of the generalized fluctuation-dissipation theorem Equation (35) clearly indicates that the theory holds for small forcing ( f (t)). This observation seems, at first sight, rather counterintuitive: Unlike standard Gaussian uncorrelated fluctuations, Lévy noises introduce bursty-like large jumps to the overall displacement (dX t ) causing divergence of the mean-squared displacement [45,52]. However, suitably chosen conjugate variable provided by the generalized fluctuation dissipation theorem (FDT) allow the response of the system to be analyzed using the concept of susceptibility (χ(t)). It should be noted that although the exact response can be derived using Equation (34) directly, the formulas are cumbersome and can be solved numerically in the simplest cases only. Therefore, the generalized FDT is a practical tool [55]. This is, however, not without a price. The conjugate variable used in the framework of FDT represents the change in the probability distribution of the system under the action of external disturbances. Close to equilibrium, such change is related to the energy absorbed by the system from external forcing. For non-equilibrium systems driven by Lévy fluctuations, this interpretation may be invalid. First, the conjugate variable is not just a displacement related to the external forcing (cf. Figure 3). Next, its proper identification requires knowledge of the stationary state which might be a demanding task for nonlinear systems. Last, but not least, divergence of the moments of dynamic variable x plague the system with infinite energy [72,74,75] (e.g., the potential energy of the system in the harmonic potential studied above).
We have checked whether the linear response function (35) can be used to predict responses of different nonlinear transformations (observables) of the dynamic variable, see Figure 4. To our surprise, some of the observables responded closer to the prediction of the linear response function than the conjugate variable. This suggests that the conjugate variable is not necessarily the best observable to monitor linear response and more detailed research is required to identify the conditions under which the linear response works.

Summary and Conclusions
We analyzed the peculiarities of jump-like Lévy noises with heavy tail statistics and demonstrated their distinctively different features from their Gaussian counterparts. For probability density functions generated by Langevin equation with non-Gaussian noise (Poisson shot noise or Lévy noises), violations of conventional steady states and transient fluctuation relations have been demonstrated [56,[74][75][76], and the athermal character of non-Gaussian random forces has been investigated [77]. Here, we showed that despite the strong nonequilibrium character of Lévy flights, the generalized version of Onsager's theory and fluctuation-dissipation theorem can be adapted to capture the dynamic responses of the system subjected to this type of noise. The linear response to external perturbation can be correctly predicted if the system is properly analyzed, i.e., one should transform the observables into the corresponding conjugate variables.
We predict that studying the fluctuations and linear responses of nonequilibrium systems with heavy-tailed noise will be useful for understanding active transport phenomena in biology. Lévy flights-like patterns have been observed in movement data of diverse organisms, including albatrosses [78], humans [79], honeybees [80], swarms of bacteria [81], and many others [82]. This is being commonly explained by the fact that this kind of behavior can constitute an advantageous random foraging strategy [20,21,24,83]. Although many models have been proposed [22,82,84], most of them apply to very specific cases and are yet to be tested experimentally. Therefore, the generative mechanisms behind Lévy patterns are still obscure. The linear response theory with nonlinear conjugate variables may be a fruitful way of looking at these systems, helping to identify the source of Lévy noise by applying perturbative interventions and inferring the underlying mechanisms.
The framework described in this paper is limited to white noise. Generalizations of the linear response theory to non-Markovian processes, such as subdiffusive continuous time random walks [47,48] and superdiffusive Lévy walks [33] are needed. Another interesting avenue for future research are systems with stochastic resetting [20,21,85,86], in which the particle is transported back to the initial position at random times. The stochastic thermodynamics of resetting has recently been analyzed [87]. However, it seems likely that the linear response cannot be directly applied in systems subject to stochastic resetting. For example, the nonequilibrium steady state of one-dimensional diffusion with stochastic resetting is given by the Laplace distribution, which leads to a noninformative conjugate variable that is proportional to the sign of the position.

Conflicts of Interest:
The authors declare no conflict of interest. The founding sponsors had no role in the design of the study; in the collection, analyses, or interpretation of data; in the writing of the manuscript, and in the decision to publish the results.