We consider stochastic differential equations with additive noise and conditions on the coefficients in those equations that allow a time singularity in the drift coefficient. Given a maximum step size, , we specify variable (adaptive) step sizes relative to which decrease as the time node points approach the singularity. We use an Euler-type numerical scheme to produce an approximate solution and estimate the error in the approximation. When the solution is restricted to a fixed closed time interval excluding the singularity, we obtain a global pointwise error of order . An order of error for any is obtained when the approximation is run up to a time within of the singularity for an appropriate choice of exponent . We apply this scheme to Brownian bridge, which is defined as the nonanticipating solution of a stochastic differential equation of the type under consideration. In this special case, we show that the global pointwise error is of order , independent of how close to the singularity the approximation is considered.

1. Introduction

Numerical approximation methods for stochastic differential equations (SDEs) are well developed for SDEs with coefficients satisfying Lipschitz conditions. One of basic methods is the Euler-Maruyama algorithm which is presented in [1] and also well explained by Kloeden and Platen [2]. It is developed using Taylor series expansions and the Itô formula and is analogous to the Euler algorithm for the deterministic case. The Euler method with constant step sizes gives the global error when drift and diffusion terms satisfy global Lipschitz and linear growth conditions. In addition to methods with constant step sizes, methods with variable (adaptive) step sizes have been developed, for example, those in [37]. In general, these methods seek to improve the efficiency of the algorithms by adapting the step sizes to the state of the numerical solution at successive stages of the algorithm.

We consider stochastic differential equations with additive noise and an endpoint singularity with respect to the time variable in the drift term. With this time singularity, the global error estimate for the Euler-Maruyama scheme can only be applied up to a fixed time before the singularity. Using fixed step sizes to generate a numerical approximation for a time closer to the singularity comes at a great cost in efficiency due to the much smaller step size required to produce the same order of global error. We seek to increase the efficiency of the algorithm by using variable step sizes adapted to the shape of the singularity. In making this adaptation, we must also ensure that the step sizes do not become so small that the sum of the steps cannot reach times arbitrarily close to the singularity.

We consider a stochastic differential equation of the form where is a constant, is standard Brownian motion, and satisfies linear growth and Lipschitz conditions with bounds for derivatives up to the second order. For example, suppose that, uniformly in on , where denotes the partial derivative of with respect to the th variable.

Note that the singularity in the drift term stops us from using the standard results of numerical analysis. To overcome this problem, we use variable step sizes and stop the approximation at where . We estimate bounds for the pointwise local error and we use a lemma [8] to estimate the global error. We fix a maximum step size and produce variable steps and node points depending on in order to produce an estimate of the global error of order for fixed . The same algorithm and analysis using fixed step sizes fail to produce an estimate of global error. We further show that there is a choice of depending on so that the global error is of order for positive and converges to 0 as converges to 0. In this analysis we also obtain an estimate of the number of node points required for a given and . The trade-off for producing our global error estimate is that this number grows faster than , the corresponding number of node points for fixed step size approximations.

In Section 4, we apply this method to the process well known as the Brownian bridge (BB). Here we use the version of Brownian bridge that is the strong, adapted solution of a stochastic differential equation of the form (1a)-(1b) [9, p. 75]. BB is a stochastic process that starts at a specific point and converges to a specific point at a given time. BB is widely used in Monte Carlo simulations and many other applications [10, 11]. Most of the applications use simpler anticipating versions of BB like to avoid the approximation problem that occurs from the singularity in the drift term. We apply our methods to numerically approximate Brownian bridge and obtain even better estimates of global error than those for our general problem. For more properties of Brownian bridge, one can refer to [12, p. 86–89].

Brownian bridge from to is given by the following SDE [9, p. 75]. For the sake of simplicity, we use and later in our analysis. For fixed , , Then the solution of the above SDE is given by Also, It is clear that BB satisfies conditions given in (2a)–(2d). It can also be proven that BB satisfies the bounded moment condition, (2e), but we will show that our proposed numerical method can be applied directly to the Brownian bridge process to obtain the improved global error estimate without using the assumption of bounded moments.

In Section 2, we discuss the Euler-Maruyama method and some lemmas used in this paper. In Section 3, we propose the numerical approximation with variable (adaptive) step sizes based on the Euler-Maruyama algorithm and we analyze the error. In Section 4, we apply these methods to approximate the specific example, Brownian bridge, and we show that the order of the global error estimate is better for BB than for the general case.

2. Euler-Maruyama Algorithm and Fundamental Lemmas

Let be standard Brownian motion, and let be the natural filtration. We use the following general stochastic differential equation and Itô formula to present the development of Euler-Maruyama algorithm: where .

The same thing can be given in integral notation assuming all the integrals exist [9, p. 22]. The second integral is the Itô integral

