A Robust Weak Taylor Approximation Scheme for Solutions of Jump-Diffusion Stochastic Delay Differential Equations

Abstract

Stochastic delay differential equations with jumps have a wide range of applications, particularly, in mathematical finance. Solution of the underlying initial value problems is important for the understanding and control of many phenomena and systems in the real world. In this paper, we construct a robust Taylor approximation scheme and then examine the convergence of the method in a weak sense. A convergence theorem for the scheme is established and proved. Our analysis and numerical examples show that the proposed scheme of high order is effective and efficient for Monte Carlo simulations for jump-diffusion stochastic delay differential equations.

1. Introduction

Stochastic delay differential equation models play a very important role in the study of many fields such as economics and finance, chemistry, biology, microelectronics, and control theories. Models of this type can more accurately describe many phenomena in the real world by taking into account the effect of time delay and random noise. By time delay, it means a certain period of time is required for the effect of an action to be observed after the moment when the action takes place. This phenomenon exists in most systems in almost any area of science; for example, a patient shows symptoms of an illness days or even weeks after he/she was infected. Similarly, random noises appear in almost all real world phenomena and systems, for example, the motion of molecules and the price of assets in financial markets. Hence, study of SDDEs is an important undertaking in order to understand real world phenomena and systems precisely.

Over the last couple of decades, a lot of work has been carried out to study differential equations with delay and/or random noises, for example, [1–5]. The best-known and well-studied theory and systems include the delay differential equations (DDEs) presented by Kolmanvskii and Myshkis [6] and their stochastic generalizations and the stochastic delay differential equations (SDDEs) established by Mohammed [7, 8], Mao [9, 10], and Mohammed and Scheutzow [11]. Other SDDEs theories of interest include, for instance, the so-called SDDEs with Markovian switching and Poisson jumps. These models have been investigated in the literature [12–14].

Analytical solutions of SDDEs can hardly be obtained. It is thus important to develop and study discrete-time approximation methods for solving SDDEs. Discrete-time approximations may be divided into two categories: weak approximations and strong approximations [15]. Some implicit and explicit numerical approximation methods for SDDEs in strong approximation sense were derived by Küchler and Platen [16]. Weak numerical methods for SDDEs have been studied by Küchler and Platen. The Monte Carlo simulation method has also been developed as a powerful simulation method for SDEs, while weak numerical approximations are required for Monte Carlo simulation [9, 17–20].

In this paper, we extend a fully weak approximation method for SDDEs to the type of jump-diffusion SDDEs: (𝑡)(𝑡)(𝑡)𝑑𝐗=𝐚𝑡,𝐗,𝐗𝑡−𝛾𝑑𝑡+𝐛𝑡,𝐗,𝐗𝑡−𝛾𝑑𝐖_𝑡+_𝜺𝐜(𝑡,𝐗(𝑡⁻),𝜐)𝑝_𝜑(𝑑𝜐,𝑑𝑡)(1) subject to the initial condition 𝑋(𝜃)=𝜒(𝜃)for𝜃∈−𝛾,0,(2) where 𝑡∈[0,𝑇], 𝛾 is the time delay which is assumed to be constant at all time, 𝐖_𝑡={(𝑊¹_𝑡,…,𝑊^𝑚_𝑡),𝑡∈[0,𝑇]} is an 𝒜-adapted 𝑚-dimensional Wiener process, and 𝑝_𝜑 denotes the Poisson random measure. Also here we denote by 𝐗(𝑡⁻) the almost sure left-hand limit of 𝐗(𝑡). The coefficient 𝐚(𝑡,𝑥,𝑥^𝑟)∶[0,𝑇]×ℝ^𝑑×ℝ^𝑑→ℝ^𝑑 and 𝐜(𝑡,𝑥,𝜐)∶[0,𝑇]×ℝ^𝑑×𝜺→ℝ^𝑑 are 𝑑-dimensional vectors of Borel measurable functions. Further, 𝐛(𝑡,𝑥,𝑥^𝑟) defined on [0,𝑇]×ℝ^𝑑×ℝ^𝑑 is a 𝑑×𝑚-matrix of Borel measurable function.

The main contributions of this work include construction of a numerical approximation method for weak solutions of SDDEs with jump-diffusion, and establishment of a theoretical result for the convergence of the scheme.

The remaining part of this paper is organised as follows. In the following section, we give the conditions for the existence and uniqueness of the solution of the jump-diffusion SDDEs (1) which is applicable to any cases where solution exists, and present various lemmas to be used later for the proof of the convergence theorem. We then introduce, in Section 3, a general weak approximation scheme, where the simplified stochastic Taylor approximation scheme with order 𝛽 is constructed; followed by a convergence theorem and its proof. In Section 4, we give a numerical example to demonstrate the application and the convergence of the numerical scheme.

2. Preliminaries

In this section, we present some basic concepts and definitions to be used in later sections, then establish the conditions for the existence and uniqueness of solution to the Jump-diffusion SDDEs (1), and then give some lemmas to be used later for the proof of the convergence theorem. In this work, we assume that the 𝑑-dimensional vector valued function for the initial condition, 𝜒={𝜒(𝑠),𝑠∈[−𝛾,0]}, is right continuous and has left-hand limits.

