Monte Carlo simulation for a doubly nonlinear problem in finance - - PowerPoint PPT Presentation

SLIDE 1

Monte Carlo simulation for a doubly nonlinear problem in finance

Lokman Abbas-Turki

First part from a joint work with M. A. Mikou Last part from a joint work with S. Graillat

UPMC, LPMA

10 June 2015

Lokman (UPMC, LPMA) Journées inaugurales 1 / 22

SLIDE 2

Plan

Introduction From linear/linear to linear/nonlinear From nonlinear/linear to nonlinear/nonlinear Simulation algorithms Without funding constraints With funding constraints Implementation issues

• n GPUs

The difference with the references Complexity vs. accuracy Some numerical results

Lokman (UPMC, LPMA) Journées inaugurales 2 / 22

SLIDE 3

Introduction

Plan

Introduction From linear/linear to linear/nonlinear From nonlinear/linear to nonlinear/nonlinear Simulation algorithms Without funding constraints With funding constraints Implementation issues

• n GPUs

The difference with the references Complexity vs. accuracy Some numerical results

Lokman (UPMC, LPMA) Journées inaugurales 3 / 22

SLIDE 4

Introduction From linear/linear to linear/nonlinear

Credit Valuation Adjustment

In a financial transaction between a party C that has to pay another party B some amount V , the CVA value is the price of the insurance contract that covers the default of party C to pay the whole sum V . CVAt,T = (1 − R)Et

• V +

τ 1t<τ≤T

• (1)

◮

R is the recovery to make if the counterparty defaults (Assume R = 0),

◮

τ is the random default time of the counterparty,

◮

T is the protection time horizon.

Numerical simulation

CVA0,T ≈

N−1

• k=0

E

• V +

tk 1τ∈(tk ,tk+1]

• ,

(2) N ≤ the number of time steps used for SDEs discretization.

Importance

◮

Hold sufficient amount of liquid assets to face the counterparty default.

◮

Basel III includes the calculation of the CVA (Credit Valuation Adjustment) as an important part of the prudential rules.

Lokman (UPMC, LPMA) Journées inaugurales 4 / 22

SLIDE 5

Introduction From linear/linear to linear/nonlinear

Various kind of contracts

Simulating assets St = (S1

t , ..., Sd t ) trajectories then contracts trajectories

to get Vt as a sum: Vt =

• ie

φexp

ie (St) +

• ii

φeui

ii (St) +

• id

φeud

id,t(St) +

• ia

φam

ia,t(St),

(3) where ie, ii, id and ia are the exposure indices and: φexp explicit function, for example: φexp(Stk ) = S1

tk − S2 tk .

φeui is a path-independent European contract, φeui(Stk ) = E(f eui(ST )|Stk ). (4) φeud

t

is a path-dependent European contract, φeud

tk (Stk ) = E(f eud tk

(Stk+1)|Stk ), (5) for example: f eud

tk

(Stk+1) = ( max

i=0,..,k S1 ti ∨ S1 tk+1 − S2 tk+1)+.

φam

t

is an American contract, involving an optimal stopping problem φam

tk (Stk ) = f (Stk ) ∨ E(φam tk+1(Stk+1)|Stk )

(6) with f an explicit payoff that generally does not depend on the asset path.

Common problem

ϕ(x) = E(f (Stk+1)|Stk = x).

Lokman (UPMC, LPMA) Journées inaugurales 5 / 22

SLIDE 6

Introduction From nonlinear/linear to nonlinear/nonlinear

TVA definition

Total valuation adjustment (>CVA+DVA+FVA), it covers:

◮

Both defaults: τ = τ c ∧ τ b, CVA and DVA.

◮

Funding our risk and the risk of the counterparty: Nonlinear BSDE part, FVA.

• S. Crépey(2012)

Ignoring the external funding and denoting βt = e−

t

0 rudu where r is the

risk-free short rate process, Θ satisfies the following BSDE on [0, τ ∧ T] βtΘt = E

• βτ1τ<T (Vτ − Rτ) +

τ∧T

t

βsgs(Vs − Θs)ds

• Gt
• (7)

where G is the extension of F by the natural filtration generated by τ c and by τ b. R is the total close-out cash-flow specified thanks to CSA (Credit Support Annex) and g is the funding coefficient.

TVA BSDE simulation

◮

Only for European contracts.

◮

Requires a good approximation of the exposure V .

◮

Practitioners usually use rough approximations.

◮

