Calculation of Probability of the Exit of a Stochastic Process from a Band by Monte-Carlo Method: A Wiener-Hopf Factorization

: This paper considers a method of the calculation of probability of the exit from a band of the solution of a stochastic di ﬀ erential equation. The method is based on the approximation of the solution of the considered equation by a process which is received as a concatenation of Gauss processes, random partition of the interval, Girsanov transform and Wiener-Hopf factorization, and the Monte-Carlo method. The errors of approximation are estimated. The proposed method is illustrated by numerical examples.


Introduction and Related Work
The paper is dedicated to the calculation of probability of the exit from a band. This problem is very important in applications of the theory of stochastic processes, such as financial mathematics (barrier options), communication lines (signals), and so on.
We consider a method of the calculation of probability of the exit from a band of the solution of a stochastic differential equation dX t = f (X t )dt + ϕ(X t )dW t with an initial condition X 0 , where W is a standard Wiener process. It is supposed that sufficient conditions of the existence and uniqueness of a strong solution of the equation with a given initial condition are satisfied. The solution is considered on a finite interval [0, T]. The desired probability is P out = 1 − EI {(X T >a)∧(X T <b)} , where The main contribution of the paper is the proposed numerical method. The method is based on (1) the approximation of the solution of the considered equation by a process which is received as a concatenation of Gauss processes; (2) random partition of the interval; (3) Girsanov transform and Wiener-Hopf factorization; and (4) Monte-Carlo method, which is the calculation core of the approach [1].
The applied random partition is proposed by Carr [2], who used it for calculation of the price in an American option. In our case, the Carr partition permits to use the Wiener-Hopf factorization [3] for the calculation of supremum and infimum processes for a Lévy process. This mode, together with the Monte-Carlo method, was applied in [4][5][6] to Lévy processes. In our case the solution of the considered equation is not a Lévy process and the direct use of the mode is impossible. A preliminary approximation of the solution is required.
Namely, the following approximation of the solution is used: Here W (i) are independent copies of a standard Wiener process W, µ i = f (H i ), σ i = ϕ(H i ), H 0 = X 0 ; or a refined approximation: In these formulas, {τ i } is an increasing sequence of the random points on the time axis. This approximation differs from a classic one. First of all, it is not a piecewise-constant or a linear approximation of the solution of a stochastic differential equation [7]. We will show later that this approach gives some advantages in the estimation of errors.
There are some references on the use of Brownian bridges in an approximation process. For example, an approximation based on the Brownian bridges was used for calculation of the prices of barrier options [8,9]. However, we think it impossible to use Brownian bridges for the calculation of the probability of the exit of a stochastic process from a band.
The approximation process permits to calculate an approximate value of the probability of the exit of a stochastic process from a band as follows: The Girsanov transform is a standard tool of the analysis of stochastic differential equations [10]. In our case, the Girsanov transform simplifies the modeling of paths and the estimation of approximation errors.
The numerical calculations were made for a geometrical Brownian motion and Ornstein-Uhlenbeck process [11]. The geometrical Brownian motion is chosen due to the possibility of an analytical calculation of the probability of the exit of a stochastic process from a half-band. Respectively, it is possible to compare the numerical results with the exact solution. The Ornstein-Uhlenbeck processes are widely used as a model in econometrics, for example, in the models of stochastic volatility describing financial indices [12,13]. If the volatility exceeds a threshold value, then the economic system is unbalanced. On the contrary, if the volatility is below a threshold value, then it signalizes about stagnation.
It is important to establish the following connection. There is a notion of homeostasis that means a dynamic stability of complex systems. A system satisfies to the condition of homeostasis if the values of all its key state variables belong to a given range during the whole time of its functioning. A comprehensive mathematical formalization of this notion is proposed by J.-P. Aubin in the viability theory [14]. An integral approach to the sustainable management based on the notion of homeostasis is presented in [15].
In this connection notice that the probability of the exit from a band gives a stochastic characteristic of the possible violation of the homeostasis. One of the possible interpretations was given for volatility.
The rest of the paper is organized as follows. Section 2 considers the Girsanov transform for the simplification of the analyzed equation, as well as Carr partition and Wiener-Hopf factorization. Section 3 gives estimates of the possible errors. In Section 4, the Wiener-Hopf factorization is used without the Girsanov transform. Section 5 contains the results of numerical calculations and their analysis. Section 6 concludes and presents a future work.

Girsanov Transform
Suppose that ϕ(x) > 0. The considered equation takes the form ϕ(X t ) dt + dW t . Let X t be a strong solution of the Equation (1). If the factor ϕ satisfies to the Lipschitz condition, then the existence and uniqueness of the strong solution of the Equation (1) are guaranteed [7]. Consider Given EZ T = 1 the following result holds.
Theorem 1 [16]. The process W * t is a Wiener process on the interval [0, T] respective to the measure P * : dP * T = Z T dP T .
The measure P T is generated by the Wiener process W t , t ∈ [0, T] on a Borel sigma-algebra B C[0, T] of the continuous functions on the interval [0, T].
As follows from the Theorem 1, if W * is a standard Wiener process respective to the measure P * , then the process W t is a standard Wiener process respective to the measure dP T = Z T dP * where A sufficient condition for the equalities EZ T = 1 and E * Z T = 1 hold is the following: If the condition is satisfied, then for the calculation of a mathematical expectation, we can use the equality EI {(X T >a)∧(X T <b)} = EI {(X T >a)∧(X T <b)} and the measure substitution, therefore: where X t is a solution of the Equation (2). Assume that τ is a Carr partition of the interval [0, T]. For the Carr partition, h i = τ i − τ i−1 is a sequence of independent and equally distributed random variables with a common exponential distribution law with a parameter λ n = n T . It is evident that the mathematical expectation Eτ n = T and the variance Dτ n = T 2 n . As lim n→∞ Dτ n = 0, then lim n→∞ τ n = T in the mean quadratic sense. The approximation process for the Equation (1) is determined by the equalities A substitution in (2) instead of the process X t , the process H t from (3) gives an approximate equality.
Now we need the random values H τ n and H τ n , which are calculated by the recurrent formula where we use denotations H i = sup As any Wiener process is a Lévy process, then the Wiener-Hopf factorization is applicable [3], according to which where ξ 1 иξ 2 are independent sequences of the independent and equally distributed by an exponential law random values with the parameter ν n = √ 2λ n . For the calculation of the multiplier connected with the Girsanov transform, we use an approximate formula.
The desired mathematical expectation is calculated approximately as follows: Note: the expression is a multiplier of the Girsanov transform for the original approximation process.

