On moments of the integrated exponential Brownian motion

We present new exact expressions for a class of moments for the geometric Brownian motion, in terms of determinants, obtained using a recurrence relation and combinatorial arguments for the case of a Ito's Wiener process. We then apply the obtained exact formulas to computing averages of the solution of the logistic stochastic differential equation via a series expansion, and compare the results to the solution obtained via Monte Carlo.


I. INTRODUCTION
The geometric Brownian motion is the stochastic process described by the differential equation where W t is a Wiener process and µ, σ are constants describing the drift and the variance of the noise, respectively. The solution can be written as Geometric Brownian motion is used for modelling many phenomena in a variety of contexts [1].
A prominent role is played in financial applications, where the distribution of returns can be approximated by a log-normal distribution [1,2], at least in specific regimes.
For the computation of certain properties, it is necessary to compute the integral of f (W, t) over a time interval The evaluation of this functional is also involved in the solution of the geometric Brownian motion with logistic corrections. In general, averages of the form are quite common. The evaluation of averages of powers of the integrated exponential Brownian motion, then, is instrumental for the computation of these observables.
Detailed studies of this functional and of its powers [3,4] are already available in the literature.
In this paper, we will derive exact formulas for the evaluation of these integrals, under the assumption of the Ito formulation for the Wiener process. Similar results have been given in [2,5].
Motivated by obtaining exact formulas for Asian options, in [5] Yor obtained an exact formula in terms of polynomials for the following moments: Using Girsanov's theorem [6], one can derive a series of identities, in which the last is Bougerol's formula where and In this work, we take a different route with the use of combinatorics. We prove a recurrence relation for the integrals involved at the k-th order in terms of integrals at the (k − 1)-th order, and after resummation, we get an identity in terms of a determinant.

II. CALCULATION OF MOMENTS
The central quantity of interest in the present paper is given by the average over the Wiener process W s : If we expand Eq. (9), we obtain By using this property, we can now prove the following fact: For the average over the Wiener process W s of Ito type, the following formula holds true: Proof. By a direct application of Eq. (11) Due to symmetry of integrand, we can order the integration variables ast i <t i+1 , obtaining whence: After rearranging carefully the terms, we arrive at the final result: Let us now expand further on Eq. (16). It is convenient to first perform the rescaling t i = tu i .
Then, by defining λ j = t(µ + (k − j)σ 2 ), we obtain Therefore, the computation of (9) reduces to the computation of a very similar formula had been obtained in [5]; the main difference between the treatment made there and our relies on what follows. An important observation will enable us to evaluate these integrals exactly by means of combinatorics: it is possible in fact to prove the following result, which establishes a recursion relation among the s k .
Proof. Let us write Integrating by parts, where Finally, the identity f (0) = 0 gives the desired result.
Eq. (19) suggests the evaluation of the averages by means of combinatorial considerations.
Indeed, the evaluation of the integral can proceed graphically, for any fixed order of the moment k, as in Fig. 1. Starting from the top and using the properties of the recurrence relation, at each order k, we organize the recurrence on a binary tree. In Fig. 1, each left branch will pull out a factor e λ /λ, where λ is the λ obtained from the previous order. One has to consider the fact that however, in each right branch one pulls out a factor −1/λ, and sums the two factors of λ's. To obtain the final formula, once the empty set has been reached, one multiplies the final term by all the factors in the branch. We now give exact formulas for all the terms obtained from the recurrence. To lighten the notation, we define µ k,j = λ j + · · · + λ k and s k = s k (λ 1 , . . . , λ k ).
By iterating the second term of the recurrence, we obtain Therefore, Proof. We proceed the proof by induction on k. We have f k = e 0 0 , which holds for k = 0. We assume that the claim holds for k and let us prove it for k + 1. By the recurrence relation, we have

By induction hypothesis, we have
which completes the induction step.
and e k 0 = , then Eq. (23) can be written as Applying Lemma 3, we obtain , which leads to the following result.

The above theorem gives
.
This expression can be elucidated by means of an example. Consider the case n = 3 and σ = 0.
In this case, we have an exact formula for λ j = µt: which provides the exact value obtained from the deterministic logistic equation.
The previous results are general enough to hold also for averages of the form It is easy to see that Eq. (24) applies. We just substituteλ k = λ k + tσ 2 and multiply by a factor: r n (µ, σ, t) = e (µ− σ 2 2 )t t n e σ 2 2 t Γ(n)s n (λ 1 , · · · ,λ n ) = e µt t n Γ(n)s n (λ 1 , · · · ,λ n ). The above expression will be used in the following section as an application to averages in the logistic stochastic differential equation. Appendix A contains the analytical values of the functions s k (λ 1 , · · · , λ k ) and s k (λ 1 , · · · ,λ k ) for the cases n ≤ 4. It is immediate that, when t ≈ 0, we have This facts implies that Eq. (26) is analytic in t = 0. Fig. 2 is a plot of s n (t) for µ = 1, σ = 0.1.

III. AVERAGES OF LOGISTIC SDE IN PERTURBATION THEORY
As an application of the formula (26), we focus on the solution of the logistic stochastic differential equation motion, given by (we follow [10,11]) We evaluate the average of the solution x(t), x = 100, and with dt = 10 −3 , solved numerically using a stochastic Euler method and averaged over 1000 simulations (solid red), versus the analytical solution obtained at k-th order, considering ξ = x0 x the perturbative parameter. We can observe that at higher order we obtain a solution closer to the one simulated. and observe that it involves the moments of (26). In the limitx → ∞, Eq. (30) reduces to the average of the geometric Brownian motion. We consider now truncations of the mean of x(t) at k-th order, x(t) k = x 0 e µt k n=0 (−1) n (t µx 0 x ) n Γ(n)s n (λ 1 , · · · ,λ n ), and compare the truncated solution to the one obtained numerically. In Fig. 3, we plot x(t) k for k = 0, 1, . . . , 4 obtained for µ = 1, σ = 0.1,x = 100 by means of a stochastic Euler method with dt = 10 −3 , and the averages then obtained using Monte Carlo over 1000 samples. We observe that the higher the order of the approximation, the closer we are to the solution obtained by Monte Carlo.

IV. CONCLUSIONS
In this paper, we have presented new exact formulas for the moments of the integrated exponential Brownian motion in terms of sums and determinants, and based on recent results obtained in [9]. We described a simple graphical method to evaluate them, based on a recurrence relation.
Exact formulas were proved in [5] in terms of polynomials. In this paper however, we have taken an alternative route based on combinatorics. After realizing that the mean can be evaluated exactly using the properties of Gaussian integrals, and after observing that these moments feature a recurrence relation, we have shown that exact expressions can be obtained via a combinatorial ar-gument. These exact expressions were then observe to be equivalent to evaluating the determinant of a specific linear operator which depends on the order of the moment.
To complete the presentation, we have applied the formulas to the exact solution of the logistic stochastic differential equation. There, the evaluation of the ensemble expectation values of certain observables can be carried out with our method, via Taylor expansion. In particular, the comparison of the mean solution obtained by means of Monte Carlo simulations 1 with our method shows that our formula permits to approximate to a higher precision some ensemble averages of properties of the solution of the stochastic differential equation.
Aknowledgements F. C. would like to thank the London Institute for Mathematical Sciences for hospitality. S. S. is supported by the Royal Society and EPSRC.