From (1), we have the following integral form of the jump-diffusion equation SDDE: 𝐗(𝑡)=𝐗(0)+^𝑡₀𝐚(𝜏)+𝜏,𝐗,𝐗𝜏−𝛾𝑑𝜏_𝑡0𝐛𝜏,𝐗(𝜏),𝐗𝜏−𝛾𝑑𝐖_𝜏+^{𝑝𝜑(𝑡)}_𝑖=1𝐜𝜏_𝑖,𝐗𝜏⁻_𝑖,𝜉_𝑖,(3) where (𝜏_𝑖,𝜉_𝑖), for 𝑖∈{1,2,3,…,𝑝_𝜑(𝑡)}, are a sequence of pairs of jump times and corresponding values generated by the Poisson random measure 𝑝_𝜑.

Definition 1. Given a filtered probability space (Ω,𝒜,𝒜,𝑃), a stochastic process given by 𝐗={𝐗(𝑡),𝑡∈[−𝛾,𝑇]} is known as a solution of (1) subject to the initial condition (2) if 𝐗 is 𝒜-adapted, the integrals in the equation are well defined, and equalities (3) and (2) hold almost surely. Moreover, if any two solution processes 𝐗^(𝑖)={𝐗^(𝑖)(𝑡),𝑖∈1,2} are indistinguishable on [−𝛾,𝑇] with the same initial segment 𝜒 and the same path on [0,𝑇], and 𝑃sup_{[]𝑡∈0,𝑇}𝑋⁽¹⁾(𝑡)−𝑋⁽²⁾(𝑡)>0=0,(4) where ‖⋅‖ is the Euclidean norm, then if (1) has a solution, it is a unique solution for this initial value problem.

To guarantee the existence of a unique solution of the jump-diffusion SDDE (1), we assume that the coefficients of (1) satisfy the following Lipschitz conditions: 𝐚𝑡,𝐲₁,𝐳₁−𝐚𝑡,𝐲₂,𝐳₂≤𝐶₁𝐲₁−𝐲₂+𝐳₁−𝐳₂,𝐛𝑡,𝐲₁,𝐳₁−𝐛𝑡,𝐲₂,𝐳₂≤𝐶₂𝐲₁−𝐲₂+𝐳₁−𝐳₂,_𝜺𝐜𝑡,𝐲₁,𝐳₁,𝜐−𝐜𝑡,𝐲₂,𝐳₂,𝜐²𝜑(𝑑𝜐)≤𝐶₃𝑦₁−𝑦₂²+𝑧₁−𝑧₂²(5) for 𝑡∈[0,𝑇] and 𝐱₁,𝐱₂,𝐲₁,𝐲₂∈ℝ^𝑑, as well as the growth conditions |𝐚(𝑡,𝐲,𝐳)|≤𝐷₁(|𝐲|+|𝐳,|𝐛(1+|)𝑡,𝐲,𝐳)|≤𝐷₂(1+|𝐲|+|𝐳|)_𝜺,|𝐜(𝑡,𝐲,𝐳,𝜐)|²𝜑(𝑑𝜐)≤𝐷₃1+𝑦²+𝑧²(6) for 𝐲,𝐳∈ℝ^𝑑, 𝑡∈[0,𝑇].

We denote by ℂ=ℂ([−𝛾,0],ℝ^𝑑) the Banach space of all 𝑑-dimensional continuous functions 𝜂 on [−𝛾,0] equipped with the supremum norm ‖𝜂‖_ℂ=sup_{𝑠∈[−𝛾,0]}|𝜂(𝑠)|. Furthermore, we suppose that the set 𝐿₂(Ω,ℂ,𝒜₀) of ℝ^𝑑-valued continuous process 𝜂={𝜂(𝑠),𝑠∈[−𝛾,0]} is 𝒜₀-measurable with 𝐸𝜂²_ℂ=𝐸sup_{𝑠∈[−𝛾,0]}𝜂(𝑠)²<∞.(7)

Following the work in [7, 10], the following theorem can be established for the existence and uniqueness of a solution to the problem defined by (1) and (2).

Theorem 2. Suppose that the Lipschitz conditions and the growth conditions are satisfied, and the initial condition 𝜒 is in 𝐿₂(Ω,ℂ,𝒜₀). Then (1) subject to the initial condition (2) admits a unique solution.

Now we present some lemmas to be used later for the proof of the convergence theorem. Consider a right continuous process 𝑌^Δ𝑙={𝑌^Δ𝑙(𝑡),𝑡∈[−𝛾,𝑇]}. 𝑌^Δ𝑙 is called a discrete-time numerical approximation with maximum step size Δ_𝑙, if it is obtained by using a time discretization 𝑡_Δ𝑙, and the random variable 𝑌^Δ𝑙_𝑡𝑛 is ℱ_𝑡𝑛-measurable for 𝑛∈{1,…,𝑁}. Further, 𝑌^Δ𝑙_𝑡𝑛+1 can be expressed as a function of 𝑌^Δ𝑙_𝑡−𝑙,𝑌^Δ𝑙_{𝑡−𝑙+1},…,𝑌^Δ𝑙_𝑡𝑛 and the discrete-time 𝑡_𝑛.