Errors Estimation
The algorithm contains two types of errors: An error of the calculation of approximation of a process X by a process H and an error of the calculation of the multiplier in the Girsanov transform. Let us estimate the latter one. First, we can estimate: Suppose that the function ϕ(x) satisfies the Lipschitz condition on an interval [a, b] with a constant L, the function Similarly, D * ∆ n ≤ (LqT) n 2 , and ∆ n tends to zero in L 2 -norm. Denote by U n and U n the random respectively, and estimate E * exp U n − exp U n I {H τn <b∧H τn >a} .
From the Cauchy-Buniakovsky inequality, we receive: Then, for the first multiplier, Consider the second multiplier E * (exp ∆ n − 1) 2 I {H τn <b∧H τn >a} . Use the Taylor series for the The series 1 + ∞ k=1 x k (k+1)! converges on the interval (−∞, ∞). As P(∆ n = ∞) = 0 ( E * ∆ 2 n → 0 ), then the sum S(∆ n ) < ∞ uniformly by ∆ n . Therefore, E * exp U n − exp U n I {H τn <b∧H τn >a} → 0 when n → ∞ . Now consider the error of approximation. The Girsanov transform permits to use the fact that both processes X and H are martingales under the measure P * in respect to the filtration F W * . As the processes are considered in a band, we can affirm that the martingales are quadratic integrated. Consider the process V = X − H, which is also a quadratic integrated martingale. Then the process V 2 is a submartingale, and for all t ∈ [0, T], the inequality holds: E * V 2 t ≤ E * V 2 T . Now use the formula from [17]. For the application of this formula an additional condition is required: a derivative of the function ϕ(x) should be bounded. If the inequalities hold: Come to the limit in both parts of the formula (8). As τ n → T in L 2 -norm and E * V 2 t is a continuous function then lim n→∞ V 2 t = 0 uniformly by t on the interval [0, T].

Wiener-Hopf Factorization
The process is a Lévy process, which is why we can apply the Wiener-Hopf factorization without the Girsanov transform directly to the process The factorization leads to the recurrent equations of the type (4), (5). The difference is in the sequences ξ 1 and ξ 2 in the formula (5). They are independent with independent elements distributed by the exponential laws with different parameters: ν 1,n for the sequences ξ 1 and ξ 2 , respectively. The approximate formula for the calculation of the conditional mathematical expectation looks simpler than the formula (6): However, its numerical complexity is almost the same as (6), because each iteration requires to calculate 2n parameters for two exponential laws. It is difficult to receive any estimations on errors of the direct method, as the classic results are connected to the uniform partition of the interval [10].

Examples and Numerical Calculations
The geometric Brownian motion is a process X t = X 0 exp − σ 2 2 t + σW t , which is a solution of the stochastic differential equation dX t = σX t dW t . Define the process Y t = − σ 2 2 t + σW t . Calculate the probability of the exit from a band: [0, b]. It is evident that the random events sup Here, the density p τ (s) = [18]). The numerical results for these examples are presented in Table 1. We used the following values of the parameters: b = 0.1; σ = 0.4. An Ornstein-Uhlenbeck process X t = c + (X 0 − c) exp(−αt) + σ t 0 exp(−α(t − s))dW s is a strong solution of the differential equation dX t = α(c − X t )dt + σdW t with an initial condition X 0 and parameters α > 0 and σ > 0. The Ornstein-Uhlenbeck process is a return one. The numerical results are given in the Tables 2 and 3. We used the following values of parameters: a = 0.1; b = 0.5; c = 0.3; Tables 2 and 3 shows more coincidences than differences. However, the research results give more trust to the Table 2.

Conclusions and Future Work
The proposed combined method of calculation of the probability of the exit of a stochastic process from a band is based on a natural approximation of the diffusion process. The Girsanov transform simplifies the approximation and makes it exact for a constant diffusion parameter. However, an approximate calculation of the Girsanov transform gives an additional error which was estimated. The methods permit to calculate successfully the probability of the exit of a stochastic process from a half-band. The Wiener-Hopf factorization, together with the Monte-Carlo method, is a popular tool for calculations in Lévy processes. A special merit in the application of Wiener-Hopf factorization in models under Levy processes belongs to S. Levendorsky and O. Kudriavtsev [19]. The method proposed in the paper is applicable for diffusion processes that are not Lévy processes. The considered examples are illustrative but quite representative.
The future research will be connected with other ways of the random partition of an interval and other methods of approximation, which us allow to avoid the negative properties of the Monte-Carlo method.
The results can be used for the solution of sustainable management problems in financial markets and other domains of application. Namely, the exit from a band which witnesses about the violation of the homeostasis in a stochastic sense can be prevented by some special control efforts.