Simulation studies on time discrete diffusion approximations

Simulation studies on time discrete diffusion approximations

253 Mathematics and Computers in Simulation 29 (1987) 253-260 North-Holland SIMULATION STUDIES ON TIME DISCRETE DIFFVSION APPROXIMATIONS Horst LISKE...

507KB Sizes 1 Downloads 31 Views

253

Mathematics and Computers in Simulation 29 (1987) 253-260 North-Holland

SIMULATION STUDIES ON TIME DISCRETE DIFFVSION APPROXIMATIONS Horst LISKE and E&hard PLATEN Karl- Weierstrass Institute of Mathematics of the Academy of Sciences of the GDR, DDR-1086 Berlin, German Dem. Rep.

The paper presents results of simulation studies carried out with several time discrete simulation for diffusion processes with regard to the mean square and the mean convergence criterion.

algorithms

1. Introduction Let a one-dimensional diffusion process defined by the It& differential equation: dX,=a(X,)

dt+b(X,)

dF,

X = { X, }t E ,O,T1,0 < T < 00 be given, which is

TV [0, T],

(1)

denotes the standard Wiener process and X,, = x is the fixed initial value. where w= vv,qo,., In [l-4,6] and other papers time discrete approximations of diffusion processes were proposed which recursively compute approximate trajectories of the diffusion at discretization points. Such time discrete approximations yield simulation algorithms for diffusion processes in a natural way. In this paper we will test some simulation algorithms with the help of an example proposed by Talay and Pardoux in [5]. The aim of the paper is also to show some effects which occur in practical simulations.

2. Simulation algorithms In the following we present the simulation algorithms which we will apply within our simulation studies. They are based upon an equidistant discretization of the time interval [0, T] with step size A=T/n,

no {1,2, . ..}.

Furthermore, t;=iA,

iE{O,

. . . . n}=N,

denotes the ith discretization point. The simulation algorithms which we will consider are divided into two groups. The first group is derived with the help of time discrete approximations constructed according to the so-called mean square convergence criterion: 0378-4754/87/$3.50

0 1987, Elsevier Science Publishers B.V. (North-Holland)

H. Liske, E. Platen / Time discrete diffusion approximations

254

We say that a time discrete approximation Y with step size A converges for A + 0, in the mean square sense with order p E {1, 2, . . . } if a constant K exists, such that E]X,-

YT12
(2)

This criterion evaluates a pathwise approximation of the diffusion at time T. General results concerning this type of convergence can be found in [l-3,6]. The second group of simulation algorithms can be used for the approximation of the expectation of functions of diffusions. It is related to the so-called mean criterion which we shall formulate in the following simplified manner. A time discrete approximation Y with step size A converges for A + 0 in the mean sense with order p E (1, 2, . . . } if a constant K exists such that: 1 EX,-

EY,I

< KAP.

(3)