Because of dealing with the approximation of solutions of jump-diffusion SDDEs, we introduce a concept of weak order convergence due to Kloeden and Platen [15].

Definition 3. A discrete-time approximation 𝑌^Δ𝑙 converges weakly towards 𝑋 at time 𝑇 with order 𝛽>0 if for each 𝑔∈𝒞_𝑝 there is a constant 𝐶, independent of Δ_𝑙, such that 𝐸𝑔(𝑋(𝑇))−𝐸𝑔𝑌^Δ𝑙(𝑇)≤𝐶Δ_𝑙^𝛽,(8) where 𝒞_𝑝 denotes the set of all polynomials 𝑔∶ℝ^𝑑→ℝ.

We now give some auxiliary results to prepare for the proof of the weak convergence theorem to be presented.

Lemma 4. For 𝑛∈{−𝑙+1,…,0,1,…,𝑁} and 𝐳∈ℝ^𝑑, we have 𝐸𝑢𝑛,𝐗_𝑛^{𝑛−1,𝑧}₋()+−𝑢𝑛−1,𝐳_{𝑛𝑛−1𝜺}𝐋_𝜐⁻¹𝑢𝜏,𝐗_𝜏^{𝑛−1,𝑧}𝜑(𝑑𝜐)𝑑𝜏∣𝒜_𝑛−1=0(9) for (𝜏,𝐳)∈[−𝛾,𝑡]×ℝ^𝑑 and 𝑢(𝜏,𝐳)=𝐸(𝑔(𝐗_𝑇^𝜏,𝑧)𝒜_𝜏).

The proof of the lemma for the case with no delay was established by Platen and Bruti-Liberati [21], and a similar procedure can be used for the proof of this lemma.

Lemma 5. Given 𝑝∈{1,2,3,…}, there is a bounded constant 𝑀 satisfying 𝐸𝐗_𝑛^{𝑛−1,𝑧}₋−𝐳^2𝑞∣𝒜_𝑛−1≤𝑀1+|𝐳|^2𝑞Δ_𝑙^𝑞(10) for 𝑞∈{1,…,𝑝} and 𝑛∈{−𝑙+1,…,0,1,…,𝑁}.

The proof of (10) can be obtained by following that of a lemma for SDEs with jumps but with no delay in [15]. The following results are similar to what was given in Mikulevičius and Platen [22].

Lemma 6. Given 𝑝∈{1,2,3,…}, there is a finite constant 𝑀 satisfying 𝐸sup_{−𝛾≤𝑡≤𝑇}|𝜻(𝑡)|^2𝑞≤𝑀1+𝐘₀^2𝑞(11) for every 𝑞∈{1,…,𝑝}.

Lemma 7. Given 𝑝∈{1,2,…}, there is 𝑟∈{1,2,3,…} and a bounded constant 𝑀 satisfying 𝐸𝐅_𝐩𝜻(𝑧)−𝐘^Δ𝑙_𝑧^2𝑞+𝐅_𝐩𝐗^{𝑧,𝑌Δ𝑙𝑧}_𝑧−𝐘^Δ𝑙_𝑧^2𝑞∣𝒜_𝑧≤𝑀1+𝐘^Δ𝑙_𝑧^2𝑟Δ_𝑙^𝑞𝑘(12) for each 𝑞∈{1,…,𝑝}, 𝑘∈{1,…,2(𝛽+1)}, 𝐩∈𝑃_𝑘={1,2,...,𝑑}^𝑘, and 𝑧∈[−𝑟,𝑇], where 𝐹_𝐩(𝐲)=∏^𝑘_ℎ=1𝑦^𝑝ℎ for all 𝐲=(𝑦¹,...,𝑦^𝑑)^𝑇∈ℜ^𝑑 and 𝐩=(𝑝₁,...,𝑝_𝑘)∈𝑃_𝑘.

The proof of estimate (12) can be established by following Itô’s formula for SDEs with jumps but with no delay as in [15].

Lemma 8. For 𝐩∈𝑃_𝑘, there exist 𝑟∈{1,2,…} and a finite constant 𝑀 satisfying 𝐸𝐅_𝐏𝜻(𝑡)−𝐘^Δ𝑙_𝑛−1−𝐅_𝐏𝐗^{𝑛−1,𝑌Δ𝑙𝑛−1}_𝑡−𝐘^Δ_𝑡𝑛−1∣𝒜_𝑡≤𝑀1+𝐘^Δ𝑙_𝑛−1^𝑟Δ_𝑙^𝛽(13) for each 𝑘∈{1,…,2𝛽+1}, 𝑛∈{−𝑙+1,…,0,1,…,𝑁}, and 𝑡∈[𝑡_𝑛−1,𝑡_𝑛).

3. The Jump-Adapted Weak Taylor Approximation Scheme

In Monte Carlo simulations for functionals of jump-diffusion SDDEs, one uses numerical approximations evaluated only at discretization time. Here, we first give a jump adapted weak approximation Taylor scheme of order 𝛽, and then study the basic properties of the discrete Taylor approximation in a weak order sense.

