Probabilistic numerical approximation for stochastic

1 downloads 0 Views 387KB Size Report
The author is grateful to Nizar Touzi, J. Frédéric Bonnans and Nicole El. Karoui for their comments ...... PR for every ν ∈ U. In resume, we have PS0 ⊂ PS ⊂ PR.
Probabilistic numerical approximation for stochastic control problems Xiaolu TAN



August 7, 2012

Abstract We give a probabilistic interpretation of the Monte Carlo scheme proposed by Fahim, Touzi and Warin [Ann. Appl. Probab. 21(4) : 1322-1364 (2011)] for fully nonlinear parabolic PDEs, and hence generalize it to the non-Markovian case for a general stochastic control problem. General convergence result is obtained by weak convergence method in spirit of Kushner and Dupuis [18]. We also get a rate of convergence using the invariance principle technique as in Dolinsky [7], which is better than that obtained by viscosity solution method. Finally, by approximating the conditional expectations arising in the numerical scheme with simulation-regression method, we obtain an implementable scheme. Key words. Numerical scheme, non-Markovian stochastic control, weak convergence, invariance principle. MSC 2010. Primary 65K99, secondary 93E20, 93E25

1

Introduction

Stochastic optimal control theory is largely applied in economics, finance, physics and management problems. Since its development, numerical methods for stochastic control problems have also been largely investigated. For the Markovian control problem, the value function can usually be characterized by Hamilton-Jacob-Bellman(HJB) equations, then many numerical methods are also given as numerical schemes for PDEs. In this context, a powerful tool to prove the convergence is the monotone convergence of viscosity solution method of Barles and Souganidis [1]. ∗

CMAP, Ecole Polytechnique, Paris, [email protected]. Research supported by the Chair Financial Risks of the Risk Foundation sponsored by Soci´et´e G´en´erale, the Chair Derivatives of the Future sponsored by the F´ed´eration Bancaire Fran¸caise, and the Chair Finance and Sustainable Development sponsored by EDF and CA-CIB. The author is grateful to Nizar Touzi, J. Fr´ed´eric Bonnans and Nicole El Karoui for their comments and suggestions.

1