This type of convergence is studied in [3-51. We abbreviate Y= Y, and

h;=h(Y,),

for all iE {O,..., n } and functions h I R -+ R. The following two normal distributed variables will be used in the (i + l)th time step of the algorithms: AK.= q,,, - I#$, and

AZ;=

random

“+I S2dF& ds,, J t, J t,

where E AW, = 0,

E(Av.)2=A,

E (AZ,)’ = A3/3,

E AZi = 0,

E AZ, AK. = A2/2.

This pair of random variables must be generated with the help of a random number generator at the digital computer. Each of our time discrete approximations starts with Y, = x at time 0. 2. I. Mean square approximations The simplest time discrete approximation

is the Euler approximation

(EULER):

Y+i = Y + aiA + biA14$,

(4)

which shows first order mean square convergence (see [l]). In [l-3,6] we find the following second order mean square Taylor approximation Y+i = Y + aiA + b,AW, +

b,b’((

Aw)’ - A)/2.

(MS-T2): (5)

H. Liske, E. Platen / Time discrete diffusion approximations

Furthermore,

255

we will also apply the third order mean square Taylor approximation

(MS-T3):

Y1+1= r, + aiA + &Aw, + b$‘(( AW,)2 - Aj/2 +Qz;AZ, + (a;b’ + b@;‘/2)( AKA - AZ,) + (b;‘b, + ( h;)2jbiAWi(( AW’i)2/3 - A j/2 + ( aia; + b;a:‘/2)

A2,‘2,

(6)

which is also proposed in [l-3,6]. For the inplementation of simulation algorithms on a computer it is desirable to avoid derivatives of the drift and diffusion coefficients. The following mean square second order derivative-free approximation (MS-DF) is proposed in [3]: yI+r = Y + aiA + biAK + (b( Y + biA’j2) - bi) AP112(( AK)’ - A)/2.

(7)

2.2. Mean approximations [3] and [4] show that the Euler approximation can be interpreted as the first order mean Taylor approximation. In [3] one also finds the following second order mean Taylor approximation (M-T2) : y+l

= x +

aiA + biAw;. f b,b;(( AK.)’ - A)/2

+b,u;AZ, + (a$; + b;b”/2)(AwA

- AZ,) + (a&

Finally, we give the second order derivative-free mean approximation

+ bfa92)

A2,‘2.

(8)

(M-DF) proposed in [3]:

~+, = ~ + (~( r, + aiA + b,AW,) + ai) A/2 +(b(Y+)+b(Y-)+2bi)

AK/4+(b(Y+)-b(q-))((A&)‘-A)

A-“2/4, (9

with Y* I = YI + a_A ,--I + b.A’12. In comparison with (8) one notes that here the random variable AZ, does not appear.

3. Example To get an impression about the numerical behaviour of different simulation algorithms we test them with the help of the following example. For this purpose it is important that the exact value of the realizations of the diffusion process at time T can be computed in our example so that one can also determine the mean square error of the approximations. The following stochastic process is given in [5] as an example for a nonlinear diffusion process: dX,=

-(sin2X,+$sin4X,)

dt+fi(cos

X,)‘dw,

tE [0, T], X,=1.

(IO)

H. Liske, E. Platen / Time discrete diffusion approximations

256

By the Ita formula it follows that we shall obtain the same process if we take X, = arctan( I$)

(11)

with v= V,

V,e-f+fi =

/0

fe”-LdW,

(12)

tan(l).

(13)

Improving the results in [5] we are able to generate the exact value of the realizations of the diffusion process at T and we propose the following method of computation. It holds that VT=e-r

b+Jz

i

e’l+l( Ay_l

+ ACJ_l + AQi_l)

(14)

i=l

with A~_, = J” (S - ti-1) dW, ?,-I

(15)

and

By the use of the Ita formula one can show: A&

= Aw_,A

- A&_,,

and we have the representation VT=ePT

b+fi (

i

07)

of V,: e’g-l( Ae_,(l

+ A) + AZ,_i + AQ,_i)

i=l

For each in {l,..., n} the random vector (Aq, mean and correlations: E(AK)‘=

A,

.

08)

i

AZi, AQi) is normally distributed with zero

E( AZ,)’ = A3/3,

E( AQi)2 = (e2* - 1)/2 - 2A e* + A + A2+ A3/3, E AQiAF = e* - (1 + A + A2/2),

E Aw.AZ, = A2/2,

E AQiAZi = e* - (1 + A + A2/2 + A3/6). Generating Aq, AZ, and AQ, at each time step we can evaluate the exact increment of X as well as the increments ot the time discrete approximations given in (4)-(9).

4. Simulation studies By the use of the FORTRAN-program SIRE (see [7]) for the simulation of Ita random equations the following simulation studies were made at a digital computer SM 4-20: 100 blocks

H. Liske, E. Platen / Time discrete diffusion approximations

251

with 100 trajectories were generated for each of the above simulation algorithms with three different step sizes. The simulated trajectories were arranged in blocks with the aim of also obtaining an impression about the dispersion of the simulated estimators. We separately present our results with respect to the mean square criterion and the mean criterion. The parameters A and T are chosen in such a way that one can observe some interesting effects. 4. I. Mean square criterion Let YTkj denote the value of the k th simulated approximate trajectory at time T of the jth loo}. XTkj denotes the value of the k th trajectory of the diffusion X at block, k, Jo {l,..., time T of the jth block. For all j E {1,. . . , lOO} we compute the value: 100 sj=dO

C k=l

(xTkj-

yZ%j)2Y

which estimates the mean square error between X, and Y, in the jth block. Because of the central limit theorem we know that Sj is nearly normally distributed estimate its mean by:

(19)

and we

100

(20)

S=&&S, j=l

and its dispersion by: 100 p=

&

c

(s,-

s)2.

j=l

(21)

To illustrate how the values S,, j E {1,. . . , lOO} are spread around the average S of the mean square error in Table 1 we will give the 0.9 confidence interval of Sj based upon the Student-t distribution in the form: (s - 50.166;

s + S.O.166).

The simulations with respect to the mean square criterion are made for the Euler approximation (EULER), the second-order mean square Taylor approximation (MS-T2), the second-order derivative free approximation (MS-DF) and the third-order mean square Taylor approximation (MS-T3). We have chosen the end of the time interval T = 0.5 and the three different step sizes = 2-4,2-5,2-4. In Table 1 the last column contains the ratio CT between the necessary computer time for the simulation of an approximate trajectory and the corresponding time for the Euler approximation with step size A = 2-4. Let us discuss the results of Table 1.

H. Liske, E. Platen / Time discrete diffusion approximations

258 Table 1 A

Algorithm

s - S.O.166

S

S + S.O.166

2-4

EULER MS-T2 MS-DF MS-T3

1.50x10-4 1.61 x 1O-4 1.43 x 10-4 2.76 x lo-’

1.79x1o-4 1.81 x 1O-4 1.62~10-~ 3.19x10-5

2.08 x 1O-4 2.01 x 10-4 1.79 x 10-4 3.62x10V5

1.00 1.40 1.64 4.58

2r5

EULER MS-T2 MS-DF MS-T3

5.17 x10-5 3.15 x10-5 3.45 x 10-5 6.35 x 1O-6

6.30x10-’ 3.46 x 10V5 3.82x10W5 7.35 x 10-6

7.43 x 10-5 3.76~10-~ 4.19x10-5 8.35 x 1O-6

2.00 2.80 3.28 9.16

2-6

EULER MS-T2 MS-DF MS-T3

2.43 x 8.36 x 7.78 x 1.16 x

2.88 x 9.15 x 8.54x 1.33 x

3.33 x 9.95 x 9.28 x 1.50x

lo-’ lO-‘j 1O-6 1O-6

lo-’ 10-6 1O-6 1o-6

1or5 10-6 1O-6 1or6

CT

4.00 5.60 6.56 18.32

For the step size A = 2-4 the first three algorithms show nearly the same relatively large error. Only the third order mean square Taylor approximation is considerably better but it needs nearly 5 times the computer time of the Euler approximation. For the step size A = 2-5 all algorithms improved their accuracy, where the second order mean square Taylor approximation and the second-order mean square derivative-free approximation show approximately the same results. For the step size A = 2-6, again all the algorithms show better results, where the third order scheme is better than the second-order schemes which have similar accuracy, and the second-order schemes are better than the Euler scheme. One observes that the length of the confidence intervals, which depends on the dispersion of the mean square error, becomes smaller by smaller step sizes A and higher order algorithms. To decide which algorithm is really better for a given problem depends on the desired accuracy and the necessary computer time. For instance, if one is interested in a mean square error close to 4 X lo-‘, then the second-order mean square Taylor approximation and the second order derivative-free approximation with step size A = 2-5 are good, because for such accuracy they need much less computer time than, for instance, the Euler scheme with A = 2-6 or the third-order mean square Taylor approximation with step size A = 2-4. In general in the above case one would prefer the second-order derivative-free approximation, because the scheme avoids derivatives of the drift and diffusion coefficients, which allows an easy implementation on the computer. The above simulation study shows that for the construction of efficient algorithms one has to choose the step size and the order of the algorithm depending on the desired accuracy. 4.2. Mean criterion To study the behaviour of the different mean approximations with respect to the mean criterion we consider the same example as above, which is given in (lo), on the time interval of

H. Liske, E. Platen / Time discrete diffusion approximations

length T= 16. Looking at the mean criterion (3) we study the approximation

259

of the functional:

EX,= 1.1 x 10-7,

(22)

that means we approximate the first moment of the diffusion X at time T. Let again YTki denote the value of the k th simulated approximate trajectory at time T of the j th block, k, j E { 1,. . . , loo}. For all j E {1,. ..,lOO}the value: 100

Rj =

c

&

Yrkj - EX,,

(23)

k=l

estimates the mean error (see (3)) of YT in the j th block of our example. Because of the central limit theorem we know that Rj is nearly normally distributed, where we can estimate its mean by:

l;

R=&

(24)

Rj

j=l

and its dispersion by: 100

C (~~4)~.

R*=&

(25)

j=l

Similarly as for the mean square criterion we illustrate in Table 2 how the values loo}, are spread around R by the 0.9 confidence interval:

)R, 1,

jG {l,...,

(IRI-R-0.166;

IR(+%0.166).

Within the simulation study, with respect to the mean criterion, we apply the Euler mation (EULER), the second-order derivative-free mean approximation (M-DF) second-order mean Taylor approximation (M-T2). We consider the three step sizes A = 2*, 2’, 2-l. In Table 2 the last column gives analogously as in Table 1 the ratio C between computer time for the simulation of an approximate trajectory and the corresponding the Euler approximation with step size A = 2.

Table 2 A

Algorithm

[RI-R.0.166

IRI

1R ) + i? .0.166

2l

EULER M-T2 M-DF

0.276 2.175 0.358

0.345 2.452 0.435

0.414 2.729 0.512

1.00 3.36 3.42

20

EULER M-T2 M-DF

0.047 0.226 0.014

0.111 0.309 0.051

0.175 0.392 0.088

2.00 6.72 6.84

2-l

EULER M-T2 M-DF

0.018 0.124 0.035

0.059 0.159 0.063

0.101 0.194 0.090

4.00 13.44 13.68

C

approxiand the

the used time for

H. Liske, E. Platen / Time discrete diffusion approximations

260

When, in the following, we discuss the above simulation results, we have to notice that we should not compare this table with Table 1, because the studied errors and the parameters are completely different. We note that the accuracy in general improves for smaller step sizes. Only the derivative-free scheme shows no improvement for the smallest step size A = 0.5. In this case the step size A = 1 was already adapted. The necessary computer times for the second-order mean Taylor approximation and the second-order derivative-free mean approximation are approximately 3.4 times that of the Euler approximation. But the Taylor approximation behaves better than the derivative-free approximation especially for larger step sizes (A = 2, A = 1). That means the Taylor approximation shows some instable behaviour for large step sizes. In any case, the derivative-free scheme can be implemented in an easier way than the corresponding Taylor approximation and shows better results in our example.

5. Conclusions In the above very special simulation studies it becomes clear that one can only construct efficient simulation methods by considering their dependence on the desired convergence criterion and that their accuracy is governed by the choice of an appropriate algorithm and an adapted step size. The numerical stability of simulation schemes is an important problem which needs further theoretical investigation. Other open problems are connected with the step size control and the choice of the sample size.

References [l] G.N. Milstein, Approximate

intergration

of stochastic differential

equations,

Theor. Prob. Appl. XIX (3) (1974)

583-588. [2] E. Platen, An approximation method for a class of Iti? processes, Lietuuos Matem. Rink. XXI (1) (1981) 121-133. [3] E. Platen, BeitrPge zur zeitdiskreten Approximation von Itaprozessen, Diss. B, AdW der DDR, Berlin, 1984. [4] D. Talay, Efficient numerical schemes for the approximation of expectations of functionals of the solution of a S.D.E. and applications, Proc. Colloquium ENST-CNET 1983, Lecture Notes in Control and Information Science (Springer, Berlin, 1984). [5] D. Talay and E. Pardoux, Discretization and simulation of stochastic differential equations, Publications de Mathematiques Appliquees. Marseille-Toulon, 1983-5. [6] W. Wagner and E. Platen, Approximation of ItS integral equations, Preprint ZIMM, AdW der DDR, Berlin, 1978. [7] FORTRAN-program “SIRE” for the simulation of Ita random equations, AdW der DDR, IMATH, Department of Probability, DDR-1086 Berlin, German Dem. Rep.