No trustable procedure in the general case.

Lokman (UPMC, LPMA) Journées inaugurales 6 / 22

SLIDE 7

Introduction From nonlinear/linear to nonlinear/nonlinear

WOLOG: V = P

For one American contract

◮

Pt is the American derivative exposition.

◮

Let τ ∗ ∈ [0, T] be the optimal stopping time associated to Pt.

Theorem

Θ = ˜ ΘJ + (1 − J)ξτ, where Jt = 1{τ>t} and the pre-default TVA ˜ Θ satisfies the following BSDE on [0, τ ∗] ˜ βt ˜ Θt = Et τ∗

t

˜ βs ˜ gs(Ps − ˜ Θs)ds

• , ˜

gt(Pt − ˜ Θt) = gt(Pt − ˜ Θt) + γt ˜ ξt. (8)

◮

˜ βt = e−

t

0(γu+ru)du. ◮

˜ ξt :=

1 Et(1{τ>t}) Et(ξt1{τ>t}) with ξt := Pt − E(R1{t=¯ τ}|Gt).

Extension

Possible to extend (8) on a portfolio of various American options.

Lokman (UPMC, LPMA) Journées inaugurales 7 / 22

SLIDE 8

Simulation algorithms

Plan

Introduction From linear/linear to linear/nonlinear From nonlinear/linear to nonlinear/nonlinear Simulation algorithms Without funding constraints With funding constraints Implementation issues

• n GPUs

The difference with the references Complexity vs. accuracy Some numerical results

Lokman (UPMC, LPMA) Journées inaugurales 8 / 22

SLIDE 9

Simulation algorithms

Without funding constraints

CVA0,T =

N−1

• k=0

E

• P+

k+11τ∈(kh,(k+1)h]

• ,

h = T N .

With funding constraints

Θk = Ek (Θk+1 + hg(k + 1, Pk+1, Θk+1)) , ΘN = 0.

The exposure

Pk(x) = E

• Φk,τk (Sτk )|Sk = x
• with

τN = N, ∀k ∈ {N − 1, ..., 0}, τk = k1Ak + τk+11Ac

k ,

(9) where Ak =

• Φk+1,k(Sk) > E
• Pk+1(Sk+1)
• Sk
• . The conditional

expectation involved in Ak is approximated using a regression on a basis

• f monomial functions where Kk is its cardinal.

Lokman (UPMC, LPMA) Journées inaugurales 9 / 22

SLIDE 10

Simulation algorithms

Without funding constraints

CVA0,T =

N−1

• k=0

E

• P+

k+11τ∈(kh,(k+1)h]

• ,

h = T N .

With funding constraints

Θk = Ek (Θk+1 + hg(k + 1, Pk+1, Θk+1)) , ΘN = 0.

An example of a two stage simulation with M0 = 2, M6 = 8 and M8 = 4

Lokman (UPMC, LPMA) Journées inaugurales 10 / 22

SLIDE 11

Simulation algorithms

CVA0,T approximation

• CVA0,T =

N−1

• k=0

1 M0

M0

• i=1

F 2

k+1

• P1(Si

1), ...,

Pk+1(Si

k+1)

• (10)

With F 1

k+1(x1, ..., xk+1) = E

• 1τ∈(kh,(k+1)h]|P1 = x1, ..., Pk+1 = xk+1
• ,

F 2

k+1(x1, ..., xk+1) = (xk+1)+F 1 k+1(x1, ..., xk+1).

(11)

Θk approximation

                   For k = 1, ..., N − 1

• Θk(x)= tψ(x)Ψ−1

k

  1 M0

M0

• j=1

ψ(Sj

k)

• Θk+1(Sj

k+1) + 1

N g

• k+1,

Θk+1(Sj

k+1),

Pk+1(Sj

k+1)

• 

 and ΘN(x) = 0,

• Θ0(S0) =

1 M0

M0

• j=1
• Θ1(Sj

1) + 1

N g

• 1,

Θ1(Sj

1),

P1(Sj

1)

• .

(12) Where Ψk = T   1 M0

M0

• i=0

ψ( Si

k)tψ(

Si

k)

  with: {Si}i∈{1,...,M0} and { Si}i∈{1,...,M0} are two independent simulations of the underlying asset S, ψ is a basis of monomial functions where K is its cardinal and T is an operator that must satisfy some desired properties.

Lokman (UPMC, LPMA) Journées inaugurales 11 / 22

SLIDE 12