First we define the jump-adapted time discretization. Let 𝑇>𝑟>0. The jump adapted time discretization used throughout this paper is (𝑡)_Δ=𝑡_𝑖∶𝑖=−𝑙,−𝑙+1,…,0,1,2,…,𝑁(14) with the maximum step size Δ_𝑙∈(0,1). We choose the time discretization in such a way that all jump times are at the nodes of the time discretization. If the discretization node 𝑡_𝑖 is not a jump time, then 𝑡_𝑖 is 𝒜_𝑡𝑖−1-measurable. Otherwise, 𝑡_𝑖 is 𝒜_𝑡𝑖−-measurable. Also, throughout the paper, we denote the set of all multiindices 𝛼 by ℳ_𝑚=𝑗₁,…,𝑗_𝑙∶𝑗_𝑖∈{}{}}∪{𝑣},0,1,2,…,𝑚,𝑖∈1,2,…,𝑙for𝑙∈𝒩(15) where the element 𝛼=(𝑗₁,𝑗₂,…,𝑗_𝑙) is called a multiindex of length 𝑙=𝑙(𝛼)∈𝒩 and 𝑣 has zero length. In the following, by a component 𝑗∈{0,1,2,…,𝑚} of a multi-index we refer to the integration with respect to the 𝑗th Wiener process in a multiple stochastic integral. A component with 𝑗=0 corresponds to integration with respect to time 𝑡.

Now define the following operators for the coefficient functions: 𝐿⁰=𝜕+𝜕𝑡^𝑑_𝑖=1𝑎^𝑖𝑡,𝑥,𝑥_𝑙𝜕𝜕𝑥^𝑖+^𝑑_𝑖=1𝑎^𝑖𝑡−𝑙,𝑥_𝑙,𝑥_2𝑙𝜕𝜕𝑥^𝑖_𝑙+12^𝑑𝑚_{𝑖,𝛾=1𝑗=1}𝑏^𝑖𝑗𝑡,𝑥,𝑥_𝑙𝑏^𝛾𝑗𝑡,𝑥,𝑥_𝑙𝜕²𝜕𝑥^𝑖𝜕𝑥^𝛾+12^𝑑𝑚_{𝑖,𝛾=1𝑗=1}𝑏^𝑖𝑗𝑡−𝑙,𝑥_𝑙,𝑥_2𝑙𝑏^𝛾𝑗𝑡−𝑙,𝑥_𝑙,𝑥_2𝑙𝜕²𝜕𝑥^𝑖_𝑙𝜕𝑥^𝛾_𝑙,𝐿^𝑘=^𝑑_𝑖=1𝑏^𝑖𝑘𝑡,𝑥,𝑥_𝑙𝜕𝜕𝑥^𝑖+^𝑑_𝑖=1𝑏^𝑖𝑘𝑡−𝑙,𝑥_𝑙,𝑥_2𝑙𝜕𝜕𝑥^𝑖_𝑙.(16)

A subset 𝒜∈ℳ is the hierarchical set, and its corresponding remainder set 𝐴(ℳ) is defined by 𝐴(ℳ)={𝛼∈ℳ_𝑚⧵𝒜∶−𝛼∈𝒜}. For each 𝛽∈1,2,3,…, we can then define the hierarchical set Γ_𝛽={𝛼∈ℳ_𝑚∶𝑙(𝛼)≤𝛽}. The weak Taylor method of order 𝛽 is then constructed as follows: 𝐘_(𝑛+1)⁻=𝐘_𝑛+_{𝛼∈Γ𝛽}𝑓_𝛼𝑛,𝐘_𝑛,𝐘_𝑛−𝑙𝐼_𝛼,𝐘_𝑛+1=𝐘_(𝑛+1)⁻+_𝜀𝑐𝑛,𝐘_(𝑛+1)⁻,𝜐𝑝_𝜑(𝑑𝜐,(𝑛+1)),(17) where 𝑓_𝛼()=𝑓()(𝛼)𝑡,𝑥,𝑢𝐿𝑡,𝑥if𝑙=0_𝑗1𝑓_−𝛼()(𝛼)𝑗𝑡,𝑥,𝑢if𝑙≥1,₁∈0,1,…,𝑚.(18)

The multiple stochastic integral is then defined recursively as follows: 𝐼_𝛼,𝑡=∫𝑡if𝑙=0_𝑡0𝐼_𝛼−,𝑧𝑑𝑧if𝑙≥1,𝑗_𝑙∫=0_𝑡0𝐼_𝛼−,𝑧𝑑𝑊^𝑗𝑙_𝑧if𝑙≥1,𝑗_𝑙∈1,…,𝑚,(19) where 𝛼− is obtained from 𝛼 by deleting its last component, while −𝛼 is obtained from 𝛼 by deleting its first component.

Now, we give the weak convergence theorem of the Taylor approximation with order 𝛽.