Theorem 1 (Itô Formula [9, p. 46]). Let be an Itô process given by where and are measurable, -adapted, and integrable (resp., square integrable) on every interval , almost surely. Let . Then is again an Itô process and

Consider a sequence of times and the corresponding step sizes . Following [13, p. 56], we use the Itô formula to expand the drift and diffusion terms in (7), where , , , , , and . Because is standard Brownian motion, forms a sequence of independent Gaussian random variables with mean zero and .

The first three terms of expansion generate the Euler-Maruyama algorithm with variable step sizes [2, p. 190]. The double integrals yield the local truncation error. So we define the numerical scheme by

With suitable conditions on the coefficient function and with constant, it can be shown that mean square local truncation error is of order as . Let . When step sizes are fixed at , it follows that [13, p. 60] It can be further shown that the global error as [13, p. 61]. In the process of estimating the global error, the following lemma is used. This is an easy extension of a fundamental lemma of numerical analysis [8, p. 189].

Lemma 2. If and are real numbers and is a sequence with s.t. for all , then for all .

3. Numerical Approximation for SDE with Singularity

We return to the SDE of the form (1a)-(1b) presented in Section 1: where is a constant and satisfies linear growth and Lipschitz conditions with bounds for derivatives up to the second order. In particular, we suppose that, uniformly in on , It follows then that These conditions follow, for example, from the conditions given in (2a)–(2e).

In this case, and , and thus Note that the derivatives of are not bounded. Nevertheless, on each subinterval where , we have the expansion from (10) and the numerical algorithm from (11), respectively: Label the two integral terms as and , respectively, and let . We have Using the relation we have We have also used the inequality and the Lipschitz condition (17) here. Since is a martingale with respect to the natural filtration, . Now we look at the second moments of the integral terms and . First we use Cauchy-Schwartz inequality twice on and the bounds imposed in (16a)–(16d) for some , , and that do not depend on or . It follows that where is independent of and . Now for term we use the Cauchy-Schwartz inequality and the Itô isometry where is independent of and .

We use Cauchy-Schwartz inequality and the inequality with to estimate :

For simplicity, we assume . For fixed step size, , there is no uniform control over the local truncation error since, close to , . If we fix and only consider , then we have the following inequality: where and are some constants which do not depend on or . Using Lemma 2, we have For fixed , this is the original case where the coefficients of the SDE have bounded derivatives. However, if , this estimate shows no decrease in global error unless decays exponentially relative to .

To produce a better estimate and a more efficient algorithm, we consider variable step sizes defined as follows. First fix , . Then define step sizes and node points using : This implies that This specific definition allows us to estimate the global error using , and therefore we can control to control the global error. Using these step sizes, we have the error estimates where , , , and are positive constants which do not depend on or . Using these estimates in (23), we have where and are constants which do not depend on or . Using Lemma 2, we have the following.

Theorem 3. Given that the SDE in (15a)-(15b) satisfies the assumptions (16a)–(16d) and , then the above constructed algorithm (20) with variable step sizes (32b) has , as , uniformly in for , and thus the global pointwise error for the above proposed algorithm is of order .

Proof. If we have steps, (32b) gives , and thus Then by using Lemma 2 on (34) and this we have where and are constants that do not depend on , , or .

Corollary 4. Given a positive constant and the SDE (15a)-(15b) satisfying the assumptions in Theorem 3, the above constructed numerical algorithm has global pointwise error , as , when , where is the constant in (34) which only depends on the bounds imposed on and the moments of the process.

Proof. Letting , we have

With a simple example, we show that we cannot expect a much better result. Consider the SDE: where . This equation has solution The variable step size algorithm generates the approximation It follows then that the maximum value of occurs at where and Thus if , as and is thus for all .

4. Numerical Approximation of Brownian Bridge with Variable Step Sizes

We apply the above approximation to the Brownian bridge introduced in (3a)-(3b). satisfies the Lipschitz and linear growth conditions assumed. We also show that we can relax the bounded moments conditions for this specific case. We use the BB starting at and converging to as . Consider the SDE: where is the standard Brownian motion. Then using (10) we get Note that term is zero and that gives us the opportunity to relax the bounded moments condition Using (11) we get Now define to get Proceeding as before we have that the second moment of the integral term is of order . Now using (22) and analyzing it in the same way as in the previous section, we get For we get which gives the global error using Lemma 2. Note that, in this case, for all , independent of . However, making small relative to increases the number of steps significantly. For example, if we set , it takes about twice as many steps to reach as in the case .

For this example, we can verify the error explicitly. Using Algorithm (45), we have Using this with the explicit solution (4), we calculate and thus the maximum value of as .

Conflict of Interests

The authors declare that they have no conflict of interests regarding the publication of this paper.


The authors would like to thank the referee for insightful comments and help in improving this paper.