In the one-dimensional case, the explicit finite difference scheme can be easily constructed and implemented, and the monotonicity is generally guaranteed under the Courant-Friedrichs-Lewy (CFL) condition. In two dimensional cases, Bonnans, Ottenwaelter and Zidani [3] proposed a numerical algorithm to construct monotone explicit schemes. Debrabant and Jakobsen [6] gave a semi-Lagrangian scheme which is easily constructed to be monotone but needs finite difference grid together with interpolation method for the implementation. In general, these methods may be relatively efficient in low dimensional cases; while in high dimensional cases, Monte Carlo method is preferred if possible. As a generalization of Feynman-Kac formula, the Backward Stochastic Differential Equation (BSDE) opens a way for the Monte Carlo method for optimal control problems (see e.g. Bouchard and Touzi [4], Zhang [27]). General speaking, the BSDE covers the controlled diffusion process problems of which only the drift part is controlled. However, it cannot include the general control problems when the volatility part is also controlled. This is one of the main motivations of recent developments of second order BSDE (2BSDE) by Cheridito, Soner, Touzi and Victoir [5] and Soner, Touzi and Zhang [21]. Motivated by the 2BSDE theory in [5], and also inspired by the numerical scheme of BSDE, Fahim, Touzi and Warin [11] proposed a probabilistic numerical scheme for fully nonlinear parabolic PDEs. In their scheme, one needs to simulate a diffusion process, and estimate the value function as well as the derivatives of the value function arising in the PDE by conditional expectations, and then compute the value function in a backward way on the discrete time grid. The efficiency of this Monte Carlo scheme has been shown by several numerical examples, we refer to Fahim et al. [11], Guyon and Henry-Labord`ere [14] and Tan [24] for the implemented examples. However, instead of probabilistic arguments, the convergence of this scheme is proved by techniques of monotone convergence of viscosity solution of Barles and Souganidis [1]. Moreover, their scheme can be only applied in the Markovian case when the value function is characterized by PDEs. The main contribution of this paper is first to give a probabilistic interpretation to the Monte Carlo scheme of [11] for fully nonlinear PDEs, which permits to generalize it to the non-Markovian case for a general stochastic optimal control problem. One of the motivations for the non-Markovian generalization comes from finance to price the path-dependent exotic derivative options in the uncertain volatility model. Our general convergence result is obtained by weak convergence techniques in spirit of Kushner and Dupuis [18]. In contrast to [18], where the authors define their controlled Markov chain in a descriptive way, we give our controlled discrete time semimartingale in an explicit way using the cumulative distribution functions. Moreover, we introduce a canonical space for the control problem following El Karoui et al. [9], which permits to explore the convergence conditions on the reward functions. We also provide a convergence rate using the invariance principle techniques of Sakhanenko [22] and Dolinsky [7]. Comparing to the scheme of [7] in the context of G−expectation (see e.g. Peng [19] for G−expectation), our scheme is implementable using simulationregression method.

2

The rest of the paper is organized as follows. In Section 2, we first introduce a general non-Markovian stochastic control problem, and propose a numerical scheme. Then we give the assumptions on the diffusion coefficients and the reward functions, as well as the main convergence results, including the general convergence and a rate of convergence. Next in Section 3, we provide a probabilistic interpretation of the numerical scheme, by showing that the numerical solution is equivalent to the value function of a controlled discrete time semimartingale problem. Then we complete the proofs of the convergence results in Section 4. Finally, in Section 5, we discuss some issues about the implementation of the numerical scheme, including the simulationregression method. Notations. We denote by Sd the space of all d × d matrix, and by Sd+ the space of all positive symmetric d × d matrix. Let Ωd := C([0, T ], Rd ) be the space of all continuous paths between 0 and T , denote |x| := sup0≤t≤T |xt | for every x ∈ Ωd . In the paper, E is a fixed compact Polish space, we denote QT

:= [0, T ] × Ωd × E.

Suppose that (Xtk )0≤k≤n is a process defined on the discrete time grid (tk )0≤k≤n of b [0, T ] with tk := kh and h := Tn , we usually write it as (Xk )0≤k≤n , and denote by X its linear interpolation path on [0, T ]. In the paper, C is a constant whose value may vary from line to line.

2 A numerical scheme for stochastic control problems 2.1

A non-Markovian stochastic control problem

Let (Ω, F, P) be a complete probability space containing a d−dimensional standard Brownian motion W , F = (Ft )0≤t≤T be the natural Brownian filtration. Denote by Ωd := C([0, T ], Rd ) the space of all continuous paths between 0 and T . Suppose that E is a compact Polish space with a complete metric dE , (σ, µ) are bounded continuous functions defined on QT := [0, T ] × Ωd × E taking value in Sd × Rd . We fix a constant x0 ∈ Rd through out the paper. Then given a F−progressive E−valued process ν = (νt )0≤t≤T , denote by X ν the controlled diffusion process which is the strong solution to Z t Z t Xtν = x0 + µ(s, X·ν , νs )ds + σ(s, X·ν , νs )dWs . (2.1) 0

0

To ensure the existence and uniqueness of the strong solution to the above equation (2.1), we suppose that for every progressive process (X, ν), the processes µ(t, X· , νt ) and σ(t, X· , νt ) are progressively measurable. In particular, µ and σ depend on the past trajectory of X. Further, we suppose that there is a constant C and a continuity

3

module ρ, which is an increasing function on R+ satisfying ρ(0+ ) = 0, such that µ(t1 , x1 , u1 ) − µ(t2 , x2 , u2 ) + σ(t1 , x1 , u1 ) − σ(t2 , x2 , u2 ) ≤

C|xt11 − xt22 | + ρ(|t1 − t2 | + dE (u1 , u2 )),

(2.2)

where for every (t, x) ∈ [0, T ] × Ωd , we denote xts := xs 1[0,t] (s) + xt 1(t,T ] (s). Let Φ : x ∈ Ωd → R and L : (t, x, u) ∈ QT → R be the continuous reward functions, and denote by U the collection of all E−valued F−progressive processes, the main purpose of this paper is to approximate numerically the following optimization problem: V

:= sup E ν∈U

hZ

T

i  L t, X·ν , νt dt + Φ(X·ν ) .

(2.3)

0

Similarly to µ and σ, we suppose that for every progressive process (X, ν), the process t 7→ L(t, X· , νt ) is progressively measurable. Moreover, to ensure that the expectation in (2.3) is well defined, we shall assume later that L and Φ are of exponential growth in x and discuss their integrability in Proposition 2.5.

2.2

The numerical scheme

In preparation of the numerical scheme, we shall fix, through out the paper, a progressive function σ0 : [0, T ] × Ωd → Sd such that σ0 (t1 , x1 ) − σ0 (t2 , x2 ) ≤ C|xt1 − xt2 | + ρ(|t1 − t2 |)), 1 2 ∀(t1 , x1 ), (t2 , x2 ) ∈ [0, T ] × Ωd , and with ε0 > 0, σ0 σ0T (t, x) ≥ ε0 Id for every (t, x) ∈ [0, T ] × Ωd . Denote σ0t,x := σ0 (t, x),

at,x := σ0t,x (σ0t,x )T , 0

at,x := σσ T (t, x, u) − at,x u 0 ,

bt,x := µ(t, x, u). u

(2.4)

Then we define a function G on [0, T ] × Ωd × Sd × Rd by G(t, x, γ, p) := sup u∈E

 1 t,x L(t, x, u) + at,x u · γ + bu · p , 2

(2.5)

which is clearly convex in (γ, p) as the supremum of a family of linear functions, and lower-semicontinuous in (t, x) as the supremum of a family of continuous functions. Let n ∈ N, denote the time discretization by h := Tn and tk := hk. Let us take the standard d−dimensional Brownian motion W in the complete probability space (Ω, F, P). For simplicity, we denote Wk := Wtk , ∆Wk := Wk −Wk−1 , W 0 FkW := σ(W0 , W1 , · · · , Wk ) and EW k [·] := E[·|Fk ]. Then we have a process X on the discrete grid (tk )0≤k≤n defined by X00 := x0 ,

0 c0 · )∆Wk+1 , Xk+1 := Xk0 + σ0 (tk , X

(2.6)

c0 denotes the linear interpolation process of (X 0 )0≤k≤n on interval [0, T ]. where X k

4

Then, our numerical scheme is given by c0 Yk := EW k [Yk+1 ] + h G(tk , X · , Γk , Zk ),

(2.7)

with the terminal condition c0 · ), Yn := Φ(X

(2.8)

where G is defined by (2.5) and T i h i h T −1 ∆Wk+1 ∆Wk+1 − hId −1 W T −1 ∆Wk+1 σ , Γk := EW Y (σ ) , Z := E Y (σ ) k+1 0,k k k+1 0,k k k 0,k h2 h c0 c0 · ). with σ0,k := σ0tk ,X = σ0 (tk , X To emphasize the dependance of Y on the time discretization h, we shall also denote by Y h the solution of the above numerical scheme.

Remark 2.1. By its definition, Yk is a measurable function of (X00 , · · · , Xk0 ). We shall show later in Proposition 2.5 that the function Yk (x0 , · · · , xk ) is of exponential growth in max0≤i≤k |xi | under appropriate conditions, and hence the conditional expectations in (2.7) are well defined. Therefore, the above scheme (2.7) should be well defined. Remark 2.2. In the Markovian case, when the function Φ(·) (resp. L(t, ·, u), µ(t, ·, u) and σ(t, ·, u)) only depends on XT (resp. Xt ), so that the function G(t, ·, γ, z) only depends on Xt and the value function of the optimization problem (2.3) can be characterized as the viscosity solution of a nonlinear PDE −∂t v −

1 a0 (t, x)D2 v − G(t, x, D2 v, Dv) = 0. 2

Then the above scheme reduces to that proposed by Fahim, Touzi and Warin [11].

2.3

The convergence results of the scheme

Our main idea to prove the convergence of the scheme (2.7), (2.8) is to interpret it as an optimization problem on a system of controlled discrete time semimartingales, which converge weakly to the controlled diffusion processes. Therefore, a reasonable assumption is that Φ and L are bounded continuous on Ωd (i.e. Φ(·), L(t, ·, u) ∈ Cb (Ωd )), or they belong to the completion space of Cb (Ωd ) under an appropriate norm. We shall suppose that Φ (resp. L) is continuous in x (resp. (t, x, u)), and there are a constant C and continuity modules ρ0 , (ρN )N ≥1 such that for every (t, x, u) ∈ QT and (t1 , x1 ), (t2 , x2 ) ∈ [0, T ] × Ωd and N ≥ 1,    |Φ(x)| + |L(t, x, u)| ≤ C exp C|x| ,   (2.9) |L(t1 , x, u) − L(t2 , x, u)| ≤ ρ0 (|t1 − t2 |),    |Φ (x ) − Φ (x )| + |L (t, x , u) − L (t, x , u)| ≤ ρ (|x − x |), N

1

N

2

N

1

N

2

where ΦN := (−N ) ∨ (Φ ∧ N ) and LN := (−N ) ∨ (L ∧ N ).

5

N

1

2

Denote mG :=

 1 T t,x w au w + bt,x u ·w , 2

min

(t,x,u)∈QT , w∈Rd

(2.10)

and h0 := 1mG =0 T + 1mG −∞ and h0 > 0. at,x u

≥ 0

and

1 −

≥ 0.

Remark 2.3. Assumption 1 is almost equivalent to Assumption F of Fahim et al. [11] in the context of the control problem, and it implies that the drift µ and σ are uniformly bounded, as assumed at the beginning of Section 2.1. In particular, it follows that when mG < 0, we have 1 t,x −1 + hmG ≥ 0, for every (t, x, u) ∈ QT and h ≤ h0 . 1 − at,x u · a0 2 Remark 2.4. When at,x u ≥ εId uniformly for some ε > 0, we get immediately mG > t,x −∞ since bu is uniformly bounded. When at,x u degenerates, mG > −∞ implies that t,x t,x bu lies in the image of au . Proposition 2.5. Suppose that the reward functions L and Φ satisfy (2.9), then the optimal value V in (2.3) is finite. Suppose in addition that Assumption 1 holds true. Then for every fixed n ∈ N (h := Tn ) and every 0 ≤ k ≤ n, as a function of (X0 , · · · , Xk ), Ykh (x0 , · · · , xk ) is also of exponential growth in max0≤i≤k |xi |. And hence Ykh is integrable in (2.7), the numerical scheme (2.7) is well defined. The proof is postponed in Section 3.1 after a technical lemma. Our main results of the paper are the following two convergence theorems, whose proofs are left in Section 4. Theorem 2.6. Suppose that L and Φ satisfy (2.9) and Assumption 1 holds true. Then Y0h → V

as h → 0.

To derive a convergence rate, we suppose further that E is a compact convex subset of Sd+ × Rd , and for every (t, x, u) = (t, x, a, b) ∈ QT , a > 0,

µ(t, x, u) = µ(t, x, a, b) = b,

σ(t, x, u) = σ(t, x, a, b) = a1/2 .

(2.12)

Moreover, we suppose that L(t, x, u) = `(t, x) · u for some continuous function ` : [0, T ] × Ωd → Sd × Rd and that there exists a constant C > 0 such that for every couple (t1 , x1 ), (t2 , x2 ) ∈ QT , |`(t1 , x1 ) − `(t2 , x2 )| + |Φ(x1 ) − Φ(x2 )|   ≤ C |t1 − t2 | + |xt11 − xt22 | + |x1 − x2 | exp C(|x1 | + |x2 | .

6

(2.13)

Theorem 2.7. Suppose that L and Φ satisfy conditions (2.9) and (2.13), the set E ⊂ Sd+ ×Rd is compact and convex, functions µ and σ satisfy (2.12), and Assumption 1 holds true. Then for every ε > 0, there is a constant Cε such that |Y0h − V | ≤ Cε h1/8−ε ,

∀h ≤ h0 .

(2.14)

If, in addition, L and Φ are bounded, then there is a constant C such that |Y0h − V | ≤ Ch1/8 ,

∀h ≤ h0 .

(2.15)

Remark 2.8. In the Markovian context as in Remark 2.2, Fahim, Touzi and Warin [11] obtained a global convergence rate h1/10 using Krylov’s shaking coefficient method. We get a rate h1/8 in the non-Markovian case under some additional constraints. When there is no control on the volatility part, the BSDE method in Bouchard and Touzi [4] and Zhang [27] gives a convergence rate of order h1/2 . Our method cannot achieve this rate in the BSDE context. Remark 2.9. When the covariance matrix σσ T is diagonal dominated, Kushner and Dupuis [18] gave a systematic way to construct a convergent finite difference scheme. However, the construction turns to be not easy when the matrix is not diagonal dominated, see e.g. Bonnans et al. [3]. Our scheme relaxes this constraint. Moreover, our scheme implies a natural Monte Carlo implementation, which may be more efficient in high dimensional cases, see numerical examples in Fahim et al. [11], Guyon et al. [14] and Tan [24].

3 A controlled discrete semimartingale interpretation Before giving the proofs of the above convergence theorems, we first provide a probabilistic interpretation to the scheme (2.7), (2.8) in spirit of Kushner and Dupuis [18]. Namely, we shall show that the numerical solution is equivalent to the value function of a controlled discrete time semimartingale problem. For finite difference schemes, the controlled Markov chain intrepration given by Kushner and Dupuis [18] is straightforward, where their construction of the Markov chain is descriptive. For our scheme, the probabilistic interpretaion is less evident as the space state is uncountable. Our main idea is to use the inverse function of the distribution functions. This question has not been evoked in the Markovian context of Fahim et al. [11] since they use the monotone convergence of viscosity solution technique, where the idea is to show that the terms Z and Γ defined below (2.8) are good approximations of the derivatives of the value function.

3.1

A technical lemma

Given a fixed (t, x, u) ∈ QT , let us simplify further the notations in (2.4), σ0 := σ0t,x , a0 := at,x 0 ,

t,x au := at,x u , bu := bu .

7

(3.1)

Denote fh (t, x, u, x) :=

1 −1 T −1  exp − h x a0 x 2 (2πh)d/2 |σ0 |1/2   1 −1 1 −1 −1 T T −1 + b · a x + h a · a xx (a ) . (3.2) 1 − au · a−1 u u 0 0 0 0 2 2 1

It follows by (2.10) that for every (t, x, u) ∈ QT and x ∈ Rd , bu · a−1 0 x +

h 1 −1 x 1 x xT T −1 i T T −1 h au · a−1 = h bu · a−1 + au · a−1 (a ) 0 xx (a0 ) 0 0 2 h 2 h h 0 ≥ h mG .

Then under Assumption 1, one can verify easily (see also Remark 2.3) that when h ≤ h0 for h0 given by (2.11), x 7→ fh (t, x, u, x) is a probability density function on Rd , i.e. Z fh (t, x, u, x) ≥ 0, ∀x ∈ Rd , and fh (t, x, u, x) dx = 1. Rd

Lemma 3.1. Let h ≤ h0 and R be a random variable with probability density x 7→ fh (t, x, u, x). Then for all functions g : Rd → R of exponential growth, we have h  Wh WhT − hI −1 i Wh 1 ,(3.3) + hau · (σ0T )−1 σ0 E[g(R)] = E g(σ0 Wh ) 1 + hbu · (σ0T )−1 h 2 h2 where Wh is a d−dimensional Gaussian random variable with distribution N (0, hId ). In particular, it follows that there exists a constant C1 independent of (h, t, x, u) ∈ (0, h0 ] × QT such that E[R] = bu h,

Var[R] = (au + a0 )h − bu bTu h2

and E[|R|3 ] < C1 h3/2 .

(3.4)

Moreover, for any c ∈ Rd , E[ec·R ] ≤ eC2 h (1 + C2 h),

(3.5)

where C2 is independent of (h, t, x, u) and is defined by C2 :=

sup (t,x,u)∈QT

1 2

cT a0 c + |bu · c| +

 1 T c au c . 2

 Proof. First, it is clear that (2πh)d/21 |σ |1/2 exp − 21 h−1 xT a−1 x is the density function 0 0 of σ0 Wh . Then by (3.2), Z E[g(R)] = fh (t, x, u, x) g(x) dx d R Z 1 1 −1 T −1  = h x a0 x exp − d/2 |σ |1/2 2 0 Rd (2πh)   1 1 −1 −1 −1 T T −1 g(x) 1 − au · a−1 h a · a + b · a x + x x(a ) u u 0 0 0 0 2 2 h  T Wh Wh − hI −1 i Wh 1 = E g(σ0 Wh ) 1 + hbu · (σ0T )−1 + hau · (σ0T )−1 σ0 . h 2 h2

8

Hence (3.3) holds true. In particular, let g(R) = R or g(R) = RRT , it follows by direct computation that the first two equalities of (3.4) hold true. Further, letting g(x) = |x|3 , we get from (3.3) that h  i √ 1 E[|R|3 ] = h3/2 E |σ0 N |3 1 + hbu · (σ0T )−1 N + au · (σ0T )−1 (N N T − I)σ0−1 , 2 where N is a Gaussian vector with distribution N (0, Id ). And hence (3.4) holds true with h  i p  1 C1 := sup E |σ0 N |3 1 + h0 bu · (σ0T )−1 |N | + au · (σ0T )−1 N N T − I σ0−1 , 2 (t,x,u)∈QT which is clearly bounded and independent of (h, t, x, u). √ Finally, to prove inequality (3.5), we denote Nh := N + hσ0 c for every h ≤ h0 . Then h T  Wh WhT − hId −1 i Wh 1 E[ec·R ] = E ec σ0 Wh 1 + hbu · (σ0T )−1 + hau · (σ0T )−1 σ0 h 2 h2 i h T √  √  1 = E ec σ0 N h 1 + hbu · (σ0T )−1 N + au · (σ0T )−1 N N T − Id σ0−1 2 i h √  cT a0 c 1 = e 2 h E 1 + hbu · (σ0T )−1 Nh + au · (σ0T )−1 Nh NhT − Id σ0−1 2  cT a0 c 1 = e 2 h 1 + (bu · c + cT au c)h 2 ≤ eC2 h (1 + C2 h),   where C2 := sup(t,x,u)∈QT 21 cT a0 c + |bu · c| + 21 cT au c is bounded and independent of (h, t, x, u). Remark 3.2. Since the random vector R does not degenerate to the Dirac mass, it follows by (3.4) that under Assumption 1, σσ T (t, x, u) > µµT (t, x, u)h,

for every (t, x, u) ∈ QT , h ≤ h0 .

With this technical lemma, we can give the Proof of Proposition 2.5. For the first assertion, it is enough to prove that supν∈U E[exp(C|X·ν |)] is bounded by condition (2.9). Note that µ and σ are uniformly bounded. When d = 1, X ν is a continuous semimartingale whose finite variation and quadratic variation are both bounded by a constant RT for every ν ∈ U. It follows by Dambis-Dubins-Schwartz’s time change theorem that sup E[exp(C|X·ν |)] ≤ eCRT E exp C sup |Bt | ν∈U



< ∞,

0≤t≤RT

where B is a standard one-dimensional Brownian motion. When d > 1, it is enough   to remark that for X = (X 1 , · · · , X d ), exp C|X· | ≤ exp C(|X·1 | + · · · + |X·d |) ; and we then conclude the proof of the first assertion applying Cauchy-Schwartz inequality.

9

For the second assertion, let us prove it by backward induction. Given 0 ≤ k ≤ n−1, x0 , · · · , xk ∈ Rd , we denote by x b the linear interpolation path of x(ti ) := xi , denote also Lk (x0 , · · · , xk , u) := L(tk , x b· , u),

∀u ∈ E.

(3.6)

For the terminal condition, it is clear that Yn (x0 , · · · , xn ) is of exponential growth in max0≤i≤n |xn | by condition (2.9). Now, suppose that  Yk+1 (x0 , · · · , xk+1 ) ≤ Ck+1 exp Ck+1 max |xi | . 0≤i≤k+1

Let Ru be a random variable of distribution density x 7→ fh (tk , x b, u, x). Then it follows by (2.7) and Lemma 3.1 that Yk (x0 , · · · , xk ) n h  = sup hLk (x0 , · · · , xk , u) + E Yk+1 x0 , · · · , xk , xk + σ0 Wh u∈E

Wh WhT − hI −1 io Wh 1 tk ,bx 1 + hbutk ,bx · (σ0T )−1 + hau · (σ0T )−1 σ0 h 2 h2 n h io  = sup hLk (x0 , · · · , xk , u) + E Yk+1 x0 , · · · , xk , xk + Ru . (3.7) 

u∈E

Therefore by (2.9) and (3.5), Yk (x0 , · · · , xk )    ≤ (Ck+1 + Ch) exp (Ck+1 + Ch) max |xi | sup E exp Ck+1 |Ru | 0≤i≤k u∈E  C2 h ≤ e (1 + C2 h)(Ck+1 + Ch) exp (Ck+1 + Ch) max |xi | , 0≤i≤k

where C is the same constant given in (2.9), and the constant C2 is from (3.5) depending on Ck+1 . We then conclude the proof.

3.2

The probabilistic interpretation

In this section, we shall interpret the numerical scheme (2.7) as the value function of a controlled discrete time semimartingale problem. In preparation, let us show how to construct the random variables with density function x 7→ fh (t, x, u, x). Let F : R → [0, 1] be the cumulative distribution function of a one-dimensional random variable, denote by F −1 : [0, 1] → R its generalized inverse function. Then given a random variable U of uniform distribution U ([0, 1]), it is clear that F −1 (U ) turns to be a random variable with distribution F . In the multi-dimensional case, we can convert the problem to the one-dimensional case since Rd is isomorphic to [0, 1], i.e. there is a one-to-one mapping κ : Rd → [0, 1] such that κ and κ−1 are both Borel measurable (see e.g. Proposition 7.16 and Corollary 7.16.1 of Bertsekas and Shreve [2]). Define Z Fh (t, x, u, x) := fh (t, x, u, y)κ(y)dy. κ(y)≤x

10

It is clear that x 7→ Fh (t, x, u, x) is the distribution function of random variable κ(R) where R is a random variable of density function x 7→ fh (t, x, u, x). Denote by Fh−1 (t, x, u, x) the inverse function of x 7→ Fh (t, x, u, x) and Hh (t, x, u, x) := κ−1 (Fh−1 (t, x, u, x)).

(3.8)

Then given a random variable U of uniform distribution on [0, 1], Fh−1 (t, x, u, U ) has the same distribution of κ(R) and Hh (t, x, u, U ) is of distribution density x 7→ fh (t, x, u, x). In particular, it follows that the expression (3.7) of numerical solution of scheme (2.7) turns to be h Yk (x0 , · · · , xk ) = sup E hLk (x0 , · · · , xk , u) u∈E i + Yk+1 x0 , · · · , xk , xk + Hh (tk , x b, u, U ) , (3.9) where x b is the linear interpolation function of (x0 , · · · , xk ) on [0, tk ]. Now, we are ready to introduce a controlled discrete time semimartingale system. Suppose that U1 , · · · , Un are i.i.d random variables with uniform distribution on [0, 1] in the probability space (Ω, F, P). Let Ah denote the collection of all strategies φ = (φk )0≤k≤n−1 , where φk is a universally measurable mapping from (Rd )k+1 to E. Given φ ∈ Ah , X h,φ is defined by  h,φ [ h,φ , φ (X h,φ , · · · , X h,φ ), U := Xkh,φ + Hh tk , X X0h,φ := x0 , Xk+1 · k+1 . (3.10) k 0 k We then also define an optimization problem by V0h

:=

sup E φ∈Ah

h n−1 X

i  [ [ h,φ , φ ) + Φ(X h,φ ) . hL tk , X · k ·

(3.11)

k=0

The main result of this section is to show that the numerical solution given by (2.7) is equivalent to the value function of optimization problem (3.11) on the controlled discrete semimartingales X h,φ . Remark 3.3. It is clear that in the discrete time case, every process is a semimartingale. When µ ≡ 0 and U is of uniform distribution on [0, 1], the random variable Hh (t, x, u, U ) is centered, and hence X h,φ turns to be a controlled martingale. This is also the main reason we choose the terminology “semimartingale” in the section title. Theorem 3.4. Suppose that L and Φ satisfy (2.9) and Assumption 1 holds true. Then for 0 < h ≤ h0 with h0 defined by (2.11), Y0h = V0h . The above theorem is similar to a dynamic programming result. Namely, it states that optimizing the criteria globally in (3.11) is equivalent to optimizing it step by step in (3.9). With this interpretation, we only need to analyze the “distance” of the controlled semimartingale X h,φ in (3.10) and the controlled diffusion process X ν in (2.1) to show this convergence of V0h to V in order to prove Theorems 2.6 and 2.7. Before providing the proof, let us give a technical lemma.

11

Lemma 3.5. For the function G defined by (2.5) and every ε > 0, there is a universally measurable mapping uε : Sd × Rd → E such that for all (γ, p) ∈ Sd × Rd , G(t, x, γ, p) ≤ L(t, x, uε (γ, p)) +

1 t,x a ε · γ + bt,x uε (γ,p) · p + ε. 2 u (γ,p)

Proof. This follows from the measurable selection theorem, see e.g. Theorem 7.50 of Bertsekas and Shreve [2] or Section 2 of El Karoui and Tan [10]. Proof of Theorem 3.4: First, following (3.9), we can rewrite Ykh as a measurable function of (X00 , · · · , Xk0 ), and h Ykh (x0 , · · · , xk ) = sup E hLk (x0 , · · · , xk , u) u∈E i h + Yk+1 x0 , · · · , xk , xk + Hh (tk , x b, u, Uk+1 ) , where x b is the linear interpolation function of (x0 , · · · , xk ) on [0, tk ], Uk+1 is of uniform distribution on [0, 1] and Lk is defined by (3.6). Next, for every control strategy φ ∈ Ah and X h,φ defined by (3.10), we denote Fkh,φ := σ(X0h,φ , · · · , Xkh,φ ) and Vkh,φ := E

h n−1 X

i [ [ h,φ , φ ) + Φ(X h,φ ) F h,φ , hL(ti , X · i · k

i=k

which is clearly a measurable function of (X0h,φ , · · · , Xkh,φ ) and satisfies Vkh,φ (x0 , · · · , xk ) = hLk (x0 , · · · , xk , φk (x0 , · · · , xk )) h i h,φ +E Vk+1 x0 , · · · , xk , xk + Hh (tk , x b, φk (x0 , · · · , xk ), Uk+1 ) . Then by comparing V h,φ with Y h and the arbitrariness of φ ∈ Ah , it follows that V0h ≤ Y0h . For the reverse inequality, it is enough to find, for any ε > 0, a strategy φε ∈ Ah with X h,ε as defined in (3.10) using φε such that Y0h ≤ E

h n−1 X

i d d h,ε , φε ) + Φ(X h,ε ) + n ε. hL(tk , X · k ·

(3.12)

k=0

Let us write Γk and Zk defined below (2.7) as a measurable function of (X00 , · · · , Xk0 ), and uε be given by Lemma 3.5, denote φεk (x0 , · · · , xk ) := uε (Γk (x0 , · · · , xk ), Zk (x0 , · · · , xk z)). Then by the tower property, the semimartingale X h,ε defined by (3.10) with φε satisfies (3.12).

12

4

Proofs of the convergence theorems

With the probabilistic interpretation of the numerical solution Y h in Theorem 3.4, we are ready to give the proofs of Theorems 2.6 and 2.7. Intuitively, we shall analyze the “convergence” of the controlled semimartingale X h,φ in (3.10) to the controlled diffusion process X ν in (2.1).

4.1

Proof of Theorem 2.6

The main tool to prove Theorem 2.6 is the weak convergence technique due to Kushner and Dupuis [18]. We adapt their idea in our context. We shall also introduce an enlarged canonical space for control problems following El Karoui, Huu Nguyen and Jeanblanc [9], in order to explore the convergence conditions. Then we study the weak convergence of probability measures on the enlarged canonical space.

4.1.1

An enlarged canonical space

In Dolinsky, Nutz and Soner [8], the authors studied a similar but simpler problem in the context of G−expectation, where they use the canonical space Ωd := C([0, T ], Rd ). We refer to Stroock and Varadhan [23] for a presentation of basic properties of canonical space Ωd . However, we shall use an enlarge canonical space introduced by El Karoui et al. [9], which is more convenient to study the control problem for the purpose of numerical analysis.

An enlarged canonical space Let M([0, T ]×E) denote the space of all probability measures on [0, T ] × E, which is a Polish space equipped with the weak convergence topology. Denote by M the collection of measures m ∈ M([0, T ] × E) such that the projection of m on [0, T ] is the Lebesgue measure, so that they admit the disintegration m(dt, du) = m(t, du)dt, i.e. M :=



m ∈ M([0, T ] × E) : m(dt, du) = m(t, du)dt .

(4.1)

In particular, (m(t, du))0≤t≤T is a measure-valued process. The measures in space M are examples of Young measures and have been largely used in deterministic control problems. We also refer to Young [26] and Valadier [25] for a presentation of Young measure as well as its applications. Clearly, M is closed under weak convergence topology and hence is also a Polish space. We define also the σ−fields on M by Mt := σ{ms (ϕ), s ≤ t, ϕ ∈ Cb ([0, T ] × Rs E)}, where ms (ϕ) := 0 ϕ(r, u)m(dr, du). Then (Mt )0≤t≤T turns to be a filtration. In particular, MT is the Borel σ−field of M. As defined above, Ωd := C([0, T ], Rd ) is the space of all continuous paths between 0 and T equipped with canonical filtration d Fd = (Ftd )0≤t≤T . We then define an enlarged canonical space by Ω := Ωd × M, as d well as the canonical process X by Xt (ω) := ωtd , ∀ω = (ω d , m) ∈ Ω = Ωd × M, and d d d d the canonical filtration F = (F t )0≤t≤T with F t := Ftd ⊗ Mt . Denote also by M(Ω ) d the collection of all probability measures on Ω .

13

Four classes of probability measures A controlled diffusion process as well as d

the control process may induce a probability measure on Ω . Further, the optimization d criterion in (2.3) can be then given as a random variable defined on Ω . Then the d optimization problem (2.3) can be studied on Ω , as the quasi-sure approach in Soner, d Touzi and Zhang [20]. In the following, we introduce four subclasses of M(Ω ). δ

Let δ > 0, we consider a particular strategy ν δ of the form νsδ = wk (Xrνk , i ≤ Ik ) i

δ

for every s ∈ (kδ, (k + 1)δ], where Ik ∈ N, 0 ≤ r0k < · · · < rIkk ≤ kδ, X ν is the controlled process given by (2.1) with strategy ν δ , and wk : RdIk → E is a function such that there are finite number of disjoint hyper-rectangles which cover RdIk and wk is constant on each rectangle. Clearly, ν δ is an adapted piecewise constant strategy. Denote by U0 the collection of all strategies of this form for all δ > 0. It is clear that U0 ⊂ U. Given ν ∈ U in the probability space (Ω, F, P), denote mν (dt, du) := δνt (du)dt ∈ ν d M. Then (X ν , mν ) can induce a probability measure P on Ω by ν   EP Υ ω d , m := EP Υ(X·ν , mν , (4.2) d

for every bounded measurable function Υ defined on Ω . In particular, for any bounded function f : RdI+IJ → R with arbitrary I, J ∈ N, si ∈ [0, T ], ψj : [0, T ] × E → R bounded, Z si ν   P ν P ψj (νr )dr, i ≤ I, j ≤ J . E f Xsi , msi (ψj ), i ≤ I, j ≤ J = E f Xsi , 0

d

Then the first and the second subsets of M(Ω ) are given by  ν  ν P S0 := P : ν ∈ U0 and P S := P : ν∈U . Pn−1 δφk (du)1(tk ,tk+1 ] (dt), Now, let 0 < h ≤ h0 and φ ∈ Ah , denote mh,φ (dt, du) := k=0 h,φ [ h,φ , mh,φ ) and X be the discrete semimartingale defined by (3.10). It follows that (X induces a probability measure P we introduce is

h,φ

d

d

on Ω as in (4.2). Then the third subsets of M(Ω )

P h :=



P

h,φ

: φ ∈ Ah . d

Finally, for the fourth subset, we introduce a martingale problem on Ω . Let Lt,x,u be a functional operator defined by Lt,x,u ϕ := µ(t, x, u) · Dϕ +

1 σσ T (t, x, u) · D2 ϕ. 2

(4.3)

d

Then for every ϕ ∈ Cb∞ (Rd ), a process M (ϕ) is defined on Ω by Z Z 1 t Mt (ϕ) := ϕ(Xt ) − ϕ(X0 ) − Lt,X· ,u ϕ(Xs ) m(s, du) ds. 2 0 E d

(4.4)

Denote by P R the collection of all probability measures on Ω under which X0 = x0 d a.s. and Mt (ϕ) is a F −martingale for every ϕ ∈ Cb∞ (Rd ). In El Karoui et al. [9], a probability measure in P R is called a relaxed control rule.

14

d

Remark 4.1. We denote, by abuse of notations, the random processes in Ω Z Z d d d σσ T (t, ω d , u)m(t, du), µ(t, ω , u)m(t, du), a(t, ω , m) := µ(t, ω , m) := E

E d

which are clearly adapted to the filtration F . It follows that Mt (ϕ) defined by (4.4) is equivalent to Z 1 t 1 Mt (ϕ) := ϕ(Xt ) − ϕ(X0 ) − µ(s, ω d , m) · Dϕ(Xs ) + a(s, ω d , m) · D2 ϕ(Xs )ds. 2 0 2 Therefore, under any probability P ∈ P R , since a(s, ω d , m) is non-degenerate, there is ˜ on Ωd such that the canonical process can be represented as a Brownian motion W Z t Z t ˜ s. Xt = x0 + µ(s, X· , m)ds + a1/2 (s, X· , m)dW 0

0 ν

Moreover, it follows by Itˆ o’s formula as well as the definition of P in (4.2) that ν P ∈ P R for every ν ∈ U. In resume, we have P S0 ⊂ P S ⊂ P R . d

A completion space of Cb (Ω ) Now, let us introduce two random variables on d

Ω by T

Z

d

 L t, ω d , u m(dt, du) + Φ(ω d )

Ψ(ω) = Ψ(ω , m) :=

(4.5)

0

and Z

d

Ψh (ω) = Ψh (ω , m) :=

T

 Lh t, ω d , u m(dt, du) + Φ(ω d ),

0

where Lh (t, x, u) := L(tk , x, u) for every tk ≤ t ≤ tk+1 given the discretization parameter h := Tn . It follows by the uniform continuity of L that sup |Ψ(ω) − Ψh (ω)| ≤ ρ0 (h),

(4.6)

ω∈Ω

where ρ0 is the continuity module of L in t given before (2.9). Moreover, the optimization problem (2.3) and (3.11) are equivalent to V

=

sup EP [Ψ],

and V0h =

sup EP [Ψh ].

P∈P S

P∈P h d

Finally, we introduce a space L1∗ of random variables on Ω . Let  ∗ P := ∪0 0, there is ξε ∈ Cb (Ω ) such that supP∈P ∗ EP [|ξ − ξε |] ≤ ε. It follows that lim sup EPn [ξ] − EP [ξ] n→∞ h  i   ≤ lim sup EPn ξ − ξε + EPn [ξε ] − EP [ξε ] + EP ξε − ξ n→∞

≤ 2ε. Therefore, (4.8) holds true by the arbitrariness of ε. The next result shows that the random variable Ψ defined by (4.5) belongs to L1∗ , when L and Φ satisfy (2.9). Lemma 4.3. Suppose that Assumption 1 holds true, Φ and L satisfy (2.9). Then the random variable Ψ defined by (4.5) lies in L1∗ . Proof. We first claim that for every C > 0   sup EP eC|XT | < ∞, P∈P

(4.9)



h,φ

which implies that E[eC|Xn | ] < ∞ is uniformly bounded in h and in φ ∈ Ah . Let C0 > 0 such that |µ(t, x, u)| ≤ C0 , ∀(t, x, u) ∈ QT . Then (|Xkh,φ − C0 tk |)0≤k≤n is a subh,φ

martingale, and hence for every C1 > 0, (eC1 |Xk Therefore, by Doob’s inequality,

−C0 tk |

)0≤k≤n is also a submartingale.

q  h i h i h,φ h,φ h,φ  E sup eC1 |Xk | ≤ edC0 T E sup eC1 |Xk −C0 tk | ≤ 2e2dC0 T E e2C1 |Xn | , 0≤k≤n

0≤k≤n

where the last term is bounded uniformly in h and φ by the claim (4.9). With the same arguments for the continuous time case in spirit of Remark 4.1, it follows by (2.9) and (4.5) that   sup EP |Ψ|2 ≤ ∞. P∈P



  Similarly, we also have supP∈P ∗ EP |Ψ0 |2 ≤ ∞ for 0

d

Z

Ψ (ω , m) :=

T

L t, ω d , u) m(dt, du) + Φ(ω d ) .

0

16

(4.10)

Let ΦN := (−N ) ∨ (Φ ∧ N ), LN := (−N ) ∨ (L ∧ N ) and Z TZ ΨN (ω d , m) := LN t, ω d , u)m(dt, du) + ΦN (ω d ). 0

E d

Then, ΨN is bounded continuous in ω = (ω d , m), i.e. ΨN ∈ Cb (Ω ). It follows by Cauchy-Schwartz inequality that r r sup EP |Ψ − ΨN |2 sup P(|Ψ0 | > N ) sup EP |Ψ − ΨN | ≤ P∈P



P∈P



r



P∈P

sup EP |Ψ − ΨN |2

P∈P



1 sup EP [|Ψ0 |] √ ∗ N P∈P

r



→ 0,

where the last inequality is from P(|Ψ0 | > N ) ≤ N1 EP [|Ψ0 |]. And hence Ψ ∈ L1∗ . Therefore, it is enough to justify the claim (4.9) to conclude the proof. By Lemma 3.1, for every random variable R of density function fh (t, x, u, x) and every c ∈ Rd , we have E[ec·R ] ≤ eC2 h (1 + C2 h), where C2 := sup(t,x,u)∈QT conditional expectation on

1 T t,x 2 |c a0 c| h ec·Xn that

 1 T t,x + |bt,x u · c| + 2 |c au c| . It follows by taking

h

E[ec·Xn ] ≤ C0 (c) := sup eC2 T (1 + C2 h)T /h < ∞. h≤h0

Let c be the vectors of the form (0, · · · , 0, ±C, 0, · · · , 0)T , we can easily conclude that h EeC|Xn | is uniformly bounded for all h ≤ h0 and X h = X h,φ with φ ∈ Ah . Further more, in spirit of Remark 4.1 and by the same arguments as in the proof of Proposition 2.5, supP∈P R EP [eC|XT | ] is bounded. And therefore, we proved the claim (4.9). Finally, we finish this section by providing two convergence lemmas, but leave their proofs in Appendix. Lemma 4.4. (i) Let (Ph )0 0, we can construct d h d d two semimartingales X and X on (Ω, F, P) as well as ∆A and ∆B such that d d d (X , ∆A , ∆B ) in (Ω, F, P) has the same distribution as that of (X d , ∆Ad , ∆B d ) in (Ω, F, P). Moreover, h

d

Xk = X0 +

k X i=1

  d CH n d d h H(∆Ai , ∆B i , U i ) and P max X k − X k > Θ ≤ C 3 . (4.19) 1≤k≤n Θ

Proof. The proof is almost the same as that of Lemma 3.2 in Dolinsky [7]. It is enough to replace φk (ΘX0 , · · · , ΘXk−1 )Yk in his proof by H(φk (ΘX0 , · · · , ΘXk−1 ), U k ) in our context, where φk is a measurable Eh −valued function such that  1 1 d ). ∆Adk , ∆Bkd = φk−1 (X0d , · · · , Xk−1 h h The main idea is to use the techniques of invariance principle of Sakhanenko [22]. h

d

Remark 4.9. The process X is defined in (4.19) by characteristics of X , with U as well as function H. Remark 4.10. With Hh given by (3.8), set H(a, b, x) := Hh (0, 0, a + bbT h, b, x), ∀(a, b, x) ∈ Eh × Rd .

(4.20)

Since (a + bbT h, b) ∈ E for every (a, b) ∈ Eh , then the above function (4.20) is well defined. Moreover, it follows by Lemma 3.1 that H satisfies (4.18) with CH ≤ C0 n−3/2 , for C0 independent of n. In particular, let Θ = n−1/8 , then k  d X 1 1 d d P max X k − H(∆Ai , ∆B i , Ui ) > 8 ≤ C1 8 , 1≤k≤n n n i=1

for another constant C1 independent of n (or equivalently of h := T /n).

20

4.2.3

Proof of Theorem 2.7

By Theorem 3.4 and Lemma 4.7, we only need to prove separately that V d ≤ V0h + Cε h1/8−ε ,

(4.21)

V0h ≤ V d + Cε h1/8−ε ,

(4.22)

and h V0 − V d ≤ Ch1/8 , if L and Φ are bounded.

(4.23)

First inequality (4.21) For every ν ∈ U as well as the discrete time semimartingale X d,ν defined at the beginning of Section 4.2.1, we can construct, following Lemma 4.8, d d d h (X , ∆A , ∆B ) and X in a probability space (Ω, F, P) with H(a, b, x) := Hh (0, 0, a+ d d d bbT h, b, x) as in Remark 4.10, such that the law of (X , ∆A , ∆B ) under P is the same as (X d,ν , ∆Ad,ν , ∆B d,ν ) in (Ω, F, P), and (4.19) holds true for every Θ > 0. Fix Θ := h1/8 and denote  d h E := max |X k − X k | > Θ . 0≤k≤n

Denote also ν d = (ν dk )1≤k≤n

with ν dk :=

1 d d ∆Ak , ∆B k . h

 cd By the same arguments in proving the claim (4.9), we know that EP exp C(|X · | + ch  |X · |) is bounded by a constant independent of ν ∈ U. It follows by the definition of Ψd in (4.12) as well as (2.9) and (2.13) that for every ε > 0, there is a constant Cε independent of ν ∈ U such that i h  ch  cd ch  cd d h EP Ψd (X , ν d ) − Ψd (X , ν d ) ≤ CEP exp C(|X · | + |X · |) X · − X ·  ≤ Cε h1/8 + P(E)1/(1−8ε) ≤ Cε h1/8−ε ,

(4.24)

where the second inequality follows from H¨older inequality and Remark 4.10. Next, we claim that h n−1 X    ch ch i h EP Ψd (X , ν d ) = EP hL tk , X · , ν dk+1 + Φ(X · ) ≤ V0h .

(4.25)

k=0

Then by the arbitrariness of ν ∈ U, it follows by the definition of V d in (4.13) that (4.21) holds true. Hence we only need to prove the claim (4.25). We can use the randomness argument as in Dolinsky, Nutz and Soner [8] for provh ing their Proposition 3.5. By the expression of X in (4.19), using regular conditional ˜ together with inde˜ F, ˜ P) probability distribution, there is another probability space (Ω, ˜k )1≤k≤n , (U ˜ 0 )1≤k≤n and measurable pendent uniformly distributed random variables (U functions Πk : [0, 1]k × [0, 1]k → Eh such that with ˜k ) := Πk (U ˜1 , · · · , U ˜k , U ˜0, · · · , U ˜0) (∆A˜k , ∆B 1 k

21

and ˜ k := x0 + X

k X

˜i , U ˜i ), H(∆A˜i , ∆B

i=1

˜ equals to (X h , ∆Ad , ∆B d )1≤k≤n ˜ k , ∆A˜k , ∆B ˜k )1≤k≤n in (Ω, ˜ F, ˜ P) the distribution of (X k k k in (Ω, F, P). Denote for every u = (u1 , · · · , un ) ∈ [0, 1]n , ˜ u := x0 + X k

k X

 ˜1 , · · · , U ˜i , u1 , · · · , ui ), U ˜i . H Πi ( U

i=1

Since H is given by (4.20), it follows by the definition of Eh in (4.14) as well as that ˜1 , · · · , U ˜k , u1 , · · · , uk ), of V0h in (3.11) that, with strategy ν˜ku := n1 Πk (U ˜ u , ν˜u )] ≤ V h . E[Ψd (X 0 And hence,     h ˜ U˜ 0 , ν˜U˜ 0 ) = E Ψd (X , ν d ) = E Ψd (X

Z [0,1]n

˜ u , ν˜u )]du ≤ V h . E[Ψd (X 0

Therefore, we proved the claim, which concludes the proof of inequality (4.21).

Second inequality (4.22) Let F˜h (a, b, x) denote the distribution function of the random variable bh + a1/2 Wh , where Wh is of Gaussian distribution N (0, hId ). Denote also by F˜h−1 (a, b, x) the inverse function of x 7→ F˜h (a, b, x). Then for every X h,φ with h φ ∈ Ah , we can construct X as in (4.19) with H(a, b, x) = F˜ −1 (a, b, x) such that its h

distribution is closed to that of X h,φ . By the same arguments as in the proof of (4.21), we can prove (4.22).

Third inequality (4.23) When Φ and L are both bounded, we can improve the estimations in (4.24) to d h EP Ψd (X , ν d ) − Ψd (X , ν d )  cd ch  ≤ 2|Ψd |∞ P(E) + CEP exp C(|X · | + |X · |) h1/8 ≤ Ch1/8 . And all the other arguments in the proof of (4.21) and (4.22) hold still true. We hence conclude the proof of (4.23).

5

The implementation of the scheme

We shall discuss some issues for the implementation of the scheme (2.7).

22

5.1

The degenerate case

The numerical scheme (2.7) demands that σ(t, x, u) ≥ σ0 (t, x) > ε0 Id for every (t, x, u) ∈ QT in Assumption 1, which implies that the volatility part should be all non-degenerate. However, in many applications, we may encounter degenerate cases. √ Example 5.1. Let d = 1, and E = [a, a], µ ≡ 0 and σ(u) := u for u ∈ E. A concrete optimization problem is given by Z T   ν sup E Φ XT , νt dt . ν∈U

˜ ν := Introducing X t

Rt 0

0

νs ds, the above problem turns to be   ˜ Tν , sup E Φ XTν , X ν∈U

which can be considered in the framework of (2.3). However, the volatility matrix of ˜ ν ) is clear degenerate. the controlled process (X ν , X The above example is the case of variance option pricing problem in uncertain volatility model in finance. Example 5.2. An typical example of variance option is the option “call sharpe”, where the payoff function is given, with constants S0 and K, by +  S0 exp(XT − VT /2) − K √ Φ XT , VT := . VT To make the numerical scheme (2.7) implementable in the degenerate case, we can perturb the volatility matrix. Concretely, given an optimization problem (2.3) with coefficients µ and σ, we set 1/2 σ ε (t, x, u) := σσ T (t, x, u) + ε2 Id , aε (t, x, u) := σ ε (σ ε )T (t, x, u). (5.1) Clearly, aε is nondegenerate. Given ν ∈ U, let X ν,ε be the solution to SDE Z t Z t ν,ε ν,ε Xt = x0 + µ(s, X· , νs )ds + σ ε (s, X·ν,ε , νs )dWs . 0

(5.2)

0

Then a new optimization problem is given by hZ T i V ε := sup E L(t, X·ν,ε , νt )dt + Φ(X·ν,ε ) , ν∈U

0

which has no more degenerate problem. A similar idea was also illustrated in Guyon and Henry-Labord`ere [14] as well as in Jakobsen [16] for degenerate PDEs. We give here in addition an error analysis on the perturbation in our context. Proposition 5.3. Suppose that L and Φ satisfy conditions (2.9) and (2.13), functions µ and σ are uniformly bounded and satisfy (2.2). Then there is a constant C independent of ε such that |V − V ε | ≤ C ε.