Theorem 9. Given 𝛽∈{1,2,…}, let 𝑌^Δ𝑙={𝑌^Δ𝑙_𝑛,𝑛∈[−𝑙,…,0,1,…,𝑁]} be the results obtained from the Taylor scheme (17) corresponding to (𝑡)_Δ with maximum step size Δ_𝑙∈(0,1). Suppose that 𝐸(|𝐗_𝜉|^𝑖)<∞ for 𝜉∈(−𝛾,0),𝑖∈{1,2,…}, and 𝐘^Δ_𝜉 converges to 𝐗_𝜉 weakly with order 𝛽∈{1,2,…}. Assume that the coefficients, 𝑎^𝑘, 𝑏^𝑘𝑗, 𝑐^𝑘, are in the space 𝒞_𝑃^2(𝛽+1)(ℝ^𝑑,ℝ), for 𝑗∈{1,2,…,𝑚} and 𝑘∈{1,2,…,𝑑}, and the coefficient functions 𝑓_𝛼, with 𝑓(𝑡,𝐲)=𝐲, satisfy the growth condition |𝑓_𝛼(𝑡,𝐲)|≤𝑀(1+|𝐲|), with 𝑀<∞, for all 𝑡<𝑇, 𝐲∈ℝ^𝑑, and 𝛼∈Γ_𝛽. Then for any 𝑔∈𝒞_𝑃^2(𝛽+1)(ℝ^𝑑,ℝ) there is a positive constant 𝐶, which does not depend on Δ, such that 𝐸𝑔(𝑋(𝑇))−𝐸𝑔𝑌^Δ𝑙(𝑇)≤𝐶Δ_𝑙^𝛽.(20)

