key: cord-0491770-gti4y0iv authors: Gao, Fei; Liu, Shuaiqiang; Oosterlee, Cornelis W.; Temme, Nico M. title: Solution of integrals with fractional Brownian motion for different Hurst indices date: 2022-03-04 journal: nan DOI: nan sha: d16123b7ee6c7240a61bc7d338f0d8abbb355122 doc_id: 491770 cord_uid: gti4y0iv In this paper, we will evaluate integrals that define the conditional expectation, variance and characteristic function of stochastic processes with respect to fractional Brownian motion (fBm) for all relevant Hurst indices, i.e. $H in (0,1)$. The fractional Ornstein-Uhlenbeck (fOU) process, for example, gives rise to highly nontrivial integration formulas that need careful analysis when considering the whole range of Hurst indices. We will show that the classical technique of analytic continuation, from complex analysis, provides a way of extending the domain of validity of an integral, from $Hin(1/2,1)$, to the larger domain, $Hin(0,1)$. Numerical experiments for different Hurst indices confirm the robustness and efficiency of the integral formulations presented here. Moreover, we provide accurate and highly efficient financial option pricing results for processes that are related to the fOU process, with the help of Fourier cosine expansions. In this paper, we will evaluate the integrals, that define the conditional expectation and variance of stochastic processes that are based on fractional Brownian motion (fBm). For this, we first define fBm, to which a Hurst index is associated, as a generalization of standard Brownian motion (sBm). Let (Ω, Σ, P) be a complete probability space, B H = {B H t , t ≥ 0} be an fBm with Hurst index H ∈ (0, 1) defined on this probability space, and F t be its natural filtration. Fractional Brownian motion, B H , is then uniquely characterized by the following properties [ 3) E[B H t ] 2 = t 2H , for t ≥ 0; 4) B H t has a Gaussian distribution, for t > 0. From the first three properties, it follows that the covariance function is given by The special case of Equation (1.1) is sBm (H = 1/2), with covariance function Cov (B 1 2 t , B 1 2 s ) = min (t, s). Particularly for fBm, and different from sBm, the increments are not independent and fBm is neither a Markov process nor a semi-martingale, for H = 1/2. As a consequence, fBm enables modeling of stochastic processes that are self-similar or exhibit long-range dependence in time. Properties of fBm have been presented, for example, in [25] and [4] . Fractional Brownian motion forms the basis for a broad class of stochastic processes, like the fractional Ornstein-Uhlenbeck (fOU), the fractional Cox-Ingersoll-Ross (fCIR) and the geometric fOU (GfOU) processes. Such processes have played an increasingly important role in diverse application fields, like hydrology [22, 20] , telecommunication [24, 28] , weather forecast [34, 31] or epidemic disease (such as the Corona virus) modeling [1] . FBm has also been used in finance for the modeling of commodities, weather (now considered as an asset) or currencies in financial products, that are based on these assets, in [30, 6, 2, 3, 37] , for bond markets in [15] and recently for rough volatility in [16, 19] . In these application areas, we require the computation of conditional expectations, which relate, next to the average behaviour under uncertainty, also to the valuation of options on these non-tradable assets in finance. For Markov processes, the accurate computation of conditional expectations and variances has been given ample attention. In the diverse literature, we encounter, for such processes, approaches based on numerical integration, Fourier or Laplace inversion, Monte Carlo simulation and recently neural network approximations. For non-Markov processes, to which fBm belongs, the literature is more scarce but some significant books and papers have been presented in the last twenty years, like [25, 4, 2, 15, 14, 8] . The conditional characteristic function for the fOU process has been derived in [15] and was evaluated in [15, 35] , where involved integration formulas had to be calculated numerically. These integration formulas are nontrivial when considering all Hurst indices H ∈ (0, 1), due to occurring singularities in the integrand. Particularly, for the range H ∈ (0, 1/2), we encounter approximation formulas in the literature that are not easily understood from a numerical point-of-view. In this paper, we therefore reconsider the expectation and variance related integration formulas for all relevant Hurst indices, and aim to approximate them numerically on the basis of accurate (classical) approximations and numerical techniques. The technique of analytic continuation [33, 32] is used in this paper to handle the integration for certain Hurst index values. Analytic continuation provides a way of extending the domain over which a complex function is defined. Here, it is used to extend the domain of validity of an integral which is defined on an interval with respect to the Hurst index, from H ∈ (1/2, 1), to the larger domain, H ∈ (0, 1). Some examples and pointers for details on analytic continuation are presented in A. Here, we calculate the conditional expectation, variance and the conditional characteristic function of fBm and fBm driven stochastic differential equations (SDEs). The accuracy of our numerical results is verified by results obtained by Monte Carlo simulation. As an application in finance, we discuss the option pricing problem where the asset price is modeled as a fractional stochastic process. Clearly, fBm is not a semi-martingale, but some assets, like those based on the weather or volatility, are non-tradable products so that the well-known financial derivative pricing theory, based on the no-arbitrage and market completeness, cannot be used anyway. Valuation of such a derivative is then based on the expected discounted value approach under the real world probability measure. With the characteristic function for the fOU process available, we can compute the conditional expectation, by means of the COS method, a Fourier-based numerical integration technique, see [10] , [11] . The COS method is based on an approximation of the appearing (conditional) probability density function by a Fourier cosine series expansion. The characteristic function of some of the fractional processes considered here is either not available or not easy to derive, however, by means of a transformation of variables, the COS method can still be employed. The remainder of this paper is organized as follows. Section 2 briefly introduces fBm and the related stochastic processes. In Section 3, we analyze the conditional variance of the fOU process in detail and present a convenient numerical method to compute the appearing integrals. Option pricing with the COS method is described in Section 4. Section 5 presents numerical results including conditional distributions and option pricing. Section 6 gives our conclusions of the paper. In this section, we discuss some well-known properties of fBm and also several stochastic processes driven by fBm. Fractional Brownian motion B H has the following three important properties. where 0 ≤ s ≤ t. The function Cov (Y t , Y s ) is presented, for various H, in Figure 1 . Cov (Y t , Y s ) = 0 for any t and s, when H = 1/2, since sBm has independent increments, whereas increments of fBm are not independent. The increments are positively correlated when H > 1/2 and negatively correlated when H < 1/2. Moreover, when H > 1/2, Cov (Y t , Y s ) decays very slowly with increasing time increment (t − s), and therefore fBm exhibits long-range dependence when H > 1/2. The fractional Ornstein-Uhlenbeck (fOU) process, X = {X t , t ≥ 0}, is based on fBm and is defined by the following SDE: with an fBm B H and given initial value X 0 . The parameter µ represents the long-term mean value of the process, λ the speed of mean reversion, and σ the volatility. It has a unique path-wise solution as For H = 1/2, using Itô's isometry property, the analytical solution of the Ornstein-Uhlenbeck process driven by standard Brownian motion reads: where W is a sample from the standard normal distribution N (0, 1). Using κ := H − 1/2, and based on the derivation of conditional characteristic function, it is proved in [15] that under the filtration F s (s ≥ 0), X t (t ≥ s) is normally distributed, with where c(r) = σe −λ(t−r) ; and, for v ∈ (0, s), while for v ∈ {0, s}, Ψ c (s, t, v) = 0. For κ = 0, the norm of a function f is defined by the fractional Riemann-Liouville integral where Γ(·) is the gamma function. For κ = 0, f 2 0,T is set as the L 2 norm, i.e., The accurate numerical evaluation of Equation (2.5) and (2.6) is a focus of the present paper. These equations also form the basis for the conditional characteristic function. Remark 1 (Conditional distribution of fBm). The fBm can be considered as a special case of the fOU process by taking λ = 0 and σ = 1 in Equation (2.2). Then, the conditional distribution of the fBm is also a normal distribution, with where, for v ∈ (0, s), There are some important processes Z t that are related to the fOU process, i.e., they may be defined based on the fOU process, as follows, where X t is the fOU process and the function g(·) is invertible. With z = g(x), the cumulative distribution functions (CDFs) of Z and X are related, as follows, By differentiation, we have the relation of their probability density functions (PDFs) as Since the PDF of the fOU process is available, the PDF of a stochastic process related to the fOU process is also available, by means of Equation (2.16). We are interested in the following three processes. 1) The geometric fractional Ornstein-Uhlenbeck (GfOU) process is defined as the process whose logarithm is the fOU process, This process is found in finance, where we tend to model the logarithm of asset prices. 2) The fractional Cox-Ingersoll-Ross (fCIR) process is an extension of the classical Cox-Ingersoll-Ross (CIR) process, which is used for short-term interest rate and stochastic volatility modeling [7, 18] . Compared to the CIR process, the fCIR process has the advantage to model the "memory phenomenon" in financial data [21] , like in rough volatility. It is shown in [14] that, for the case κ > 0, an fCIR process with zero mean, is the square of the fOU process before hitting zero for the first time, i.e., 3) Equation (2.18) is an example of a polynomial process. Polynomial processes form a class of processes that are obtained by a polynomial map of the underlying process [12] . This class plays an increasingly important role in finance as it provides a tractable relationship between underlying factors and resulting prices. For example, a second-order polynomial process is employed for the valuation of electricity storage contracts in [36, 5] . Using the same increasing polynomial map construction method as in [36] , we here consider a polynomial process with a third-order polynomial map, i.e., In the previous section, we reviewed the conditional distribution of the fOU process, which is a conditional normal distribution with expectation as in Equation (2.5) and variance in Equation (2.6) . In this section, we calculate the relevant integrals by which the conditional distributions of the fOU process and several related processes can be derived for all relevant Hurst indices. With a discretization in the time interval [0, s], i.e., 0 = s 0 ≤ s 1 ≤ · · · ≤ s n = s, and the corresponding discrete realizations of B H , according to the introduction of integrals with respect to fBm in [29] , the conditional expectation can be calculated by The variance, in Equation (2.6), is expressed as the difference of two terms, each being the square of a norm that is not easily calculated, since the expression involves a nontrivial integration with multiple variables and singularities in its integral function. We first simplify the expression of the conditional variance in Section 3.1 and obtain a solution for all H ∈ (0, 1) in Section 3.2, which is summarized in Algorithm 1. Noting that both terms in Equation (2.6) contain indicator functions, we rewrite them as Using a new integration variable, u ∈ (0, ∞), with , (3.5) and interchanging the order of integration, yield Here, we recognise the Mellin transform (see [27, Table 1 This gives us, We have derived this relation for 0 < κ < 1, as assumed in Equation (3.7), however, in later formulas we need this expression, or a similar result, also for κ ≤ 0. For this, we use the classical technique of analytic continuation. For details on this topic, we refer to A. Substituting Equation (3.8) into Equation (3.3) and comparing with Equation κ,T and the expression of the conditional variance in Equation (2.6) can be simplified, from two terms to one term only: (3.9) Knowing that Γ(κ + 1) = κΓ(κ) and πκ the front factor A κ of Equation (3.9) can also be written as To further illustrate how to calculate the variance, we write Equation (3.9) in the form, and c(r) = σe −λ(t−r) . (3.14) A set of transformations is now proposed for the calculation of H(s, t, κ) in Equation (3.12): 1) First, we use the transformation x = z 1−2κ , and obtain 3) The final transformation is given by y = 1 2 erfc(−w), and we get where erfc(·) is the complementary error function with the following properties The final integral in Equation (3.17) can be approximated efficiently by the Trapezoidal Rule, where m > 0 is the step size. Because of the fast convergence of the series, the range of summation can be cut off to a smaller interval t ∈ [−a, a], a > 0: For details about the Trapezoidal Rule, see [17, §5.4] . In our case, f (w) = h 2 (κ, z), with w a function of z. Here we will derive convenient computable expressions for the h(κ, z) term from Equation (3.13), which is the only involved part remaining. We expand c(r) from Equation (3.14) in powers of r and get an expansion with fast convergence, i.e., Then, we obtain Initially, we consider κ > 0, because for r = z we need convergence of the integral. Integrating by parts gives us, In this representation, we can accept κ > −1, so the integration by parts gives us again an analytic continuation 2 of R n (κ, z) with respect to κ. By using the fact that (r − z) κ = (r − z) κ−1 (r − z), a simple recursion formula follows The functions R n can be expressed in terms of Gauss hypergeometric functions. For details about these functions, see [26] . In particular, we have the integral representation and power series, In this integral, we take u = (1 − z/t)v, which gives us, by invoking the integral in (3.26), (3.28) The representation in the second line of Equation (3.28) is valid for κ > −1, and it gives us again the analytic continuation to the domain of h(κ, z) initially defined for κ > 0. Remark 2. From Equation (3.23), we have observed that we do not need to handle a singularity at the endpoint z = s when s > 0. However, this will be needed when s = 0, in which case the variable of integration, z in Equation (3.12), will reach the value z = 0, but the 2 F 1 -function in Equation (3.28) is not defined at z = 0. For this case, a different form of R n (κ, z) will be required, which is given by R n (κ, z) = Γ(κ + 1)Γ(n + 1 + κ) 2 cos(πκ)Γ(n + 1 + 2κ) The way this form is derived is presented in B. The 2 F 1 -function here is not defined for z → t when κ > 0. In that case, Equation (3.28) should be used. In this section, we will consider the computation of conditional expectations based on the calculation of the expectation and variance of the fOU process described in the previous sections. Usually, option valuation starts from the risk-neutral valuation formula, which is a discounted expectation under the risk-neutral asset price process. Risk neutrality is related to complete markets, no arbitrage and the asset price process being a martingale, i.e., the expected fair value of the discounted price process at a future time point is the present asset price. Processes based on fractional Brownian motion are not (semi-)martingales, so the risk-neutral pricing formula does not apply in this case. Financial option valuation related to weather [2, 6] , or the asset price's volatility, is based on a non-tradable asset for which the real world probability measure is often employed. Under the assumption that such a measure has been defined, option pricing can be cast in an integral form: where ∆t := T − t, V denotes the option value, x the state variable at time t and y at time T , f (y|T, F t ) the transition probability density, and r the interest rate. Here, we focus on the European option pricing problem when the asset prices are modeled as the fOU related processes, as described in Section 2.3. For these processes, the transition probability, f (y|T, F t ) in Equation (4.1), is the PDF of a normal distribution with expectation µ y = E[X T |F t ] and variance σ 2 y = Var[X T |F t ], as in Equation (2.5) and Equation (2.6), respectively. It is therefore possible to obtain a closed-form solution of Equation (4.1) for some specific processes. For example, if the asset price follows a GfOU process as Equation (2.17), and the payoff reads where K is the strike price, the closed-form solution of Equation (4.1) for a call option can be derived as and for a put option, Here, Φ(·) denotes the CDF of the standard normal distribution. For more general processes, e.g., polynomial processes, a closed-form option value expression is difficult or even impossible to obtain. It is however possible to employ the COS method in such cases. The COS method, [10, 11] , can be applied to compute expectations of (underlying) processes for which the characteristic function is available. For Markov processes with affine dynamics, the characteristic function can be easily derived by solving the Ricatti differential equations [9] . Affinity is not invariant under polynomial transformations [13] and therefore the characteristic function of a stochastic model which is based on a polynomial process, Z t = g(X t ), generally does not exist. However, by using a transformation z = g(x) in the definition of the COS method, the state variables can be conveniently chosen so that the characteristic function of the underlying process X t , can be employed within the pricing formula, with, again, µ y = E[X T |F t ] and σ 2 y = Var[X T |F t ], as in Equation (2.5) and Equation (2.6) . The conditional probability density function in Equation (4.1) is approximated by a truncated Fourier cosine expansion, which uses the characteristic function to recover the Fourier coefficients, as follows [10] : where [b, d] is the truncated integration interval, L is the number of terms in the Fourier-cosine expansion and implies that the first term of the summation is multiplied by 1 2 . By replacing the conditional density function f (y|T, F t ) in Equation (4.1) by its Fourier cosine expansion approximation Equation (4.6) and interchanging integration and summation, with the help of Fubini's theorem, the COS formula is obtained: where coefficients V l are defined by: where, importantly, the state variables x and y can be any function of respectively the asset prices Z t and Z T , e.g., in our case, x = g −1 (Z t ) = X t and y = g −1 (Z T ) = X T , with g(·), invertible, as in Equation (2.14) . This means that we can simply work with the characteristic function of the fOU process X t , but price options under the more involved process Z t . A closed-form solution of the coefficients V l is available for various payoff functions and several choices of the state variables x and y [10] . We present several numerical results for the methods presented, starting with the calculation of the conditional variances of the fBm and fOU processes, followed by the computation of probability density functions of the related processes, i.e., the GfOU, fCIR and polynomial processes, and some option pricing results using the COS method. These results show the accuracy and efficiency of the proposed computational techniques. First of all, we focus on the accurate representation of the conditional variance by means of the integrals that were derived in Section 3. As explained in the previous sections and also shown in Algorithm 1, there are three hyperparameters in our numerical method to compute the integrals in the approximation of the conditional variance for the fOU process, i.e., the step size m, summation range a for the Trapezoidal Rule in Equation (3.21) , and the number of expansion terms N in Equation (3.23) . To determine suitable values for these parameters, we perform the following experiment. For two values of H, and the value of N sufficiently large (i.e. N = 20), the computed results with varying m and a are presented in Figure 2a . The figure shows that the summation range value a has a bigger impact than the step size m. The numerical approximation has converged when a is larger than 3 and m is smaller than 0.6. We therefore set a = 5 and m = 0.5 in the subsequent numerical experiments. Next, we consider the effect of the number of expansion terms N on the computed results, also for the two values of H, and plot it in Figure 2b . It shows that our numerical scheme converges rapidly as N increases, and setting N to be larger than 10 already provides a reliable computational result. We may define a stopping criterion by checking when, e.g., |c n R n /h| is less than 10 −8 , then stop the summation. The results presented are representative for a large range of model parameters. To analyze the accuracy of our numerical technique, we compare the conditional (s = 3) and unconditional (s = 0) variances obtained by our method with Monte Carlo simulation. FBm and fOU, with σ = 0.3, λ = 0.5, are considered, with different Hurst indices H. For the simulation of fBm, we use a Python package named "fbm" 3 . Two sets of parameters are considered, one with a step size of 0.01 and the number of Monte Carlo paths of 10 4 , and another with the step size halved and the number of paths multiplied by 10. We report the relative errors in Table 1 , which are defined as where the standard deviation of the Monte Carlo results isσ 1 and the computed standard deviation by numerical integration isσ 2 . As shown in Table 1 , our numerical technique provides highly accurate results very close to those with the Monte Carlo simulations, with all average relative errors less than 1%. By using smaller time step and simulating more paths within the Monte Carlo simulation, the differences are even smaller. The CPU times 4 are also compared in Table 2 . Clearly, a numerical integration technique requires much less time than a Monte Carlo sampling method, because Monte Carlo simulation requires the generation of many paths with a small time step to obtain reliable results. We further present the computed variances for different H, time point s and with increasing time intervals (i.e., t − s) in Figure 3 (fBm) and Figure 4 (the fOU process with σ = 0.3, λ = 0.5). From these figures, we can conclude that our computational technique converges very well for all relevant H indices, in a robust way. The curves for the different values of H are very accurate, as compared to the literature [15] . As the time interval t − s increases, the variance increases. With a sufficiently large time interval t − s, large values of the Hurst parameter H results in higher variances. This makes sense as the Hurst index indicates the roughness of a process, with a higher value leading to smoother paths. Comparing the presented variances for s = 0 and s = 3, different results are observed, for the same time interval t − s, which confirms that the variance is not a function t − s, but of s and t separately. As described in Section 2.3, due to the relation between the PDFs of the fOU process and the related processes in Equation (2.16), we will compute the PDFs of the related processes by means of the PDF of the fOU process, which is normally distributed with the expectation in Equation (3.1) and the variance computed as in Algorithm 1, and the transformation g(·) in Equation (2.14). Figure 5 plots the PDFs of the GfOU, fCIR and polynomial process, as in Equation (2.17), Equation (2.18) and Equation (2.19) , respectively. We calculate the PDFs on the condition of either initial value (s = 0, t = 3) or realized historical paths (s = 3, t = 6), and compare them with the results by Monte Carlo simulation. It is shown that the calculated PDFs resemble the histograms of Monte Carlo simulation very well, which is an indication of the accuracy of our technique. Moreover, by comparing the first and third columns of Figure 5 , we can see the impact of the memory of the processes (the second column of Figure 5 ) on the PDFs. The COS method enables us to price options if we have the characteristic function available for a function of the asset price process, like X t in Z t = g(X t ). When the asset prices are modeled as one of the fOU related processes, the COS method will be based on the characteristic function of the fOU process X t . To test the performance of the COS method, we price a European call option under the GfOU process. This is just a reference example, since a closed-form solution is available in this case, as in Equation (4.3) . The errors of the COS method are presented in Table 3 with the strike price K = 10, the interval of integration [b, d] = [µ y − 10 σ 2 y , µ y + 10 σ 2 y ] and a varying number of cosine terms L, see Equations (4.6) and (4.7). It is shown that the COS method converges rapidly, as with only L = 16 Fourier-cosine expansion terms, the errors are of the order 10 −8 . We further present the prices of European call and put options under the GfOU process in Figure 6 and the polynomial process in Figure 7 , respectively, for different Hurst indices H and strike prices K. Obviously, the Hurst index H has a significant impact on the option price, which indicates the importance of modeling the asset price as a stochastic process with a suitable H. 2.8e-08 2.5e-08 2.1e-08 1.6e-08 9.1e-09 1.0e-09 8.6e-09 2.0e-08 3.2e-08 In this paper, we have shown that, on the basis of several transformations of variables and the analytic continuation technique, integrals that define the conditional expectation, variance and characteristic function of stochastic processes, with respect to fractional Brownian motion (fBm), can be reformulated so that they exist for all relevant Hurst indices. Numerical experiments for the conditional variance of the fOU process, and to recover the probability density function of certain stochastic processes that can be derived from the fOU process, confirm the robustness and efficiency of the integral formulations and of the numerical technique, regarding different values of the Hurst index. Moreover, by means of the COS method, we have presented accurate and highly efficient option pricing results for processes connected to the fOU process, on the basis of the conditional characteristic function of the fOU process. The methodology presented here can be used as an alternative to Monte Carlo simulation when appropriate, or as a validation technique for implementations of other numerical techniques. We first take f (t) = 1 for t ≥ 0. In that case, F (x) = 1, and we see that the combination of the front factor and the integral, considered as one quantity, becomes a function which is defined for all real and complex values of x. A second step is writing f (t) = f (0) + (f (t) − f (0)), or introducing the functions The function G(x) is an analytic function for {x} > −1, and F (x) = G(x) for {x} > 0, and we conclude that F (x) which was initially defined for {x} > 0, has G(x) as its analytic continuation for {x} > −1, and F (x) = G(x) for {x} > −1. For further details on analytic continuation, we refer to the references [32, Section 1.10(ii)] and [33, Chapter IV] . A second method is based on integration by parts. We have, using xΓ(x) = Γ(x + 1), for (x) > 0, The new integral is defined for {x} > −1, and, again, we have obtained the analytic continuation of F (x) from {x} > 0 to {x} > −1. We can repeat these two methods to obtain for F (x) the analytic continuation for {x} > −2, {x} > −3, and so on. In (3.12), we see that we do not need to handle a singularity at the endpoint z = s when s > 0. However, this will be needed when s = 0, in which case the variable of integration, z in (3.13), will take the value z = 0. To deal with this, first of all, we write the representation of R n (κ, z) in (3.28) in a different form, because the 2 F 1 -function is not defined at z = 0 (although the front factor z 2κ+n will control the product of the two quantities). and w 4 (z) is given by This gives us, for (3.28), R n (κ, z) = Γ(κ + 1)Γ(n + 1 + κ) 2 cos(πκ)Γ(n + 1 + 2κ) This 2 F 1 -function is not defined for z → t when κ > 0. In that case Equation (3.28) should be used. Solutions of a disease model with fractional white noise On arbitrage-free pricing of weather derivatives based on fractional Brownian motion Modeling and pricing in financial markets for weather derivatives Stochastic calculus for fractional Brownian motion and applications Valuation of electricity storage contracts using the COS method Dynamical pricing of weather derivatives A theory of the term structure of interest rates Simulation of fractional Brownian motion Transform analysis and asset pricing for affine jump-diffusions A novel pricing method for European options based on Fourier-cosine series expansions Pricing early-exercise and discrete barrier options by Fourier-cosine series expansions Polynomial diffusions and applications in finance Polynomial jump-diffusion models Fractional Lévy-driven Ornstein-Uhlenbeck processes and stochastic differential equations Conditional distributions of processes related to fractional Brownian motion Volatility is rough. Quantitative finance Numerical methods for special functions A closed-form solution for options with stochastic volatility with applications to bond and currency options. The review of financial studies Rough volatility: evidence from option prices An efficient, three-dimensional, anisotropic, fractional Brownian motion and truncated fractional Lévy motion simulation algorithm based on successive random additions Fractional Cox-Ingersoll-Ross process with non-zero mean Fractional Brownian motion and fractional Gaussian noise in subsurface hydrology: A review, presentation of fundamental properties, and extensions A storage model with self-similar input. Queueing systems On the use of fractional Brownian motion in the theory of connectionless networks Fractional Brownian motion: stochastic calculus and applications Chapter 15, Hypergeometric function Nist handbook of mathematical functions, us department of commerce national, national institute of standards and technology Simulation of telecommunication traffic using statistical models of fractional Brownian motion Are classes of deterministic integrands for fractional Brownian motion on an interval complete? Bernoulli Long-range dependence in the volatility of commodity futures prices: Wavelet-based evidence A new approach for time series forecasting: Bayesian enhanced by fractional Brownian motion with application to rainfall series Chapter 1, Algebraic and analytic methods The theory of functions Long-range correlations in the extratropical atmospheric circulation: Origins and implications Modeling and forecasting realized volatility with the fractional Ornstein-Uhlenbeck process Polynomial processes for power prices Pricing currency options in a fractional Brownian motion with jumps In complex analysis, analytic continuation is a classical technique to extend the domain of definition of a given analytic function. Analytic continuation often succeeds in defining further values of a function, for example, in a new region where an infinite series representation, in terms of which it is initially defined, becomes divergent.As a standard example to explain analytic continuation for a class of integrals, we consider the function F (x) defined bywhere the function f (t) is analytic in [0, ∞), and in a complex domain around this interval. F (x) is an analytic function in the domain {x} > 0. Next, we assume that f (0) = 0 and, in that case, the integral does not converge if x = 0. The reciprocal gamma function in front of the integral will allow us to give an interpretation to F (x) for x ≤ 0.