Simulation algorithms Without funding constraints

Theorem

E

• CVA0,T − CVA0,T

2 ≤ N2 M0 max

k∈{0,...,N−1} Var

• F 2

k+1

• P1(Si

1), ...,

Pk+1(Si

k+1)

• +

N

• j=1

1 4NM2

j

• E
• Vj(Si

j )f ′′ j (Pj(Si j ))F 3 j (P1(Si 1), ..., Pj(Si j ))

2 +

N

• j=1

1 4NM2

j

  E   Vj(Si

j )f ′′ j (Pj(Si j )) N−1

• k=j

F 4

k+1(P1(Si 1), ..., Pj(Si j ), Si j )

   

2

+

N

• j=1

N 4M2

j

• E
• Vj(Si

j )F 1 j (P1(Si 1), ..., Pj(Si j ))|Pj(Si j ) = 0

• ϕj(0)

2 + N

N

• j=1

(N − j + 1)2O

• 1

M4

j

• Where

ϕj is the density of Pj(Si

j ), Vj(x) = Var

Mj

• Pj(x) − Pj(x)
• .

Good choices

If ϕj(0) is big then take Mj ∼ √M0, otherwise Mj ∼ √M0/N. In both cases, N must be small when compared to √M0. When American options are involved, make sure that K 3

j /Mj is small enough.

For example

Mj = N − j N − 1 M1 with either M1 = √M0 N

• r M1 =
• M0.

(13)

Lokman (UPMC, LPMA) Journées inaugurales 12 / 22

SLIDE 13

Simulation algorithms With funding constraints

Theorem

As long as {Θi(x)}0≤i≤N−1 are of class Cs on the support of S ∈ Rd, there exists a positive constant C such that for each 0 ≤ k ≤ N − 1 E

• Θk(Si

k) − Θk(Si k)

2 ≤ CK N2

N−1

• l=k

  E   Vl+1(Sj

l+1)∂2 Pg

• l + 1, Θl+1(Sj

l+1), Pl+1(Sj l+1)

• 2Ml

   

2

+O

• K

M0 + K 2 N2M0 + K N4M2

l

+ K 1−2s/d N2 + K −2s/d

• .

Good choice

Take Ml ∼

• M0/N, N must be sufficiently small N ∼ 10. When

American options are involved, make sure that K 3

l /Ml is small enough.

For example

Ml = N − l N − 1 M1 with M1 =

• M0

N . (14)

Lokman (UPMC, LPMA) Journées inaugurales 13 / 22

SLIDE 14

Implementation issues on GPUs

Plan

Introduction From linear/linear to linear/nonlinear From nonlinear/linear to nonlinear/nonlinear Simulation algorithms Without funding constraints With funding constraints Implementation issues

• n GPUs

The difference with the references Complexity vs. accuracy Some numerical results

Lokman (UPMC, LPMA) Journées inaugurales 14 / 22

SLIDE 15

Implementation issues on GPUs

An example of a two stage simulation with M0 = 2, M6 = 8 and M8 = 4 Inner dynamic programming

τN = N, ∀k ∈ {N − 1, ..., 0}, τk = k1Ak + τk+11Ac

k ,

(15) where Ak =

• Φk+1,k(Sk) > R · ψ(Sk)
• . The vector R minimizes the

quadratic error ||Pk+1(Sk+1) − R · ψ(Sk)||L2 (16) and denoting Ψ = E (ψ(Sk)tψ(Sk)), we get R = Ψ−1E (ψ(Sk)Pk+1(Sk+1)) . (17)

Lokman (UPMC, LPMA) Journées inaugurales 15 / 22

SLIDE 16

Implementation issues on GPUs The difference with the references

Three main methods for symmetric big matrices

Cholesky factorization

◮

• V. Volkov and J. Demmel. LU, QR and Cholesky Factorizations using

Vector Capabilities of GPUs. Berkeley Technical Report. 2008.

◮

• G. Ballard, J. Demmel, O. Holtz and O. Schwartz,

Communication-Optimal Parallel and Sequential Cholesky Decomposition. SIAM J. SCI. COMPUT. 32(6), 3495–3523. 2010.

Tridiagonal form + cyclic reduction

◮

• Y. Zhang , J. Cohen and J. D. Owens. 15th ACM SIGPLAN Symposium
• n Principles and Practice of Parallel Programming, 127–136. 2010.

◮

• D. Goddeke and R. Strzodka. Cyclic Reduction Tridiagonal Solvers on