Proof. For 𝛽∈{1,2,3,…} and 𝑔∈𝒞_𝑝^2(𝛽+1)(ℝ^𝑑,ℝ), consider the Itô process below: 𝐗_𝑡^𝑧,𝑦=𝐲+^𝑡_𝑧𝐚𝐗_𝑢^𝑧,𝑦𝑑𝑢+^𝑡_𝑧𝐛𝐗_𝑢^𝑧,𝑦𝑑𝐖_𝑢+^𝑡_𝑧𝜺𝐜𝐗^𝑧,𝑦_𝑢−,𝜐𝑝_𝜑(𝑑𝜐,𝑑𝑢).(21)
Then we can get 𝑢(0,𝑋₀)=𝐸(𝑔(𝑋^0,𝑋0_𝑇))=𝐸(𝑔(𝑋_𝑇)).
Define also the process 𝜁=𝜁(𝑡),𝑡∈(−𝛾,𝑇), by 𝜁(𝑡)=𝜁𝑡_𝑛+_{𝛼∈Γ𝛽}𝐼_𝛼𝑓_𝛼𝑛,𝜻_𝑛,𝜻_𝑛−𝑙+^𝑡_𝑛𝜀𝑐𝑛,𝜻_(𝑛+1)⁻,𝜐𝑝_𝜑(𝑑𝜐,(𝑛+1))(22) for 𝑛∈{−𝑙,−𝑙+1,…,0,1,…,𝑁−1},𝑡∈(𝑡_𝑛,𝑡_𝑛+1], 𝜁(0)=𝑌₀, and 𝜁(𝑡_𝑛)=𝑌_𝑡𝑛 for 𝑛∈{−𝑙,…,0,1,2,…,𝑁}.
By the definition of the functional 𝑢 and the terminal condition of the stochastic process 𝑋, we have 𝐻=𝐸𝑔𝐘^Δ𝑙_𝑇−𝐸𝑔𝐗_𝑇=𝐸𝑢𝑇,𝐘^Δ𝑙_𝑇−𝑢0,𝐗₀.(23)
Since 𝐘₀ converges towards 𝐗₀ weakly with order 𝛽, one has 𝐻≤𝐸^𝑁_{𝑛=−𝑙+1}𝑢𝑛,𝐘_𝑛−𝑢𝑛,𝐘_𝑛−+𝑢𝑛,𝐘_𝑛−−𝑢𝑛−1,𝐘_𝑛−1+𝐾Δ_𝑙^𝛽.(24)
By Lemma 4, we can write 𝐻≤𝐸^𝑁_{𝑛=−𝑙+1}𝑢𝑛,𝐘_𝑛−𝑢𝑛,𝐘_𝑛−+𝑢𝑛,𝐘_𝑛−−𝑢𝑛−1,𝐘_𝑛−1−𝑢𝑛,𝐗^{𝑛−1,𝑌𝑛−1}_𝑛⁻−𝑢𝑛−1,𝐘_𝑛−1+^𝑛_𝑛−1𝜺𝐋_𝜐⁻¹𝑢𝑧,𝐗^{𝑛−1,𝑌𝑛−1}_𝑧×𝜑(𝑑𝜐)𝑑𝑧+𝐾Δ_𝑙^𝛽.(25)
From the properties of stochastic integrals, we obtain 𝐸^𝑁_𝑛=−𝑙𝑢𝑛,𝐘_𝑛−𝑢𝑛,𝐘_𝑛−=𝐸^𝑇_−𝑟𝜺𝐋_𝜐⁻¹𝑢(𝑧,𝜻(𝑧))𝜑(𝑑𝜐)𝑑𝑧.(26)
Therefore, we have 𝐻≤𝐻₁+𝐻₂+𝐾Δ_𝑙^𝛽,(27) where 𝐻₁=𝐸^𝑁_{𝑛=−𝑙+1}𝑢𝑛,𝐘_𝑛−−𝑢𝑛,𝐘_𝑛−1−𝑢𝑛,𝐗^{𝑛−1,𝑌𝑛−1}_𝑛⁻−𝑢𝑛,𝐘_𝑛−1,𝐻(28)₂=𝐸^𝑇_−𝑟𝜺𝐿_𝜐⁻¹𝑢(𝑧,𝜻(𝑧))−𝐿_𝜐⁻¹𝑢𝑧,𝐘_𝑧−𝐿_𝜐⁻¹𝑢𝑧,𝐗^{𝑧,𝑌𝑧}_𝑧−𝐿_𝜐⁻¹𝑢𝑧,𝐘_𝑧×𝜑(𝑑𝜐)𝑑𝑧.(29)
In the following, we proceed to estimate 𝐻₁ and 𝐻₂ in Steps 1 and 2, respectively, and then complete the proof in Step 3.
Step 1. Let us assume that 𝑢 is so smooth that the deterministic Taylor expansion may be applied. Hence, by expanding 𝑑𝑢 in 𝐻₁, we get 𝐻₁=𝐸^𝑁_{𝑛=−𝑙+1}^2𝛽+1_𝑘=11𝑘!_{𝐏∈𝑃𝑘}𝜕^𝐏_𝑦𝑢𝑛,𝐘_𝑛−1×𝐹_𝐏𝐘_𝑛−−𝐘_𝑛−1+𝑅_𝑛𝐘_𝑛−−^2𝛽+1_𝑘=11𝑘!_{𝐏∈𝑃𝑘}𝜕^𝐏_𝑦𝑢𝑛,𝐘_𝑛−1×𝐹_𝐏𝐗^{𝑛−1,𝑌𝑛−1}_𝑛⁻−𝐘_𝑛−1+𝑅_𝑛𝐗^{𝑛−1,𝑌𝑛−1}_𝑛⁻,(30) where the remainder term is 𝑅_𝑛(𝐙)=12𝛽+1!×_{𝐩∈𝑃2(𝛽+1)}𝜕^𝐩_𝑦𝑢𝑛,𝐘_𝑛−1+𝜽_𝐩,𝑛𝐙−𝐘_𝑛−1×𝐹_𝐩𝐙−𝐘_𝑛−1,(31) where 𝜽_𝐩,𝑛(𝐙) is a 𝑑×𝑑 diagonal matrix with 𝜃^𝑘,𝑘_𝐩,𝑛(𝐙)∈(0,1) for 𝑘∈{1,2,3,…,𝑑} and 𝐙=𝐘_𝑛− and 𝐗^{𝑛−1,𝑌𝑛−1}_𝑛⁻, respectively.
Therefore, according to the properties of expectation and absolute value, we get 𝐻₁≤𝐸^𝑁_{𝑛=−𝑙+1}^2𝛽+1_𝑘=11𝑘!_{𝐩∈𝑃𝑘}𝜕^𝐩_𝑦𝑢𝑛,𝐘_𝑛−1×𝐹_𝐩𝐘_𝑛−−𝐘_𝑛−1−𝐹_𝐩𝐗^{𝑛−1,𝑌𝑛−1}_𝑛⁻−𝐘_𝑛−1×𝑅_𝑛𝐘_𝑛−+𝑅_𝑛𝐗^{𝑛−1,𝑌𝑛−1}_𝑛⁻≤𝐸^𝑁_{𝑛=−𝑙+1}^2𝛽+1_𝑘=11𝑘!_{𝐩∈𝑃𝑘}𝜕^𝐩_𝑦𝑢𝑛,𝐘_𝑛−1×𝐸𝐹_𝐩𝐘_𝑛−−𝐘_𝑛−1−𝐹_𝐩𝐗^{𝑛−1,𝑌𝑛−1}_𝑛⁻−𝐘_𝑛−1∣𝒜_𝑛−1+𝐸𝑅_𝑛𝐘_𝑛−∣𝒜_𝑛−1+𝐸𝑅_𝑛𝐗^{𝑛−1,𝑌𝑛−1}_𝑛⁻∣𝒜_𝑛−1.(32)
By (31), the Hölder inequality, and Lemma 7, we get 𝐸𝑅_𝑛𝐘_𝑛−∣𝒜_𝑛−1≤𝑀_{𝐩∈𝑃2(𝛽+1)}𝐸𝜕^𝐩_𝑦𝑢𝑛,𝐘_𝑛−1+𝜽_𝐩,𝑛𝐘_𝑛−×𝐘_𝑛−−𝐘_𝑛−1²∣𝒜_𝑛−1^1/2×𝐸𝐹_𝐩𝐘_𝑛−−𝐘_𝑛−1²∣𝒜_𝑛−1^1/2≤𝑀𝐸1+𝐘_𝑛−1^2𝑟+𝐘_𝑛−−𝐘_𝑛−1^2𝑟∣𝒜_𝑛−1^1/2×𝐸𝐘_𝑛−−𝐘_𝑛−1^4(𝛽+1)∣𝒜_𝑛−1^1/2≤𝑀1+𝐘_𝑛−1^2𝑟Δ_𝑙^𝛽+1.(33)
Similarly, by Lemma 5 and the Cauchy-Schwarz inequality, we have 𝐸𝑅_𝑛𝐗^{𝑛−1,𝑌𝑛−1}_𝑛⁻∣𝒜_𝑛−1≤𝑀_{𝐩∈𝑃2(𝛽+1)}𝐸𝜕^𝐩_𝑦𝑢𝑛,𝐘_𝑛−+𝜽_𝐩,𝑛𝐗^{𝑛−1,𝑌𝑛−1}_𝑛⁻×𝐗^{𝑛−1,𝑌𝑛−1}_𝑛⁻−𝐘_𝑛−1²∣𝒜_𝑛−1^1/2×𝐸𝐹_𝐩𝐗^{𝑛−1,𝑌𝑛−1}_𝑛⁻−𝐘_𝑛−1²∣𝒜_𝑛−1^1/2≤𝑀𝐸1+𝐘_𝑛−1^2𝑟+𝐗^{𝑛−1,𝑌𝑛−1}_𝑛⁻−𝐘_𝑛−1^2𝛾∣𝒜_𝑛−1^1/2×𝐸𝐗^{𝑛−1,𝑌𝑛−1}_𝑛⁻−𝐘_𝑛−1^4(𝛽+1)∣𝒜_𝑛−1^1/2≤𝑀1+𝐘_𝑛−1^2𝑟Δ_𝑙^𝛽+1.(34)
Now, from the Cauchy-Schwarz inequality, Lemmas 8 and 6, and inequalities (33) and (34), we obtain 𝐻₁≤𝐸𝐾^𝑁_{𝑛=−𝑙+1}1+𝐘_𝑛−1^2𝛾Δ_𝑙^𝛽≤𝑀Δ_𝑙^𝛽1+𝐸max_{−𝑙≤𝑛≤𝑁}𝐘_𝑛^2𝛾≤𝑀Δ_𝑙^𝛽1+𝐘₀^2𝛾≤𝐾Δ_𝑙^𝛽.(35)
Step 2. Now we estimate the term 𝐻₂ in inequality (27). By the jump coefficient 𝑐 and the smooth function 𝑢, applying the Taylor expansion yields 𝐻₂=𝐸^𝑇_−𝑟𝜀^2𝛽+1_𝑘=11𝑘!_{𝐩∈𝑃𝑘}𝜕^𝐩_𝑦𝐿_𝜐⁻¹𝑢𝑧,𝐘_𝑧×𝐹_𝐩𝜻(𝑧)−𝐘_𝑡𝑧+𝑅_𝑛(𝜻(𝑧−))_{2𝛽+1𝑘=1}1𝑘!_{𝐩∈𝑃𝑘}𝜕^𝐩_𝑦𝐿_𝜐⁻¹𝑢𝑧,𝐘_𝑧×𝐹_𝐩𝐗^{𝑧,𝑌𝑧}_𝑧−𝐘_𝑡𝑧+𝑅_𝑛𝐗^{𝑧,𝑌𝑧}_𝑧𝜑(𝑑𝜐)𝑑𝑧≤^𝑇_−𝑟𝜀𝐸^2𝛽+1_𝑘=11𝑘!_{𝐩∈𝑃𝑘}𝜕^𝐩_𝑦𝐿_𝜐⁻¹𝑢𝑧,𝐘_𝑧×𝐸𝐹_𝐩𝜻(𝑧)−𝐘_𝑧−𝐹_𝐩𝐗^{𝑧,𝑌𝑧}_𝑧−𝐘_𝑧∣𝒜_𝑧+𝐸𝑅_𝑛(𝜻(𝑧))∣𝒜_𝑧+𝐸𝑅_𝑛𝐗^{𝑧,𝑌𝑧}_𝑧∣𝒜_𝑧𝜑(𝑑𝜐)𝑑𝑧.(36)
Similarly, we can estimate the reminders as follows: 𝐸𝑅_𝑛(𝜻(𝑧))∣𝒜_𝑧≤𝑀1+𝐘_𝑡𝑧^2𝑟𝑧−𝑡_𝑧^𝛽+1,𝐸𝑅_𝑛𝐗^{𝑧,𝑌𝑧}_𝑧∣𝒜_𝑧≤𝑀1+𝐘_𝑡𝑧^2𝑟𝑧−𝑡_𝑧^𝛽+1.(37)
Then, by applying the Hölder inequality, Lemmas 8 and 6, inequalities (37) to estimate the inequality above, we get 𝐻₂≤𝑀^𝑇_−𝑟𝜀𝐸1+𝐘_𝑡𝑧^2𝑟Δ_𝑙^𝛽𝜑()Δ𝑑𝜐𝑑𝑧≤𝑀_𝑙^𝛽𝑇₀𝐸1+max_{0≤𝑛≤𝑛𝑇}𝐘_𝑡𝑛^2𝑟𝑧−𝑡_𝑧Δ𝑑𝑧≤𝑀_𝑙^𝛽.(38)
Step 3. Finally, by inequalities (27) and (35) as well as (38), we have 𝐸𝑔(𝑋(𝑇))−𝐸𝑔𝑌^Δ𝑙(𝑇)≤𝑀Δ_𝑙^𝛽.(39)