23

˜ be two independent standard d−dimensional Brownian motions. Proof. Let W and W ˜ ν,ε which is the solution Then it is clear that (5.2) is equivalent, in distribution, to X to Z t Z t ν,ε ν,ε ˜ ν,ε , νs )dWs + εW ˜ s. ˜ ˜ σ(s, X µ(s, X· , νs )ds + Xt = x0 + · 0

0

˜ ν,ε |2 , then taking expectations and Applying Ito’s formula on the process |Xtν − X t using classical method with Gronwall lemma, we can get that there is a constant C independent of ε such that ˜ ν,ε 2 ≤ Cε2 . E sup Xtν − X t 0≤t≤T

Then using conditions (2.9) and (2.13), with Cauchy-Schwartz inequality, we conclude the proof.

5.2

The simulation-regression method

To make scheme (2.7) implementable, a natural technique is to use the simulationregression method to estimate the conditional expectations arising in the scheme. First, given a function basis, we propose a projection version of the scheme (2.7). Next, replacing the L2 −projection by least-square regression with empirical simulations of X 0 , it follows an implementable scheme. An error analysis of the simulation-regression method has been achieved by Gobet and Turkedjiev [13] in the context of BSDE numerical schemes. In this paper, we shall just describe the simulation-regression method for our scheme and leave the error analysis for further study.

5.2.1

The Markovian setting

In practice, we usually add new variables in the optimization problem and make the dynamic of X 0 (given by (2.6)) Markovian. Suppose that for d0 > 0, there are functions 0 0 0 µ0,k : Rd × Rd → Rd , σ 0,k : Rd × Rd → Sd and sk : Rd × Rd → Rd for every 1 ≤ k ≤ n such that c0 · ) = µ (X 0 , S 0 ), µ0 (tk , X 0,k k k

c0 · ) = σ 0,k (X 0 , S 0 ) σ0 (tk , X k k

0 0 ). Then (X 0 , S 0 ) and Sk+1 := sk+1 (Sk0 , Xk+1 k k 0≤k≤n is a Markovian process from (2.6). 0 Suppose further that there are functions (µk , σ k , Lk ) : Rd × Rd × E → Rd × Sd × R 0 and Φ : Rd × Rd → R such that

and

c0 · , u) = µ (X 0 , S 0 , u), µ(tk , X k k k 0 c 0 L(tk , X · , u) = Lk (Xk , Sk0 , u),

c0 · , u) = σ k (X 0 , S 0 , u), σ(tk , X k k   0 0 c 0 Φ X · = Φ X ,S . n

n

Then it is clear that the numerical solution Ykh of (2.7) can be represented as a measurable function of (Xk0 , Sk0 ), where the function G in (2.5) turns to be G(tk , x, s, γ, p)    1 T := sup Lk (x, s, u) + σ k σ Tk (x, s, u) − σ0,k σ0,k (x, s) · γ + µk (x, s, u) · p (. 5.3) 2 u∈E

24

Remark 5.4. In finance, when we consider the payoff functions of exotic options such as Asian option, lookback option, we can usually add the cumulative average, or cumulative maximum (minimum) to make the system Markovian.

5.2.2

The projection scheme

To simplify the notations, let us just give the scheme for the case d = d0 = 1, while that in general case can be easily deduced. Let (pYk,i )1≤i≤I, 0≤k≤n−1 be a family of basis functions where every pYk,i is function defined on R2 , so that SkY

:=

I nX

o αi pYk,i (Xk0 , Sk0 ), α ∈ RI ,

i=1

A projection operator PkY is defined by 2 PkY (U ) := arg min E U − S , ∀U ∈ L2 (Ω, FT ).

is a convex subclass of

L2 (Ω, F

T ).

S∈SkY

(5.4)

Γ Z Γ Similarly, with basis functions pZ k,i and pk,i , we can define Sk , Sk as well as the Γ Z projections operators Pk , Pk . And inspired by Gobet and Turkedjiev [13], we propose the following two projection schemes: with the same terminal condition

Yˆn = Φ(XT0 , ST0 ),

First scheme    ˆ k , Zˆk ) , ˆk = P Y Yˆk+1 + hG(tk , X 0 , S 0 , Γ  Y  k k    k  ∆Wk+1 T −1 Z ˆ ˆ Zk = Pk Yk+1 (σ0 ) , h     T ∆W ∆W  k+1 −hId −1 k+1 Γ ˆ k = P Γ Yˆk+1 (σ T )−1 . σ 2 0 0 k h

Second scheme    Pn−1 Y Y 0, S0, Γ ˆ ˆ ˆ ˆ  hG(t + , Z ) Y = P , X  i i , T k k i=k k k k  h i   P T )−1 ∆Wk+1 , 0, S0, Γ ˆ ˆ Zˆk = PkZ YˆT + n−1 hG(t , Z ) (σ , X i i k 0 i=k+1 k k h  i  h  T P ∆W  n−1 k+1 ∆Wk+1 −hId −1 0 0 ˆ ˆ T −1 Γ ˆ k = P Γ YˆT + σ . 2 0 i=k+1 hG(tk , Xk , Sk , Γi , Zi ) (σ0 ) k h We note that the numerical solutions are of the form Yˆk = yk (Xk0 , Sk0 ), Zˆk = ˆ k = γk (X 0 , S 0 ) with functions yk , zk , γk . zk (Xk0 , Sk0 ) and Γ k k

5.2.3

Empirical regression scheme

The simulation-regression scheme consists in simulating M empirical processes by (2.6), denoted by (X 0,m , S 0,m )1≤m≤M , then replacing the projection of (5.4) by empirical least square method to estimate functions yk , zk and γk . Concretely, with the simulation-regression method, the first scheme turns to be yk =

I X

y Y α bk,i pk,i ,

zk =

PI

z pZ , bk,i i=1 α k,i

i=1

γk =

I X i=1

25

γ Γ α bk,i pk,i ,

(5.5)

where α bky := arg min

α∈RI

M X I X m=1

αi pYk,i (Xk0,m , Sk0,m )

i=1

2 0,m 0,m − yk+1 (Xk+1 , Sk+1 ) − hG(tk , ·, γ, z)(Xk0,m , Sk0,m ) . and α by , α bγ are also given by the corresponding least square method. Similarly, we can easily get an empirical regression scheme for the second projection scheme. Finally, we finish by remarking that in error analysis as well as in practice, we usually need to use truncation method in formula (5.5) with the a priori estimations ˆ k ). of (Yˆk , Zˆk , Γ

A

Appendix

We shall give here the proofs of Lemmas 4.4 and 4.5. The arguments are mainly due to Section 8 of Kushner [17], we adapt his idea of proving his Theorems 8.1, 8.2 and 8.3 in our context. We first recall that given φh ∈ Ah , X h,φh is given by (3.10) and mh,φh (t, du) := h δ(φh )k (du) for t ∈ (tk , tk+1 ]. Denote φhs := (φh )k for s ∈ (tk , tk+1 ] and Eh,φ [·] := k h,φh h,φh h,φh h,φh d ∞ := σ(X0 , · · · , Xk ). Then for every ϕ ∈ Cb (R ), it follows E[·|Fk ] for Fk by Taylor expansion that i h Z tk+1 \ h,φ h h,φh h,φh h,φh  h,φh  \ h,φh ) ds + ε , (A.1) Ls,X h ,φs ϕ(X Ek ϕ(Xk+1 ) − ϕ(Xk ) = Ek s h tk

where |εh | ≤ C(h3/2 + hρ(h)) and ρ is continuity module of µ and σ given by (2.2), with a constant C depending on ϕ but independent of (h, φh ). Proof of Lemma 4.4. (i) First, let (Ph )h≤h0 be the sequence of probability measures d \ h,φh , mh,φh ) with φ ∈ A , on Ω given in the lemma. Suppose that Ph is induced by (X h h then by (3.4), it is clear that there is a constant C3 such that for all 0 ≤ s ≤ t ≤ T , sup 0 0, we consider the strategy νtδ = vk (Xs , s ≤ kδ) for t ∈ (δk, δ(k + 1)], where vk are measurable functions defined on C([0, δk], Rd ). Denote by P Sc the d ν δ , mν δ ) as in (4.2), with ν δ of this collection of all probability measures induced by (X form. Then it is clear that P S0 ⊂ P Sc ⊂ P S ⊂ P R . Proof of Lemma 4.5. First, by almost the same arguments as in Section 4 of El Karoui et al. [9] (especially that of Theorem 4.10) that for every P ∈ P R , there is a sequence of probability measures Pn in P Sc such that Pn → P, where the main idea is using Fleming’s chattering method to approximate a measure on [0, T ] × E by piecewise constant processes. We just remark that the uniform continuity of µ and σ w.r.t. u in (2.2) is needed here, and the “weak uniqueness” assumption in their paper is guaranteed by Lipschitz conditions on µ and σ. Then we conclude by the fact that we can approximate a measurable function vk (x) defined on C([0, δk], Rd ) by functions wk (xti , i ≤ I) which are constant on rectangles as shown in Theorem 7.1 of Kushner [17].

27

References [1] G. Barles and P.E. Souganidis, Convergence of approximation schemes for fully nonlinear second order equations, Asymptotic Anal., 4(3):271-283, 1991. [2] D.P. Bertsekas, and S.E. Shreve, Stochastic optimal control, the discrete time case, volume 139 of Mathematics in Science and Engineering, Academic Press Inc. [Harcourt Brace Jovanovich Publishers], New York, 1978. [3] J.F. Bonnans, E. Ottenwaelter and H. Zidani, A fast algorithm for the two dimensional HJB equation of stochastic control, M2AN Math. Model. Numer. Anal., 38(4):723–735, 2004. [4] B. Bouchard and N. Touzi, Discrete-time approximation and Monte-Carlo simulation of backward stochastic differential equations, Stochastic Process. Appl., 111(2):175-206, 2004. [5] P. Cheridito, H.M. Soner, N. Touzi, and N. Victoir, Second-order backward stochastic differential equations and fully nonlinear parabolic PDEs, Comm. Pure Appl. Math., 60(7):1081-1110, 2007. [6] K. Debrabant and E.R. Jakobsen, Semi-Lagrangian schemes for linear and fully nonlinear diffusion equations, preprint, 2009. [7] Y. Dolinsky, Numerical Schemes for G–Expectations, preprint, 2011. [8] Y. Dolinsky, M. Nutz and H.M. Soner, Weak approximations of G−expectations, Stochastic. Processes. Appl. 2012. [9] N. El Karoui, D. Huu Nguyen, and M. Jeanblanc-Picqu´e, Compactification methods in the control of degenerate diffusions: existence of an optimal control, Stochastics, 20:169–219, 1987. [10] N. El Karoui, X. Tan, Capacities and measurable selection: Applications in stochastic control, in preparation. [11] A. Fahim, N. Touzi, and X. Warin, A probabilitic numerical method for fully nonlinear parabolic PDEs, Ann. Appl. Probab., 21(4):1322-1364, 2011. [12] E. Gobet, J.P. Lemor, and X. Warin, A regression-based Monte-Carlo method to solve backward stochastic differential equations, Ann. Appl. Probab., 15(3):21722202, 2005. [13] E. Gobet and P. Turkedjiev, Approximation of discrete BSDE using least-squares regression, preprint, 2011. [14] J. Guyon and P. Henry-Labord`ere, Uncertain Volatility Model: A Monte-Carlo Approach, Journal of Computational Finance, 14(3), 2011. [15] J. Jacod, and J. M´emin, Sur un type de convergence interm´ediaire entre la convergence en loi et la convergence en probabilit´e, S´eminaire de Probabilit´e XV 1979/80, Lecture Notes in Mathematics, Vol 850, 529-546, 1981. [16] E. R. Jakobsen, On error bounds for approximation schemes for non-convex degenerate elliptic equations, BIT 44(2): 269-285, 2004.

28

[17] H.J. Kushner, Numerical methods for stochastic control problems in continuous time, SIAM J. Control Optim., 28(5):999-1084, 1990. [18] H.J. Kushner and P. Dupuis, Numerical methods for stochastic control problems in continuous time, Vol 24 of Applications of Mathematics, Springer-Verlag, 1992. [19] S. Peng, G-Expectation, G-Brownian motion, and related stochastic calculus of Itˆ o type, Stochastic Analysis and Applications, volume 2 of Abel Symp., Springer Berlin, 541-567, 2007. [20] M. Soner, N. Touzi and J. Zhang, Quasi-sure Stochastic Analysis through Aggregation, Electronic Journal of Probability, 2012. [21] H.M. Soner, N. Touzi and J. Zhang, Wellposedness of second order backward SDEs, Probability Theory and Related Fields, to appear. [22] A.I. Sakhanenko, A new way to obtain estimates in the invariance principle, High Dimensional Probability II, 221-243, 2000. [23] D.W. Stroock, and S.R.S. Varadhan, Multidimensional diffusion processes, volume 233 of Fundamental Principles of Mathematical Sciences, Springer-Verlag, Berlin, 1979. [24] X. Tan, A splitting method for fully nonlinear degenerate parabolic PDEs, preprint, 2012. [25] M. Valadier, A course on Young measures, Rend. Istit. Mat. Univ. Trieste, 26(suppl.):349:394 (1995), 1994. Workshop on measure theory and real analysis (Italian) (Grado, 1993). [26] L. C. Young, Lectures on the calculus of variations and optimal control theory, Foreword by Wendell H. Fleming. W.B. Saunders Co. Philadelphia, 1969. [27] J. Zhang, A numerical scheme for backward stochastic differential equations, Annals of Applied Probability, 14(1), 459-488, 2004.

29