GPUs Applied to Mixed Precision Multigrid. Parallel and Distributed Systems, IEEE Trans. 22(1), 22–32. 2010.

Tridiagonal form + eigenproblem

◮

• J. W. Demmel, O. A. Marques, B. N. Parlett and C. Vomel. Performance

and Accuracy of Lapack’s Symmetric Tridiagonal Eigensolvers. SIAM J.

• SCI. COMPUT. 30(3), 1508–1526. 2008.

◮

• C. Vomel, S. Tomov and J. Dongarra. Divide & Conquer on Hybrid

Gpu-Accelerated Multicore Systems. SIAM J. SCI. COMPUT. 34(2), 70–82. 2012.

Lokman (UPMC, LPMA) Journées inaugurales 16 / 22

SLIDE 17

Implementation issues on GPUs Complexity vs. accuracy

None of the previous works can be used directly

The reason

◮

Large number of small systems, the size is at most of the order of the warp size. Basically, the communication is much more reduced.

◮

Some of these systems could be ill-conditioned.

Parallel complexity = serial complexity

◮

Cholesky factorization operations less adapted to parallel architectures than Householder tridiagonalization.

◮

Cholesky factorization involves less cache memory than Householder tridiagonalization.

◮

Cyclic reduction and parallel cyclic reduction complexities can be neglected when compared to the Householder tridiagonalization complexity.

◮

The Householder tridiagonalization less adapted to parallel architectures than divide and conquer algorithm for symmetric tridiagonal eigenproblem.

Numerical accuracy

◮

Check the residual of the solutions for each procedure.

◮

Use CADNA to test each procedure.

Lokman (UPMC, LPMA) Journées inaugurales 17 / 22

SLIDE 18

Some numerical results

Plan

Introduction From linear/linear to linear/nonlinear From nonlinear/linear to nonlinear/nonlinear Simulation algorithms Without funding constraints With funding constraints Implementation issues

• n GPUs

The difference with the references Complexity vs. accuracy Some numerical results

Lokman (UPMC, LPMA) Journées inaugurales 18 / 22

SLIDE 19

Some numerical results

Within less than 1 minute simulation on GPU: M0 = 131K, N = 10, Neds = 50

European Path-dependent

• ption

Φ(ST ) =

• S1

T

2 + S2

T

2 − S

3 T

• +

M1 Θ0 Θ0 std CVA0,T CVA0,T std √M0 N 0.01364 4 ∗ 10−5 0.0296 2 ∗ 10−4 √M0 √ N 0.01307 4 ∗ 10−5 0.0294 2 ∗ 10−4 √M0 0.01265 3 ∗ 10−5 0.0291 2 ∗ 10−4

Lokman (UPMC, LPMA) Journées inaugurales 19 / 22

SLIDE 20

Some numerical results

Within less than 1 minute simulation on GPU: M0 = 131K, N = 10, Neds = 50

European Path-dependent

• ption

Φ(ST ) =

• 3S1

T

10 + 7S2

T

10 − S

3 T

• +

−

• 7S1

T

10 + 3S2

T

10 − S

3 T

• +

M1 Θ0 Θ0 std CVA0,T CVA0,T std √M0 N 2.72 ∗ 10−3 10−5 0.0365 8 ∗ 10−4 √M0 √ N 2.44 ∗ 10−3 10−5 0.0453 8 ∗ 10−4 √M0 2.28 ∗ 10−3 10−5 0.0520 8 ∗ 10−4 √ N√M0 2.24 × 10−3 10−5 0.0528 8 × 10−4

Lokman (UPMC, LPMA) Journées inaugurales 20 / 22

SLIDE 21

Some numerical results

Within less than 1 minute simulation on GPU: M0 = 131K, N = 10, Neds = 50

American option

Φ(ST ) =

• K − S1

T

3 − S2

T

3 − S3

T

3

• +

M1 Θ0 Θ0 std CVA0,T CVA0,T std √M0 √ N 0.0242 10−4 0.0356 2 ∗ 10−4 √M0 0.0229 10−4 0.0351 2 ∗ 10−4

◮

L.A. Abbas-Turki and M.A. Mikou. TVA on American Derivatives: https://hal.archives-ouvertes.fr/hal01142874

Lokman (UPMC, LPMA) Journées inaugurales 21 / 22

SLIDE 22

The End

Thank you Questions?

Lokman (UPMC, LPMA) Journées inaugurales 22 / 22