4. A Numerical Example

Here we give an illustrative example to demonstrate the application and the convergence of the proposed numerical scheme. We consider the following linear SDDE with Poisson jumps: 𝑑𝑋_𝑡=𝜇−𝜈𝜆𝑋_𝑡+𝛼𝑋_𝑡−1+𝑑𝑡𝜎𝑋_𝑡+𝛽𝑋_𝑡−1𝑑𝑊_𝑡+𝜈𝑋_𝑡𝑋(𝑡)[],𝑑𝑁,=𝑡+1,𝑡∈−1,0(40) where 𝜇, 𝜎, and 𝜈 are, respectively, the drift coefficient, the diffusion coefficient and the jump coefficient; 𝜆 is the jump intensity; and 𝛼 and 𝛽 are the delay coefficients.

By the method for solving linear stochastic differential equations in [23], the analytical solution for 𝑡∈[0,1] is derived as follows: 𝑋(𝑡)=Φ(𝑡)1+^𝑡₀Φ(𝑠)⁻¹+𝛼−𝜎𝛽𝑠𝑑𝑠_𝑡0Φ(𝑠)⁻¹𝛽𝑠𝑑𝑊_𝑠,(41) where Φ(𝑡)=(𝜈+1)^𝑁(𝑡)𝜎exp𝜇−𝜈𝜆−²2𝑡+𝜎𝑊_𝑡.(42)

According to the weak Taylor approximation scheme (17)–(19) proposed in Section 3, we now expand it with weak order 1 (well known as Euler scheme): 𝑌_(𝑛+1)⁻=𝑌_𝑛+𝜇−𝜈𝜆𝑌_𝑛+𝛼ℎ𝑛ℎ+𝜎𝑌_𝑛+𝛽ℎ𝑛Δ𝑊_𝑛,𝑌_𝑛+1=𝑌_(𝑛+1)⁻+𝜈𝑌_𝑛Δ𝑁_𝑛.(43) Here we have used the jump adapted time discretization, and ℎ is the maximum step size.

For higher accuracy and efficiency, one needs to construct higher order numerical schemes. We now give a Taylor scheme of weak order two below: 𝑌_(𝑛+1)⁻=𝑌_𝑛+𝜇−𝜈𝜆𝑌_𝑛ℎ++𝛼ℎ𝑛𝜎𝑌_𝑛+𝛽ℎ𝑛Δ𝑊_𝑛+𝜇−𝜈𝜆𝜎𝑌_𝑛𝑌+𝛽ℎ𝑛+𝜎𝜇−𝜈𝜆_𝑛ℎ+𝛼ℎ𝑛2Δ𝑊_𝑛+𝜇−𝜈𝜆𝜇−𝜈𝜆𝑌_𝑛ℎ+𝛼ℎ𝑛²2+𝜎𝜎𝑌_𝑛+𝛽ℎ𝑛Δ𝑊_𝑛²−ℎ2,𝑌_𝑛+1=𝑌_(𝑛+1)⁻+𝜈𝑌_𝑛Δ𝑁_𝑛.(44)

Next, we study the convergence of the two numerical schemes presented above by using the weak errors measured by 𝜀(ℎ)=|𝐸(𝑋(𝑇))−𝐸(𝑌(𝑇))|(45) and compare the results obtained from these two schemes to the explicit exact solution. We estimate the weak errors 𝜀(ℎ) by running a very large number of simulations. The exact number depends on the implemented scheme. We use the following parameters: 𝛼=0.01, 𝛽=0.01, 𝜇=0.001, 𝜎=0.6, 𝜈=0.002, and 𝜆=0.001.

In Figure 1, we give the sample paths under the two approximation schemes and the numerical explicit solution of (40). We can see from the figure that the weak Taylor scheme path is closer to the analytical solution line than the Euler scheme.

Figure 1 

Sample paths of linear SDDE with jumps.

Now we present the numerical errors generated by the two numerical schemes presented above. From Table 1, we notice that, for all the step sizes used in the numerical experiments, the weak Taylor method is more accurate. Moreover, the errors of the weak order two Taylor method decrease faster than the Euler scheme.

Table 1 

Convergence results for the linear SDDE with jumps.

5. Conclusions

In this work, we have extended previous research on weak convergence to a more general class of stochastic differential equations involving both jumps and time delay. We proved that under the Poisson random measure and a fixed time delay, a simplified Taylor method gives weak convergence rate arbitrarily close to order 𝛽. There is much scope for further work in the context of weak solution of jump-diffusion SDDEs. For example, it is clearly of great importance to extend the weak convergence theory to the case where coefficients in the equations are not globally Lipschitz, and to develop and analyze new methods that maintain good properties of convergence and stability.

Acknowledgments

The authors thank the anonymous referees and the editors for their constructive comments and suggestions. Yanli Zhou acknowledges financial support from Curtin International Postgraduate Research Scholarship (CIPRS) and Chinese Scholarship Council (CSC). Acknowledgement is also made to the China National Social Science Foundation (SSF) for the research support (10BJY104).