Transcription

1 (CRR ) Richardson

2 Crank-Nicolson : Black-Scholes Euler Milstein

3 0 3 0 Options, Futures, and Other Derivatives, John Hull, ( ) Investment Science, David G. Luenberger, Oxford Univ Pr, (, ) Stochastic Calculus for Finance I: The Binomial Asset Pricing Model, S.E.Shreve, Springer, ( I, ) Stochastic Calculus for Finance II: Continuous-Time Models, S.E.Shreve, Springer, ( II, ) Arbitrage Theory in Continuous Time,T. Bjork, Oxford Univ Pr.,, Dynamic Asset Pricing Theory J.Darrell Duffie Princeton Univ Pr..,,,,,,, Stochastic Differential Equations, An introduction with application, B.Øksendal, ( ) Brownian Motion and Stochastic Calculus, I.Karazats and S.E.Shreve, Springer, ( ), III, Implementing Derivatives Models, L. Clewlow and C. Strickland, Wiley&Sons Options, Futures, and Other Derivatives, John Hull, ( ) Dynamic Asset Pricing Theory J.Darrell Duffie Princeton Univ Pr Quantitative Methods in Derivatives Pricing, D. Tavella, John Wiley& Sons Monte Carlo Methods in Financial Engineering, P. Glasserman, Springer

4 0 4 Glasserman Glasserman Numerical Recipes in C H. P. William ( C ) (C++ ) Numerical Recipes R RjpWiki ( R R-Tips ( Web PDF R Rmetrics ( RQuantLib ( RQuantLib QuantLib C++ R QuantLib Frequently Asked Questions in Quantitative Finance P. Wilmott Wilmott Wilmott ( wilmott.com/)

5 ( ) ,2 p, q (p, q > 0, p + q = 1) , 2 5 ( (15, 5) ) , 2 9 (3, 2) 3 (3, 2) ? ,2 3 1 x 1 2 x 2 3 ( ) ( ) ( x 1 + x 5 2 = 10 2) x 1, x 2 (x 1, x 2 ) = (1/10, 3/20) (x 1, x 2 ) 3 1,2 (x 1, x 2 ) = (1/10, 3/20) x 1 + 9x 2 = = p, q

6 (us, ds), S 2 (1 + r, 1 + r), 1 2 r 0 < d < 1 + r < u (1.1) (V 1, V 2 ) 3 v 1,2 us 1 + r V 1 S 1 ds 1 + r V x 1, x 2 3 ( ) ( ) ( ) us 1 + r V1 x 1 + x ds 2 =. 1 + r V 2 x 1 = V 1 V 2 S(u d), x 2 = dv 1 + uv 2 (1 + r)(u d) 3 v v = x 1 S + x 2 = 1 + r d u d V r + 1 r + u u d V r (1.2) 3 p, q (1.2) p = 1 + r d u d, q = 1 r + u u d (1.3) (1.1) 0 < p, q < 1, p + q = 1 ( p, q) Ẽ[ ] (1.2) [ ] V v = Ẽ 1 + r V P (V = V 1 ) = p, P (V = V 2 ) = q ( p, q) 1 S 1 Ẽ[S 1 ] = pus + qds = (1 + r)s

7 1 7 ( p, q) 1 Ẽ[(S 1 S)/S] 2 r r 1 (Ẽ[S 1/(1 + r)] = S) ( p, q) nikkei_l nikkei_r Index Index 1.1 ( ) ( ) (Ω, P, F) Ω X(ω) (ω Ω) X t ( ) {X(t, ω) t T T T t ω ) X(t, ω) ω Ω t X(t, ω) t X(t, ω) X X {X(t) t 1 < t 2 < < t n (X(t 1 ), X(t 2 ),..., X(t n ))

8 {W (t) (1) W (0) = 0 (2) {W (t) ω Ω W (t, ω) t (3) 0 < t 1 < t 2 < < t n n W (t 1 ) W (t 0 ), W (t 2 ) W (t 1 ),..., W (t n ) W (t n 1 ) ( ) (4) s t W (t) W (s) 0, t s (1) (4) Wiener X(t) = µt + σw (t) (µ, σ)- X(t) (µ, σ)- X(t) BM(µ, σ 2 ) (1.4) (µ, σ)- S(t) = S exp {µt + σw (t), S S > dx(t) = µ(x(t))dt + σ(x(t))dw (t), X(0) = x. (1.5) (1.5) {X(t) (1.5) X(t) = x + t 0 µ(x(s))ds + 3 (1.5) σ 0 t dx(t) = µ(x(t))dt, X(0) = x 0 σ(x(s))dw (s) (1.6) d X(t) = µ(x(t)), X(0) = x (1.7) dt µ(x(t)) = µx(t) (µ ) X(t) = xe µt

9 1 9 W Index X (1.5) (1.7) σ(x(t))dw (t) (1.5) (1.5) dx(t) X(t) t dx(t) X(t + t) X(t) dw (t) (1.5) X X(t + t) X(t) = µ(x(t)) t + σ(x(t))(w (t + t) W (t)) (1.8) W (t + t) W (t)

10 1 10 0, t (1.8) X(0) = x X( t) = X(0) + µ(x(0)) t + σ(x(0)) tz 1 X(2 t) = X( t) + µ(x( t)) t + σ(x( t)) tz 2. (1.9) X(n t) = X((n 1) t) + µ(x((n 1) t)) t + σ(x((n 1) t)) tz n X( t), X(2 t),..., X(n t),... Z 1,..., Z n {X(k t) (k = 0, 1,... ) (1.5) X(0), X( t), X(2 t),... {X(t) (1.5) {X(t) 1.2 µ(x(t)) 0, σ(x(t)) 1 dx(t) = dw (t) 1.3 µ(x(t)) µ( ), σ(x(t)) σ( ) dx(t) = µdt + σdw (t) (1.10) (µ, σ)- σ = 0 X(t) = µt σdw (t) µt 1.4 µ(x(t)) = µx(t), σ(x(t)) = σx(t) (µ, σ ) dx(t) = µx(t)dt + σx(t)dw (t) (1.11) Black-Scholes - (1.11) dx(t) X(t) = µdt + σdw (t) X(t + t) X(t) X(t) = µ t + σ(w (t + t) W (t)) {X(t) (1.11) N(µ t, σ 2 t) (1.11) X(t) = x exp {(µ 12 ) σ2 t + σw (t)

11 dx(t) = λ( X X(t))dt + σdw (t) (Ornstein-Uhlenbeck ) (1.12) dx(t) = λ( X X(t)) + σ X(t)dW (t) (Feller ) (1.13) *1 (1.10) X(t) (1.10) X((n + 1) t) = X(n t) + µ t + σ(w ((n + 1) t) W (n t)) = X(n t) + µ t + σ tz n, n = 0, 1,... [0, 1] N t = 1/N X(n t) = X n X n+1 = X n + µ t + σ tz n, n = 0, 1,... Z 1,..., Z n R 1.2 Bmotion <- function(mu, sigma){ T <- 1 N <- 250 Delta <- 1/N X0 <- 0 # #[0,T] #1 # X <- numeric(n) X[1] <- X0 #X for (i in 1:(N-1)){ X[i+1] <- X[i] + mu*delta + sigma*sqrt(delta)*rnorm(1) plot(x, type="l") # X matplot() 1.3 Bmotion2 <- function(mu, sigma){ *1

12 1 12 T <- 1 # N <- 250 #[0,T] Delta <- 1/N #1 Npath <- 30 # X0 <- 0 # X <- matrix(0, nrow = N, ncol =Npath) #X for (j in 1:Npath) { X[1,j] <- X0 for (i in 1:(N-1)){ X[i+1,j] <- X[i,j] + mu*delta + sigma*sqrt(delta)*rnorm(1) matplot(x, type="l") # X x(t) f : R R d dt f(x(t)) = f (x(t)) dx(t) dt {X(t) f(x) {f(x(t)) (1) {X(t) : dx(t) = µ(x(t))dt + σ(x(t))dw (t), X(0) = x. f : R R {f(x(t)) df(x(t)) = f (X(t))dX(t) f (X(t))(dX(t)) 2 [ = µ(x(t))f (X(t)) + 1 ] 2 σ2 (X(t))f (X(t)) dt + σ(x(t))f (X(t))dW (t) (1.14) (1.14) (dx(t)) 2 dw (t) dt (1) dt dt = 0 (2) dt dw (t) = 0 (3) dw (t) dw (t) = dt

13 1 13 (dx(t)) 2 (dx(t)) 2 = [µ(x(t))dt + σ(x(t))dw (t)] 2 = µ 2 (X(t))(dt) 2 + 2µX(t)σ(X(t))dtdW (t) + σ 2 (X(t))(dW (t)) 2 = σ 2 (X(t))dt (1.15) (2) {X(t) 2 f(t, x) {f(t, X(t)) df(t, X(t)) = f t (t, X(t))dt + f x (t, X(t))dX(t) f xx(t, X(t))(dX(t)) 2 [ = f t (t, X(t)) + f x (t, X(t))µ(X(t)) + 1 ] 2 f xx(t, X(t))σ 2 (X(t)) dt (1.16) + f x (t, X(t))σ(X(t))dW (t). f t = f/ t, f x = f/ x, f xx = 2 f/ x 2 (1.16) (dt) 2 dt dx(t) Y (t) = y 0 exp{µt + σw (t) (µ, σ ) X(t) = µt + σw (t), f(x) = y 0 exp(x) dy (t) = df(x(t)) = f (X(t))dX(t) f (X(t))(dX(t)) 2 = y 0 exp(x(t))(µdt + σdw (t)) y 0 exp(x(t))(µdt + σdw (t)) 2 ( = Y (t) µdt + σdw (t) + 1 ) 2 σ2 dt = Y (t) [(µ + 12 ) ] σ2 dt + σdw (t) (1.17) {X(t) dy (t) = (µ + 12 σ2 ) Y (t)dt + σy (t)dw (t) ( 1.4 ) T ds(t) = µs(t)dt + σs(t)dw (t), S(0) = S 0 > 0 (1.18)

14 1 14 (µ, σ ( ) ) T, K T K (S(T ) K) + x + = max{0, x T, K T K (K S(T )) + f : R R f(s(t )) (1.18) db(t) = rb(t)dt, B(0) = 1 B(t) = e rt (r 0) 1.1 T (S(T ) K) + ϕ(t) η(t) V (t) V (t) = ϕ(t)s(t) + η(t)b(t) (1.19) ϕ(t), η(t) self-financing V (t) dv (t) = ϕ(t)ds(t) + η(t)db(t) = [µϕ(t)s(t) + rη(t)b(t)] dt + σϕ(t)s(t)dw (t) (1.20) V (T ) = (S(T ) K) + ϕ, η t [0, T ] C(t) C(T ) = (S(T ) K) + C(t) 2 f C(t) = f(t, S(t)) df(t, S(t)) = f t (t, S(t))dt + f s (t, S(t))dS(t) f ss(t, S(t))(dS(t)) 2 [ = f t (t, S(t)) + µs(t)f s (t, S(t)) + 1 ] 2 σ2 S 2 (t)f ss (t, S(t)) dt (1.21) + σf s (t, S(t))S(t)dW (t)

15 1 15 V (T ) = (S(T ) K) + (= C(T )) (ϕ, η) V (t) = C(t), t [0, T ] (1.22) (1.22) (ϕ, η) (1.21) (1.20) µϕ(t)s(t) + rη(t)b(t) = f t (t, S(t)) + µs(t)f s (t, S(t)) σ2 S 2 (t)f ss (t, S(t)) (1.23) σϕ(t) = σf s (t, S(t)) (1.24) (1.19), (1.23), (1.24) ϕ, η, f ϕ(t) = f s (t, S(t)) (1.25) η(t) = f(t, S(t)) f s(t, S(t))S(t) B(t) (1.26) f t (t, S(t)) + rs(t)f s (t, S(t)) σ2 S 2 (t)f ss (t, S(t)) rf(t, S(t)) = 0 (1.27) f t (t, s) + rsf s (t, s) σ2 s 2 f ss (t, s) rf(t, s) = 0, f(t, s) = (s K) +, s > 0 (1.28) {f(t, S(t)) (ϕ, η) V (0) = C(0) = f(0, S 0 ) (1.28) (1.28) f(t, s) = sn(d + (T t, s)) Ke r(t t) N(d (T t, s)) (1.29) N d ± (t, x) = 1 [log σ xk ) ] (r t + ± σ2 t. (1.30) 2 (1.30) R BScall <- function(r, sigma, T, K, S0){ f(0, S 0 ) = S 0 N(d + (T, S 0 )) Ke rt N(d (T, S 0 )) #r: #sigma: #t: #T: #K: #S0:

16 1 16 dplus <- 1/(sigma*sqrt(T))*(log(S0/K)+(r+0.5*sigma^2)*T) dminus <- 1/(sigma*sqrt(T))*(log(S0/K)+(r-0.5*sigma^2)*T) price <- S0*pnorm(dplus) - K*exp(-r*T)*pnorm(dminus) return(price) r = 0.1, S 0 = 62, K = 60, σ = 0.2, T = 5 > BScall(0.1,0.2,5/12,60,62) [1] > ( ) T (S(T ) K) + C(0) = E Q [ e rt (S(T ) K) +] (1.31) E Q [ ] Q {S(t) ds(t) = rs(t)dt + σs(t)dw Q (t), S(0) = S 0 > 0 (1.32) {W Q (t) Q (1.18) S(t)dt ( ) r {W (t) {S(T ) (1.31) (1.32) Girsanov-, Radon-Nikodým 0 t < T C(t) C(t, x) = e r(t t) E Q [ (S(T ) K) + S(t) = x ] C(t, x) (1.28) Feynman-Kac T F F = f(s(t )) F = g({s(t) 0 t T ) E Q [ e rt F ]

17 (1.18) 1 1 δ ds(t) = µs(t)dt + σs(t)dw (t) δs(t)dt = (µ δ)s(t)dt + σs(t)dw (t) (1.33) (1.20) dv (t) = ϕ(t)ds(t) + η(t)db(t) + ϕ(t)δs(t)dt = [µϕ(t)s(t) + rη(t)b(t)]dt + σϕ(t)s(t)dw (t). (1.34) (1.20) C(t) = f(t, S(t)) df(t, S(t)) = [ f t (t, S(t)) + (µ δ)s(t)f s (t, S(t)) + 1 ] 2 σ2 S 2 (t)f ss (t, S(t)) dt + σf s (t, S(t))S(t)dW (t) (1.35) f t (t, s) + (r δ)sf s (t, s) σ2 s 2 f ss (t, s) rf(t, s) = 0, f(t, s) = (s K) +, s > 0 (1.36) sf s r δ f(t, s) = se δ(t t) N(d + (T t, s)) Ke r(t t) N(d (T t, s)) (1.37) d ± (t, x) = 1 σ t [log xk ) ] (r + δ ± σ2 t. (1.38) 2 {S(t) ds(t) = (r δ)s(t)dt + σs(t)dw Q (t), S(0) = S 0 > 0 f(s(t )) E Q [e rt f(s(t ))] (1.39)

18 S(t) ds(t) = µs(t)dt + σs(t)dw (t), S(0) = s 0 > 0 S(t) 1 B d (t) db d (t) = r d B d (t)dt, B d (0) = 1 r d S(t) (S(T ) K) + S(t) B f (t) S(t)B f (t) B f (t) db f (t) = r f B d (t)dt, B f (0) = 1 r f B f (t) ϕ(t) B d (t) η(t) V (t) = ϕ(t)s(t)b f (t) + η(t)b d (t) f t (t, s) + (r d r f )sf s (t, s) σ2 s 2 f ss (t, s) r d f(t, s) = 0, f(t, s) = (s K) +, s > 0 (1.40) sf s r d r f f(t, s) = se rf (T t) N(d + (T t, s)) Ke rd (T t) N(d (T t, s)) d ± (t, x) = 1 σ t [log xk ) ] (r + d r f ± σ2 t. 2 {S(t) ds(t) = (r d r f )S(t)dt + σs(t)dw Q (t), S(0) = S 0 > 0 f(s(t )) E Q [e rdt f(s(t ))] (1.41)

19 T t > 0 τ {τ < t t τ τ t f(s(t)) τ f(s(τ )) f(s(τ )) [ ] E e rτ f(s(τ )) τ τ V 0 = sup E [ e rτ f(s(τ)) ] (1.42) τ T 0,T T 0,T [0, T ] (1.42) (1.42) Black-Scholes * f ( ) 1 T ( ) { S(t)dt, f max T S(t), f (S(T )) I max S(t) < B 0 0 t T 0 t T (1.43) f I *2

20 Black-Scholes Black-Scholes Heston ds(t) = µs(t)dt + Y (t)s(t)dw 1 (t) dy (t) = α(m Y (t))dt + β Y (t)dw 2 (t) dw 1 (t) dw 2 (t) = ρdt, 1 ρ (1) Tree models (2) Finite-difference methods (3) Monte Carlo methods (1) (2) Tree & FD MC

21 (CRR ) Cox-Ross-Rubinstein ( Option Pricing: A Simplified Approach, J.Financial Economics, 7(1979) ) 2 ds(t) = rs(t)dt + σs(t)dw (t) (2.1) u, d, p (2.1) Black-Scholes : db(t) = rb(t)dt, B(0) = 1 (B(t) = e rt ) : ds(t) = rs(t)dt + σs(t)dw (t), S(0) = S 0 (2.2) : T (S(T ) K) + E[e rt (S(T ) K) + ] 2 *3 [0, T ] N t := T/N t i = i t (i = 0, 1,..., N) t 0, t 1,..., t N *3 Black-Scholes

22 2 22 S(t 0 ), S(t i ),..., S(t N ) S(t i+1 ) = D i+1 S(ti ), i = 0,..., N (2.3) D i P (D i = u) = p, P (D i = d) = 1 p, 0 < p < 1 S i = S(t i ) S i i 1 i S i = D j S 0 j=1 S i u k d i j S 0, j = 0,..., i P (S i = u j d i j S 0 ) = ( i j) p j (1 p) i j {S i i=0,1,...,n (2.2) u, d, p (2.2) S(t + t) = S(t) exp {(r 12 ) σ2 t + σ(w (t + t) W (t)) (2.4) W (t + t) W (t) N(0, t) S(t) S(t + t) (2.4) u, d, p (2.4) u, d, p us i S i ds i t i t i + t 2 t i S i s S i+1 E[S i+1 S i = s ] = pus + (1 p)ds (2.5) (2.1) S(t i ) = s S(t i + t) [ E [S(t i + t) S(t i ) = s ] = se exp {(r 12 ) ] σ2 t + σ(w (t + t) W (t)) = se r t. (2.6)

23 2 23 *4 (2.5) (2.6) e r t = pu + (1 p)d p = er t d u d (2.7) (1.3) 2 V [S i+1 S i = s] = E[Si+1 S 2 i = s] (E[S i+1 S i = s]) 2 = p(us) 2 + (1 p)(ds) 2 (pus + (1 p)ds) 2 = s 2 ( pu 2 + (1 p)d 2 (pu + (1 p)d) 2) = s 2 ( pu 2 + (1 p)d 2 e 2r t) (2.1) V [S(t i + t) S i = s] = s 2 e 2r t (e σ2 t 1) (2.8) pu 2 + (1 p)d 2 e 2r t = e 2r t (e σ2 t 1) e (2r+σ2) t = pu 2 + (1 p)d 2 (2.9) u, d, p (2.7) (2.9) CRR ud = 1 (2.10) (2.7) (2.9) e r t ( ) d u d u2 + 1 er t d d 2 = e (2r+σ2) t e r t (u + d) ud = e (2r+σ2 ) t u d (2.11) u = 1/d d u 2 (e r + e (r+σ2) t )u + 1 = 0 u > d u = β + β 2 1 *4 Z N(µ, σ 2 ) e Z [ E e Z] ) = exp (µ + σ2 2 V E [(e Z ) 2] = exp ( 2µ + 2σ 2) [ e Z] [ = E (e Z ) 2] ( [ E e Z]) 2 = e 2µ+2σ 2 e 2µ+σ2 = e 2µ+σ2 (e σ2 1).

24 2 24 β = (e r t + e (r+σ2 ) t )/2 d, p (2.7) (2.10) u = β + β 2 1 d = β β 2 1 p = er t d u d u, d β = e r t + e (r+σ2 ) t. 2 u = e σ t, d = e σ t t 0 u, d (2.7),(2.9),(2.10) *5 (2.10) p = 1 2 (2.7),(2.9),(2.12) u, d, t 0 ( ) u = e r t 1 + e σ2 t 1 ( ) d = e r t 1 e σ2 t 1 u e (r 1 2 σ2 ) t+σ t, d e (r 1 2 σ2 ) t σ t, (2.12) T ( ) t i j (j i) ( (i, j) ) C i,j (i, j) S i,j = S 0 u j d i j C N,j = (S N,j K) +, j = 0, 1,..., N (2.13) N j=0 ( N j ) p j (1 p) N j e rt (S N,j K) + (2.14) *5 u 1 + σ δt σ2 t, d 1 σ δt σ2 t, e r t 1 + r t (2.11) t e r t (u + d) ud (1 + r t)(2 + σ 2 t) r t + σ 2 t (2.11) e (2r+σ2) t 1 + 2r t + σ 2 t

25 2 25 (5.15) i N 1 C i,j = e r t (pc i+1,j+1 + (1 p)c i+1,j ), i = N 1, N 2,..., 0, j = 0, 1,..., i C 0,0 *6 (1) u, d, p (i, j) S i,j S i,j = S 0 u j d i j (2) C N,j = (S N,j K) +, j = 0, 1,..., N (3) (5.16) C i,j = e r t (pc i+1,j+1 + (1 p)c i+1,j ), C 0,0 R 2 CRREcall <- function(r, sigma, T, K, S0, N) { S <- matrix(0, nrow=n+1, ncol=n+1) # C <- matrix(0, nrow=n+1, ncol=n+1) # delta <- T/N # u <- exp(sigma*sqrt(delta)) d <- exp(-sigma*sqrt(delta)) p <- (exp(r*delta) - d)/(u - d) # # # # for (i in 0:N){ for (j in 0:i){ S[i+1,j+1] <- S0*u^j*d^(i-j) *6 R 1

26 2 26 # for (j in 0:N){ C[N+1, j+1] <- max(s[n+1,j+1] - K, 0) # # for (i in N:1) { for (j in 1:i) { C[i,j] <- exp(-r*delta)*(p*c[i+1,j+1] + (1-p)*C[i+1,j]) return(c[1,1]) # ( ) S[i,j] (i N) r = 0.1, σ = 0.2, T = 5, S 0 = 62, K = 60, N = 300 > CRREcall(0.1,0.2,5/12,60,62,300) [1] > 2.2 (2.3) h(x) *7 V 0 V 0 = sup E[e rτ h(s τ )] τ T 0,N T 0,N {0,..., N {S 1,..., S N τ V i (x) = sup E[e r(τ ti) h(s τ ) S i = x] τ T i,n τ = min{i {0,..., N V i (S i ) h(s i ) V i (x) V i (x) = max{h(x), E[e r t V i+1 (D i+1 x)]. (2.15) E[e r t V i+1 (D i+1 x)] i (2.15) *7 K h(x) = (x K) +, h(x) = (K x) +

27 2 27 V 0 (S 0 ) 2 V i,j = max { h(s i,j ), e r t (pv i+1,j+1 + (1 p)v i+1,j ) R 2 CRRAput <- function(r, sigma, T, K, S0, N) { S <- matrix(0, nrow=n+1, ncol=n+1) V <- matrix(0, nrow=n+1, ncol=n+1) # # delta <- T/N # u <- exp(sigma*sqrt(delta)) d <- exp(-sigma*sqrt(delta)) p <- (exp(r*delta) - d)/(u - d) # # # # for (i in 0:N){ for (j in 0:i){ S[i+1,j+1] <- S0*u^j*d^(i-j) # for (j in 0:N){ V[N+1, j+1] <- max(k-s[n+1,j+1], 0) # # for (i in N:1) { for (j in 1:i) { V[i,j] <- max(k-s[i,j],exp(-r*delta)*(p*v[i+1,j+1] + (1-p)*V[i+1,j])) # return(v[1,1])

28 2 28 r = 0.1, σ = 0.2, T = 5, S 0 = 62, K = 60, N = 100 > CRRAput(0.1,0.2,5/12,60,62,100) [1] > N CRR C(N) Black-Scholes C BS ( ) 1 C(N) C BS = O N 1/N 1/N 2 L-B. Chang and K. Palmer, Smooth convergence in the binomial model, Finance and Stochastics (2007) 11, M. S. Joshi, Achieving higher order convergence for the prices of European options in binomial tress, Mathematical Finance (2010) 20, N C(N) C BS C(N) K 3 Clewlow and Strickland Hull 2.4 Y. S. Tian, A flexible binomial option pricing model The Journal of Futures Markets, 19(1999), u, d, p K CRR Tian u, d. u = e σ t+λσ 2 t, d = e σ t+λσ 2 t. λ t 0 u, d CRR u = e σ t, d = e σ t λ K

29 log(error) 2 29 Price N logn 2.3 (error) log 10 ( ) 1

30 2 30 λ p = er t d u d 0 < p < 1, d < e r t < u λ r σ 2 1 σ t *8 λ λ (2.16) (1) CRR (λ = 0 ) K CRR S N,j = S 0 u j d N j j : S 0 u j d N j = K. j j = log (K/S 0) N log(d). (2.17) log (u/d) j j : [ ] log j (K/S0 ) N log(d) =. log (u/d) [x] x. j N (2) λ λ S 0 (u ) j (d ) N j = K λ = log (K/S 0) (2j N)σ t Nσ 2 t (2.17) log(k/s 0 ) (2.18) λ = 2(j j ) t σt j j 0.5 t 0 λ 0 t (2.16) u = e σ t+λσ 2 t, d = e σ t+λσ 2 t [ ] log j (K/S0 ) N log(d) = log (u/d) λ = log (K/S 0) (2j N)σ t Nσ 2 t *8 CRR CRR t

31 2 31 CRR Tian 2.4 Tian CRR N N price N 2.4 CRR Tian 2.1 Tian 2.5 V (t, S) t ( t S ) V S S f(x) (1 ) : : : : f(x + h) f(x) h f(x) f(x h) h f(x + h) f(x h) 2h

32 (0 ) V (S 0 ) ( ) h V (S 0 + h) V (S 0 ) = V (S 0 + h) V (S 0 ) (S 0 + h) S 0 h V (S 0 ) V (S 0 + h) (i, j) V i,j t V 1,1 V 1,0 S 1,1 S 1,0 V 0,1, V 0, 1, S 0,0, S 0,1, S 0, 1 ( 2.5 ) : V 0,1 V 0,0 S 0,1 S 0,0, : V 0,0 V 0, 1 S 0,0 S 0, 1, : V 0,1 V 0, 1 S 0,1 S 0, 1 2 V 0,1 S 0,1 S 1,2 S 1,1 S 1,1 S 2,0 S 0,0 S 1,0 S 1,0 S 0, 1 S 1, CRR S 0,1 = u d S 0, S 0, 1 = d u S 0, S 1,1 = 1 d S 0, S 1,0 = 1 u S 0, S 2,0 = 1 ud S 0 S 2,0 R 2 CRREcalldelta <- function(r, sigma, T, K, S0, N) { #

33 2 33 C <- matrix(0, nrow=n+1, ncol=n+3) # delta <- T/N # u <- exp(sigma*sqrt(delta)) d <- exp(-sigma*sqrt(delta)) p <- (exp(r*delta) - d)/(u - d) # # # # for (j in 0:N){ C[N+1, j+2] <- max(s0*u^j*d^(n-j) - K, 0) # C[N+1, N+3] <- max(s0*u^(n+1)/d -K, 0) C[N+1, 1] <- max(s0*d^(n+1)/u -K,0) # for (i in N:1) { for (j in 1:(i+2)) { C[i,j] <- exp(-r*delta)*(p*c[i+1,j+1] + (1-p)*C[i+1,j]) GreeksDelta <- (C[1,3]-C[1,1])/(S0*(u/d-d/u)) return(greeksdelta) # (1.29) V s = N(d +(T t, s)) r = 0.1, σ = 0.2, T = 5, S 0 = 62, K = > CRREcalldelta(0.1,0.2,5/12,60,62,100) [1] >

34 (2.2) S(t) = S(0) exp {(r 12 ) σ2 t + σw (t) (2.19) (5.8) {W (t) {W (t) {W i,j {S(t) {S i,j {W (t) {W i,j i=0,...,n,j=0,...,i CRR [0, T ] N t = T/N, t i = i t (i = 0, 1,..., N) W i,j = w W i+1,j+1 = w + u, W i+1,j = w + d (d < 0) W i,j w + u p w 1 p w + d W i+1,j+1, W i+1,j p, 1 p p u 1 p d ( d ) u, d, p 2 t i W i w W i+1 E[W i+1 W i = w] = p(w + u) + (1 p)(w + d) = w + pu + (1 p)d (2.20)

35 2 35 W (t i ) = w W (t i + t) E[W (t i + t) W (t i ) = w] = E[W (t i + t) W (t i ) + W (t i ) W (t i ) = w] = w + E[W (t i + t) W (t i )] = w (2.21) (5.9) (5.10) w + pu + (1 p)d = w pu + (1 p)d = 0 (2.22) 2 V [W i+1 W i = w] = p(w + u) 2 + (1 p)(w + d) 2 (p(w + u) + (1 p)(w + d)) 2 (2.23) = p(w + u) 2 + (1 p)(w + d) 2 w 2 V [W (t i + t) W (t i ) = w] = V [W (t i + t) W (t i ) + W (t i ) W (t i ) = w] = V [W (t i + t) W (t i )] = t (2.24) (5.12) (5.13) p(w + u) 2 + (1 p)(w + d) 2 w 2 = t (2.25) (5.11) (5.14) u, d, p 1 d = u (5.11) (5.14) pu + (1 p)( u) = 0, p(w + u) 2 + (1 p)(w u) 2 w 2 = t p = 1 2, u = t, d = t {W i,j S i,j = S(0) exp {(r 12 ) σ2 t i + σw i,j BtreeEcall <- function(r, sigma, T, K, S0, N) { W <- matrix(0, nrow=n+1, ncol=n+1) S <- matrix(0, nrow=n+1, ncol=n+1) C <- matrix(0, nrow=n+1, ncol=n+1) # # # delta <- T/N #

36 2 36 u <- sqrt(delta) d <- -sqrt(delta) p <- 0.5 # # # # mu <- r - 0.5*sigma^2 for (i in 0:N){ for (j in 0:i){ W[i+1,j+1] <- j*u + (i - j)*d #(2*j-i)*sqrt(delta) S[i+1,j+1] <- S0*exp(mu*i*delta + sigma*w[i+1,j+1]) # for (j in 0:N){ C[N+1, j+1] <- max(s[n+1,j+1] - K, 0) # for (i in N:1) { for (j in 1:i) { C[i,j] <- exp(-r*delta)*(p*c[i+1,j+1] + (1-p)*C[i+1,j]) return(c[1,1]) 2.7 Richardson [0, T ] N 1 h = T/N C(h) C(h), C(h/2), C(h/3),..., C(h/n) C(0) C(0) n i=1 ( 1) n i i n i!(n i)! C ( ) h i

37 2 37 n = 2 C(h), C(h/2) C(0) C( ) 2 C(0), C (0) C(h) C(0) + C (0)h ( ) h C C(0) + C (0) h 2 2 ( ) ( ) ( ) C(0) 1 2 C(h) C = (0) 2/h 2/h C(h/2) C(0) C(h) + 2C (h/2) (2.26) (2.26) C(0) = C(h) + ah + O(h 2 ) (2.27) h C(0) = C(h/2) + ah/2 + O(h 2 ) (2.28) (2.28) 2 (2.27) 2C(0) C(0) = 2C(h/2) C(h) + O(h 2 ) C(0) = C(h) + 2C(h/2) + O(h 2 ) h 2 C(h) + 2C(h/2) C(h) n = 3 C(0) 1 2 C(h) 4C (h/2) C (h/3) C(h) C(0) ( ) r = 0.1, S 0 = 62, K = 60, σ = 0.2, T = 5/12 Black-Scholes N = 300, 600 > CRREcallsmooth(0.1, 0.2, 5/12, 60, 62, 300) [1] > CRREcallsmooth(0.1, 0.2, 5/12, 60, 62, 600) [1] > (2.26) CRR ( )

38 Black-Scholes ,2 ds 1 (t) = rs 1 (t)dt + σ 1 S 1 (t)dw 1 (t) [ ds 2 (t) = rs 2 (t)dt + σ 2 S 2 (t)d ρw 1 (t) + ] 1 ρ 2 W 2 (t) 1 ρ 1, σ 1 > 0, σ 2 > 0 B(t) = e rt (r 0) (W 1, W 2 ) 2 S 2 W 3 = ρw ρ 2 W 2 W 1, W 3 ( ) ρ dw 1 dw 3 = dw 1 (ρdw ρ 2 dw 2 ) = ρdw 1 dw 1 = ρdt. d(w 1 (t)w 3 (t)) = W 3 (t)dw 1 (t) + W 1 (t)dw 3 (t) + dw 1 dw 3 (t) = W 3 (t)dw 1 (t) + W 1 (t)dw 3 (t) + ρdt. W 1 (t)w 3 (t) = ( 0 ) t 0 W 3 (u)dw 1 (u) + t E[W 1 (t)w 3 (t)] = ρt 0 W 1 (u)dw 3 (u) + ρt. W 1 (t), W 3 (t) t ρ W 1 (t), W 3 (t) E t [ (W 1 (t + t) W 1 (t))(w 3 (t + t) W 3 (t)) ] = ρ t ρw ρ 2 W 2 W 2 ds 1 (t) = rs 1 (t)dt + σ 1 S 1 (t)dw 1 (t) ds 2 (t) = rs 2 (t)dt + σ 2 S 2 (t)dw 2 (t) W 1, W 2 1 ρ 1 ( dw 1 dw 2 = ρdt ) (S 1, S 2 ) 2 [0, T ] N ρ W 1, W 2 2 {W 1 i,j j=0,1,...,i, i=0,1,...,n, {W 2 i,k k=0,1,...,i, i=0,1,...,n

39 2 39 Wi,j 1 i t j + 1 W i,k 2 Wi 1 = {Wi,j 1 j=0,...,i, Wi 2 = {Wi,k 2 k=0,...,i Wi,j 1 = w1, Wi,k 2 = w2 Wi+1,j+1 1 = w 1 + w 1, Wi+1,j 1 = w 1 w 1, Wi+1,k+1 2 = w 2 + w 2, Wi+1,k 2 = w 2 w 2 w 1, w 2 > 0 (w 1, w 2 ), (w 1 + w 1, w 2 + w 2 ) (w 1 w 1, w 2 + w 2 ) p uu p du (w 1, w 2 ) p dd pud (w 1 + w 1, w 2 w 2 ) (w 1 w 1, w 2 w 2 ) (w 1 + w 1, w 2 + w 2 ), (w 1 + w 1, w 2 w 2 ), (w 1 w 1, w 2 + w 2 ), (w 1 w 1, w 2 w 2 ) p uu, p ud, p du, p dd ( 2.7) w 1, w 2, p uu, p ud, p du, p dd 1 ( ) 2 ( ) 1 ( ) 2 2 Wi 1 = w 1 Wi+1 1 W i 1 E[Wi+1 1 Wi 1 Wi 1 = w 1 ] = p uu w 1 + p ud w 1 + p du ( w 1 ) + p dd ( w 1 ) = (p uu + p ud ) w 1 (p du + p dd ) w 1 W 1 (t i ) = w 1 W 1 (t i + t) W 1 (t i ) 0 p uu + p ud p du p dd = 0. (2.29) W 2 p uu p ud + p du p dd = 0. (2.30) 2 W 1 i = w 1 W 1 i+1 W 1 i 2 E[(W 1 i+1 W 1 i ) 2 W 1 i = w 1 ] = p uu ( w 1 ) 2 + p ud ( w 1 ) 2 + p du ( w 1 ) 2 + p dd ( w 1 ) 2 = (p uu + p ud + p du + p dd )( w 1 ) 2.

40 2 40 W 1 (t i ) = w 1 W 1 (t i + t) W 1 (t 1 ) 2 t W 2 p uu + p ud + p du + p dd = p uu + p ud + p du + p dd = ( ) W 1 i = w 1, W 2 i = w 2 E[(W 1 i+1 W 1 i )(W 2 i+1 W 2 i ) W 1 i = w 1, W 2 i = w 2 ] t ( w 1 ) 2. (2.31) t ( w 2 ) 2. (2.32) = p uu w 1 w 2 + p ud w 1 ( w 2 ) + p du ( w 1 ) w 2 + p dd ( w 1 )( w 2 ) = (p uu p ud p du + p dd ) w 1 w 2 E[(W 1 (t i + t) W 1 i )(W 2 (t i + t) W 2 i ) W 1 (t i ) = w 1, W 2 (t i ) = w 2 ] = ρ t (2.29)-(2.33) p uu p ud p du + p dd = p uu + p ud + p du + p dd = 1 ρ t w 1 w 2. (2.33) w 1 = w 2 = t, p uu = p dd = 1 4 (1 + ρ), p ud = p du = 1 (1 ρ) 4 h(s 1 (t), S 2 (t)) = (S 1 (t) S 2 (t) K) + (1) w 1, w 2, p uu, p ud, p du, p dd ρ 2 {W 1 i,j j=0,...,i,i=0,...,n, {W 2 i,k k=0,...,i,i=0,...,n W 1 i,j = (2j i) w 1, W 2 i,k = (2k i) w 2 {Si,j 1, {S2 i,k Si,j 1 = S 1 (0) exp {(r 12 ) (σ1 ) 2 i t + σ 1 Wi,j 1, Si,k 2 = S 2 (0) exp {(r 12 ) (σ2 ) 2 i t + σ 2 Wi,k 2 (2) V N,j,k = (S 1 N,j S 2 N,k K) +, j, k = 0, 1,..., N (3) V i,j,k = e r t (p uu V i+1,j+1,k+1 + p ud V i+1,j+1,k + p du V i+1,j,k+1 + p dd V i+1,j,k ) (4) V 0,0,0

41 dx(t) = µ(x(t))dt + σ(x(t))dw (t), X(0) = x 0. (2.34) [0, T ] N t = T/N t i = i t (i, j) (i = 0, 1,..., N, j = i, i + 1,..., 0, 1,..., i) X i,j X i = {X i,j i j i (2, 2) (1, 1) (2, 1) (0, 0) (1, 0) (2, 0) (1, 1) (2, 1) (2, 2) X i,j = x, X i+1,j+1 = x u, X i+1,j = x m, X i+1,j 1 = x d x x u, x m, x d p u, p m, p d (2.34) p u, p m, p d, x u, x m, x d i, x (2.34) X(t + t) X(t) = µ(x(t)) t + σ(x(t)) tz, z N(0, 1) X(t) = x X(t + t) X(t) E[X(t + t) X(t) X(t) = x] = µ(x) t, V [X(t + t) X(t) X(t) = x] = σ 2 (x) t (2.35)

42 2 42 p u x u x p m x m p d x d t i 2.9 t i + t 3 X i = x X i+1 X i E[X i+1 X i X i = x] = p u (x u x) + p m (x m x) + p d (x d x) = p u x u + p m x m + p d x d x. (2.35) E[X(t i + t) X(t i ) X(t i ) = x] = µ(x) t p u x u + p m x m + p d x d x = µ(x) t (2.36) 2 ( ) X i+1 X i X i = x 2 E [ (X i+1 X i ) 2 Xi = x] = p u (x u x) 2 + p m (x m x) 2 + p d (x d x) 2. (2.35) X(t i ) = x E [ (X(t i + t) X(t i )) 2 X(t i ) = x ] = V [X(t i + t) X(t i ) X(t i ) = x] + (E[X(t i + t) X(t i ) X(t i ) = x]) 2 = σ 2 (x) t + (µ(x) t) 2. p u (x u x) 2 + p m (x m x) 2 + p d (x d x) 2 = σ 2 (x) t + (µ(x) t) 2 (2.37) p u, p m, p d p u + p m + p d = 1 (2.38) (2.36), (2.37), (2.38) x u, x m, x u, p u, p m, p d 0 < p u, p m p d < 1 (2.39) p u, p m, p d (µ(x), σ(x) ) x (2.39)

43 (Black-Scholes ) ds(t) = rs(t)dt + σs(t)dw (t) (2.40) 3 2 (2.40) W (t) 3 {W i,j S i,j = S(0) exp {(r 12 ) σ2 t i + σw i,j (2.41) (2.36), (2.37) µ = 0, σ = 1 p u x u + p m x m + p d x d x = 0 p u (x u x) 2 + p m (x m x) 2 + p d (x d x) 2 = t (2.42) ( ) x (2.42) (2.38) (2.43) x u = x + x, x m = x, x d = x x p u x p d x = 0 p u = p d p u ( x) 2 + p d ( x) 2 = t p u + p d = p u = p d =, 0 < p u, p m, p d < 1 t 2( x) 2, p m = 1 t ( x) 2 0 < t ( x) 2 < 1 t (2.43) ( x) 2 x = 3 t p u = p d = 1 6, p m = 2 3 (1) x u = x + x, x m = x, x d = x x, p u = p d = 1 6, p m = {W i,j i=0,...,n,j= i,...,i W 0,0 = W (0) = 0 W i,j = j x ( i j i)

44 2 44 (2) {B i,j {S i,j S i,j = S(0) exp {(r 12 ) σ2 t i + σw i,j (3) V N,j = (S N,j K) +, j = N,..., 0,..., N (4) : V i,j = e r t (p u V i+1,j+1 + p m V i+1,j + p d V i+1,j 1 ). (5) V 0,0 3 TriEcall <- function(r, sigma, T, K, S0, N) { S <- matrix(0, nrow=n+1, ncol=2*n+1) V <- matrix(0, nrow=n+1, ncol=2*n+1) # # delta <- T/N # dx <- sqrt(3*delta) pu <- 1/6 pm <- 2/3 pd <- 1/6 # # # # # W <- matrix(0, nrow=n+1, ncol=2*n+1) W[1,1] <- 0 # S[1,1] <- S0 # # mu <- r - sigma^2/2 for (i in 2:(N+1)){ muti <- mu*delta*(i-1) W[i,i] <- 0 # S[i,i] <- S0*exp(muti + sigma*w[i,i]) # for (j in 1:(i-1)){

45 2 45 W[i,i+j] <- j*dx W[i,i-j] <- -j*dx # # S[i,i+j] <- S0*exp(muti + sigma*w[i,i+j]) # S[i,i-j] <- S0*exp(muti + sigma*w[i,i-j]) # # for (j in 1:(2*N+1)){ V[N+1, j] <- max(s[n+1,j] - K, 0) # # for (i in N:1) { for (j in 1:(2*i-1)) { V[i,j] <- exp(-r*delta)*(pu*v[i+1,j+2] + pm*v[i+1,j+1] + pd*v[i+1,j]) return(v[1,1]) > TriEcall(0.1,0.2,5/12,60,62,100) [1] > (CRR ) D. Brigo, F. Mercurio, Interest Rate Models Theory and Practice, Springer Finance ( ) dr(t) = (ϕ(t) a(t)r(t))dt + σ(r(t), t)dw (t) (2.44) Vasicek Ornstein-Uhlenbeck dr(t) = λ( r r(t))dt + σdw (t),

46 2 46 price N CIR Cox, Ingersoll, Ross Feller dr(t) = λ( r r(t)) + σ r(t)dw (t) ϕ(t), a(t), σ(r(t), t) (2.44) r(t + t) r(t) = (ϕ ar(t)) t + σ tz, z N(0, 1). t ϕ ar(t) > 0 r(t) < ϕ/a r(t + t) > r(t) r t ϕ ar(t) < 0 r(t) > ϕ/a r(t + t) < r(t) r r(t) ϕ/a ϕ/a (2.44) α(t) X(t) = r(t) α(t) X dx(t) = dr(t) α (t)dt = (ϕ(t) α (t) + a(t)r(t))dt + σ(r(t), t)dw (t) = (ϕ(t) α (t) a(t)α(t) a(t)x(t)) + σ(x(t) + α(t))dw (t) α α (t) = ϕ(t) a(t)α(t), α(0) = r(0) dx(t) = a(t)x(t) + ˆσ(X(t), t)dw (t), X(0) = 0 (2.45)

47 2 47 ˆα(X(t), t) = σ(x(t) + α(t)). (2.45) (2.44) r(t) = α(t) + X(t) α Hull-White dr(t) = ar(t)dt + σdw (t) (2.46) 3 (2.46) r(t) = r(0)e at + σ r(t) t 0 e a(t u) dw (u) E [r(t + t) r(t) = r] = re a t, V [r(t + t) r(t) = r] = σ2 ( 1 e 2a t ) 2a (2.46) E[r(t i + t) r(t i ) r(t i ) = r] = ar t E[(r(t i + t) r(t i )) 2 r(t i ) = r] = σ 2 t + ( ar t) 2 (2.47) r u = r + r, r m = r, r d = r r (2.47) 2 p u + p m + p d = 1 p u p d = ar t r p u + p d = σ2 t + (ar t) 2 ( r) 2 p u = σ2 t + (ar t) 2 ar t r 2( r) 2 p m = 1 σ2 t + (ar t) 2 ( r) 2 p d = σ2 t + (ar t) 2 + ar t r 2( r) 2 r = σ 3 t *9 r i,j = j r p u = σ2 j 2 ( t) 2 aj t 2 p m = 2 3 a2 j 2 ( t) 2 p d = σ2 j 2 ( t) 2 + aj t 2 * 10 0 < p u, p m, p d < a t < j < 3a t *9 0 < p u, p m, p d < 1 *10 p u, p m, p d j p u(j) p u

48 2 48 j j j (j < 6/(3a t)) 2.11 r u p u r m p m r p d r d t i 2.11 t i + t r u = r + 2 r, r m = r + r, r d = r 0 < p u, p m, p d < 1 p u = σ2 j 2 ( t) 2 + aj t 2 p m = 1 3 a2 j 2 ( t) 2 + 2aj t p d = σ2 j 2 ( t) 2 + 3aj t a t < j < 3 6 3a t j < 6/(3a t) ( j ) 2.12 r u = r, r m = r r, r d = r 2 r p u = σ2 j 2 ( t) 2 3aj t 2 p m = 1 3 a2 j 2 ( t) 2 + 2aj t p d = σ2 j 2 ( t) 2 aj t 2

49 2 49 r p u r u p m r m p d r d t i 2.12 t i + t 0 < p u, p m, p d < a t < j < a t j > 6/(3a t) ( j ) x, p u, p m, p d 3 (2.46) Vasicek Vasicek dr(t) = λ( r r(t))dt + σdw (t) (2.48) (2.48) r(t) = r + (r(0) r)e λt + σ t 0 e λ(t u) dw (u) E [r(t + t) r(t) = r] = r + (r r)e λ t, V [r(t + t) r(t) = r] = σ2 2λ (1 e 2λ t ) Hull-White α (t) = α(t) + λ r, α(0) = r(0) r(t) α(t) = e λt (λ r t 0 ) e λs ds + r(0) (2.49) dx(t) = λx(t) + σdw (t) r(t) = X(t) α(t) r(t)

50 CIR CIR dr(t) = λ( r r(t))dt + σ r(t)dw (t) (2.50) Nelson and Ramaswamy(1990) * 11 r(t) f(t, r(t)) f(t, r(t)) f f(t, r(t)) f r(t) f df(t, r(t)) = f f dt + t r dr(t) f [ 2 r 2 f = r λ( r r(t)) f r 2 σ2 r(t) + f t ] dt + f r σ r(t)dw (t) f r = 1 σr(t) r dr f(t, r) = 0 σ r (2.51) f(t, r(t)) 1 (2.51) f(t, r(t)) λ df(t, r(t)) = σ r(t) f(t, r) = 2 r σ ) ( r σ2 4λ r(t) dt + dw (t) f(t, r(t)) r(t) = σ 2 f(t, r(t)) 2 /4 r(t) /N (3 ) ( ) *11 Simple binomial processes as diffusion approximations in financial models, Review of Financial Studies 1990, vol.3,

51 Black-Scholes f (t, s) + rs f t s (t, s) σ2 s 2 2 f (t, s) rf(t, s) = 0, f(t, s) = h(s), s > 0 (3.1) s2 3.0 D. Tavella, Quantitative Methods in Derivatives Pricing, John Wiley& Sons D. Tavella, C. Randall, Pricing Financial Instruments: The Finite Difference Method, John Wiley& Sons R. Seydel, Tools for Computational Finance, Springer 3.1 (1) (3.1) f(t, s) (t, s) t [0, T ] M s s s max s min [s min, s max ] N N M (N + 1) (N + 1) (t i, s j ) (i, j) f i,j f(t i, s j ) (2) ( ) s = s max, s = s min, t = T f i,j f M,j = h(s j ) (j = 0, 1,..., N) (3) f(x) 1 2 : f(x + x) f(x) x : f(x) f(x x) x : f(x + x) f(x x) 2 x 2 f(x) x 2 = x ( ) f x

52 3 52 s max = s N s N 1. (t i, s j ) s 1 s min = s 0 0 = t 0 t 1... t M 1 T = t M 3.1 ( f(x + x) f(x) x ) f(x) f(x x) f(x + x) 2f(x) + f(x x) / x = x ( x) 2 (i, j) f/ t f i+1,j f i,j t t = T/M. f/ s f i,j+1 f i,j 1 2 s s = (s max s min )/N. t = t i t s f i,j+1 f i,j f i+1,j f i,j 1 {f i,j j=0,...,n t = t i+1 {f i+1,j j=0,...,n (4) ( ) 1 [0, T ], [s min, s max ] M, N t = T/M, s = (s max smin )/N t i = i t, s j = s min + j s (i = 0, 1,..., M, j = 0, 1,..., N)

53 3 53 s s min = t i+1 t i+1 t i+1 t i 1 s 1 s min = 0 (i + 1, j) f i+1,j f i,j t + rs j f i+1,j+1 f i+1,j 1 2 s σ2 s 2 f i+1,j+1 2f i+1,j + f i+1,j 1 j ( s) 2 rf i+1,j = 0. (3.2) f i,j = a j f i+1,j+1 + b j f i+1,j + c j f i+1,j 1 (3.3) a j = 1 2 (rj + (σj)2 ) t, b j = 1 ((σj) 2 + r) t, c j = 1 2 ( rj + (σj)2 ) t (3.4) f i+1,j+1 f i,j f i+1,j 3.2 f i+1,j 1 K f(t M, s j ) = { s j K, s j K 0, s j < K j = 0, 1,..., N 0 f(t i, s 0 ) = 0, i = 0, 1,..., N * 12 f(t, s) s e r(t t) K *12 t s p(s, t), c(s, t) c(t, s) p(t, s) = s e r(t t) K (p(t, s) = 0) c(t, s) = s e r(t t) K

54 3 54 f(t i, s N ) = s N e r(t ti) K, i = 0, 1,..., N. f(t, s) s * 13 R execall <- function(r, sigma, T, K, smax, M, N){ Deltat <- T/M Deltas <- smax/n # # f <- matrix(0, nrow=m+1, ncol=n+1) # # for (j in 1:(N+1)){ f[m+1, j] <- max((j-1)*deltas - K, 0) # # for (i in M:1){ f[i, 1] <- 0 # f[i, N+1] <- smax - exp(-r*(t - (i-1)*deltat))*k for (j in 2:N){ a <- (r*(j-1) + (sigma*(j-1))^2)*deltat/2 b <- 1 - ((sigma*(j-1))^2 + r)*deltat c <- (-r*(j-1) + (sigma*(j-1))^2)*deltat/2 f[i,j] <- a*f[i+1,j+1] + b*f[i+1,j] + c*f[i+1,j-1] return(f[1,]) r = 0.1, σ = 0.2, K = 60, S 0 = 62, T = 5 s max = 300, M = 5000, N = 300 > x <- execall(0.1,0.2,5/12,60,300,5000,300) > x[63] [1] *13 f(t, s)/ s = 1 s max

55 3 55 M = 300 > x <- execall(0.1,0.2,5/12,60,300,300,300) > x[63] [1] e t s ( ) (3.1) h(s) = (s K) + ( ) s = Ke x, t = T 2τ σ 2, q = 2r σ 2, (3.5) ( f(t, s) = f Ke x, T 2τ ) σ 2 = v(τ, x) (3.6) (3.7) g(τ, x) { v(τ, x) = K exp 1 ( ) 1 2 (q 1)x 2 (q 1)2 + q τ g(τ, x) (3.1) g τ = 2 g { x 2, y(0, x) = max e x 2 (q+1) e x 2 (q 1), 0. (3.8) t = 0 (3.8) y(τ, a) = y(τ, b) = 0 (a < b) (3.8) g i+1,j = λg i,j+1 + (1 2λ)g i,j + λg i,j 1, λ = τ ( x) 2 (3.9) g i,0 = g i,n = 0, 1 i M (3.9) g (i) = (g i,1, g i,2,..., g i,n 1 ) (i = 1,..., N) g (i+1) = Ag (i) 1 2λ λ 0 0 λ 1 2λ λ A = ((N 1) (N 1) ) λ 0 1 2λ

56 3 56 g (i) e (i) = ḡ (i) g (i) g (i) ḡ (i) ḡ (1) = Ag (0) + e (1) ḡ (2) = Aḡ (1) + e (2) = A(Ag (0) + e (1) ) + e (2) = A 2 g (0) + Ae (1) + e (2) ḡ (3) = Aḡ (2) + e (3) = = A 3 g (0) + A 2 e (1) + Ae (2) + e (3). ḡ (n) = A n g (0) + A n 1 e (1) + A n 2 e (2) + + e (n) n A n 1 e (1) 0 lim n An z = 0, z max µ A i < 1. i µ A i A A µ A i = 1 4λ sin 2 iπ 2N, i = 1,..., N 1 (Seydel Lemma 4.3 ) 0 0 < λ sin 2 iπ 2N < < λ < < τ ( x) 2 < 1 2 (3.10) 0 (3.10) τ, x t, s t, s 2 s t (3.2) f i+1,j f i,j f i+1,j f i,j t + rs j f i+1,j+1 f i+1,j 1 2 s σ2 s 2 j f i+1,j+1 2f i+1,j + f i+1,j 1 ( s) 2 rf i,j = 0

57 3 57 f i,j = r t (a jf i+1,j+1 + b j f i+1,j + c j f i+1,j 1 ) (3.11) a j = 1 2 (rj + (σj)2 ) t, b j = 1 (σj) 2 t, c j = 1 2 ( rj + (σj)2 ) t (3.12) * 14 a i + b i + c i = 1 0 < a i, b i, c i < 1 (3.13) a i, b i, c i (3.11) 1/(1 + r t) (3.11) 3 (3.13) s j = j s (3.14) (σj) 2 t < 1 (3.14) ( σ 2 s ) 2 j 1 t < 1 (j = 1,..., N) t < s (σn) 2 (3.15) (3.2) r = 0.1, σ = 0.2, K = 60, S 0 = 62, T = 5 s max = 300, N = 300 (3.15) M < 1 ( ) 2 M > 5 12 ( )2 = 1500 > x <- execall(0.1,0.2,5/12,60,300,1500,300) > x[63] [1] > > x <- execall(0.1,0.2,5/12,60,300,1000,300) > x[63] [1] e+68 > b j = 0 1 (σj) 2 t = 0 s = s j σ t (3.2) 2 (CRR ) *14

58 t i+1 t i t i t i+1 1 f i+1,j f i,j t + rs j f i,j+1 f i,j 1 2 s σ2 s 2 f i,j+1 2f i,j + f i,j 1 j ( s) 2 rf i,j = 0. (3.16) s min = 0 a j f i,j+1 + b j f i,j + c j f i,j 1 = f i+1,j (3.17) a j = 1 2 ( rj (σj)2 ) t, b j = 1 + ((σj) 2 + r) t, c j = 1 2 (rj (σj)2 ) t (3.18) (3.17) f i,j+1, f i,j, f i,j 1 f i+1,j (3.17) f i,j+1 f i,j f i+1,j f i,j a 1 f i,2 + b 1 f i,1 + c 1 f i,0 = f i+1,1 a 2 f i,3 + b 2 f i,2 + c 2 f i,1 = f i+1,2. (3.19) a N 2 f i,n 1 + b N 2 f i,n 2 + c N 2 f i,n 2 = f i+1,n 2 a N 1 f i,n + b M 1 f i,n 1 + c M 1 f i,n 2 = f i+1,n 1 b 1 a c 2 b 2 a c 3 b 3 a c N 2 b N 2 a N c N 1 b N 1 f i,1 f i,2 f i,3. f i,n 2 f i,n 1 = f i+1,1 c 1 f i,0 f i+1,2 f i+1,3. f i+1,n 2 f i+1,n 1 a N 1 f i,n. (3.20) {f i,j j=1,...,n R R Ax = y

59 3 59 solve(a,y) imecall <- function(r, sigma, T, K, smax, M, N){ Deltat <- T/M Deltas <- smax/n # # f <- matrix(0, nrow=m+1, ncol=n+1) # # A <- matrix(0, nrow=n-1, ncol=n-1) # a1 <- (-r - sigma^2)*deltat/2 b1 <- 1 + (sigma^2 + r)*deltat c1 <- (r - sigma^2)*deltat/2 A[1,1] <- b1 A[1,2] <- a1 #a_1 #b_1 #c_1 for (j in 2:(N-2)){ a <- (-r*j - (sigma*j)^2)*deltat/2 b <- 1 + ((sigma*j)^2 + r)*deltat c <- (r*j - (sigma*j)^2)*deltat/2 A[j,j+1] <- a A[j,j] <- b A[j,j-1] <- c anm1 <- (-r*(n-1) - (sigma*(n-1))^2)*deltat/2 #a_n - 1 bnm1 <- 1 + ((sigma*(n-1))^2+r)*deltat #b_n - 1 cnm1 <- (r*(n-1) - (sigma*(n-1))^2)*deltat/2 #c_n - 1 A[N-1,N-1] <- bnm1 A[N-1,N-2] <- cnm1 # for (j in 0:N){ f[m+1, j+1] <- max(j*deltas - K, 0) #

60 3 60 # x <- numeric(n-1) # y <- numeric(n-1) # for (i in M:1){ bl <- 0 # bu <- smax - exp(-r*(t - (i-1)*deltat))*k f[i,1] <- bl f[i,n+1] <- bu y[1:(n-1)] <- f[i+1,2:n] y[1] <- y[1] - c1*bl y[n-1] <- y[n-1] - anm1*bu x <- solve(a,y) # f[i,2:n] <- x return(f[1,]) > x <- imecall(0.1,0.2,5/12,60,300,300,300) > x[63] [1] s, t 3.4 Crank-Nicolson f(u) f(u ± u) = f(u) ± f u (u) u f 2 u 2 (u)( u)2 ± 1 3 f 6 u 3 (u)( u)3 + O ( ( u) 4). f f(u + u) f(u) (u) = + O ( u). u u

61 3 61 O( u) O( u) f(u + u) f(u u) 2 = f u (u) u + O ( ( u) 3) f f(u + u) f(u u) = + O ( ( u) 2). u 2 u O(( u) 2 ) f(u + u) + f(u u) = 2f(u) + 2 f u 2 (u)( u2 ) + O ( ( u) 4) 2 f f(u + u) 2f(u) + f(u ) = u2 ( u) 2 + O ( ( u) 2). 2 O(( u) 2 ) O ( t) + O ( ( s) 2) O(( t) 2 ) f i+1,j f i,j t + rs j f i+1,j+1 f i+1,j 1 2 s σ2 s 2 f i+1,j+1 2f i+1,j + f i+1,j 1 j ( s) 2 rf i+1,j = 0 (3.21) f i,j = a jf i+1,j+1 + b jf i+1,j + c jf i+1,j 1 (3.22) a j = 1 2 (rj + (σj)2 ) t, b j = 1 ((σj) 2 + r) t, c j = 1 2 ( rj + (σj)2 ) t. f i+1,j f i,j t + rs j f i,j+1 f i,j 1 2 s σ2 s 2 f i,j+1 2f i,j + f i,j 1 j ( s) 2 rf i,j = 0 (3.23) a j f i,j+1 + b j f i,j + c j f i,j 1 = f i+1,j (3.24) a j = 1 2 ( rj (σj)2 ) t, b j = 1 + ((σj) 2 + r) t, c j = 1 2 (rj (σj)2 ) t. (3.21) (3.23) 2 f i+1,j f i,j t σ2 s 2 j + 1 ( 2 rs fi+1,j+1 f i+1,j 1 j 2 s ( fi+1,j+1 2f i+1,j + f i+1,j 1 ( s) 2 + f ) i,j+1 f i,j 1 2 s + f i,j+1 2f i,j + f i,j 1 ( s) 2 ) 1 2 r(f i+1,j + f i,j ) = 0. (3.25) (3.25) t i+1 t i + t/2 (3.25) 1 f f (t + t/2 + t/2, s) f (t + t/2 t/2, s) f(t + t, s) f(t, s) (t + t/2, s) = t 2 t/2 t

62 3 62 f i,j+1 f i+1,j+1 f i,j f i+1,j f i,j t i t i + t 2 Crank-Nicolson t i+1 f i+1,j 1 t i + t/2 1 O(( t) 2 ). (3.25) 2 s 1 f s ( t + t 2 ± t 2, s ) ( f t + t ) s 2, s = f ( t + t ) s 2, s ± ( ( f t + t )) t t s 2, s 2 + O(( t)2 ) ( f t + t ) s 2, s = 1 ( ) f f (t + t, s) + (t, s) + O(( t) 2 ) 2 s s (3.25) O ( ( t) 2) + O ( ( s) 2) (3.25) ( (3.22) (3.24) 2 ) a j f i,j+1 + b j f i,j + c j f i,j 1 = a j f i+1,j+1 d j f i+1,j c j f i+1,j 1, j = 1, 2,..., N 1 a j = 1 4 ( rj (σj)2 ) t, b j = ((σj)2 + r) t, c j = 1 4 (rj (σj)2 ) t, d j = b j 2 b 1 a c 2 b 2 a c 3 b 3 a c N 2 b N 2 a N c N 1 b N 1 = f i,1 f i,2 f i,3. f i,n 2 f i,n 1 a 1 f i+1,2 d 1 f i+1,1 c 1 f i+1,0 c 1 f i,0 a 2 f i+1,3 d 2 f i+1,2 c 2 f i+1,1 a 3 f i+1,4 d 3 f i+1,3 c 3 f i+1,2. a N 2 f i+1,n 1 d N 2 f i+1,n 2 c N 2 f i+1,n 3 a N 1 f i+1,n d N 1 f i+1,n 1 c N 2 f i+1,n 2 a N 1 f i,n. (3.26)

63 3 63 R CrankNicolson <- function(r, sigma, T, K, smax, M, N){ Deltat <- T/M Deltas <- smax/n # # f <- matrix(0, nrow=m+1, ncol=n+1) # # a <- function(j){ return((-r*j - (sigma*j)^2)*deltat/4) b <- function(j){ return(1 + ((sigma*j)^2 + r)*deltat/2) c <- function(j){ return((r*j - (sigma*j)^2)*deltat/4) d <- function(j){ return(b(j)-2) # A <- matrix(0, nrow=n-1, ncol=n-1) # A[1,1] <- b(1) A[1,2] <- a(1) for (j in 2:(N-2)){ A[j,j+1] <- a(j) A[j,j] <- b(j) A[j,j-1] <- c(j) A[N-1,N-1] <- b(n-1) A[N-1,N-2] <- c(n-1) # for (j in 0:N){

64 3 64 f[m+1, j+1] <- max(j*deltas - K, 0) # # x <- numeric(n-1) # y <- numeric(n-1) # for (i in (M-1):0){ bl <- 0 # bu <- smax - exp(-r*(t - i*deltat))*k #---- y ---- for (j in 1:(N-1)){ y[j] <- -a(j)*f[i+2,j+2] - d(j)*f[i+2,j+1] - c(j)*f[i+2,j] y[1] <- y[1] - c(1)*bl y[n-1] <- y[n-1] - a(n-1)*bu # f[i+1,1] <- bl f[i+1,n+1] <- bu x <- solve(a,y) f[i+1,2:n] <- x return(f[1,]) # > x <- CrankNicolson(0.1,0.2,5/12,60,200,200,200) > x[63] [1] Crank-Nicolson 3.5 f/ S 0 f(s + s, t) f(s s, t) 2 s f 0,i f 0,i+1 f 0,i 1 2 s

65 3 65 Black-Scholes r = 0.1, σ = 0.2, T = 5, S 0 = 62, K = > x <- CrankNicolson(0.1,0.2,5/12,60,200,200,200) > x[63] [1] > (x[64]-x[62])/2 [1] > i <- 2:100 > delta <- (x[i+1] - x[i-1])/2 > plot(delta,type="l",xlab="stock Price", ylab="delta") 3.5 Crank-Nicolson 3.6 f(t, s) = [ ] sup E e r(τ t) h(s(τ)) S(t) = s τ T [t,t ] f f (t, s) + rs f t s (t, s) σ2 s 2 2 f (t, s) rf(t, s) 0, f(t, s) h(s), (t, s) [0, T ) [0, ) s2 ( f (t, s) + rs f t s (t, s) + 1 ) 2 σ2 s 2 2 f (t, s) rf(t, s) (f(t, s) h(s)) = 0, (t, s) [0, T ) [0, ) s2 f(t, s) = h(s) (3.27)

66 3 66 Lamberton and Lapeyre f i,j max(h(s j ), f i,j ) (3.28) R pmax(x,y) x, y max imaput <- function(r, sigma, T, K, smax, M, N){ Deltat <- T/M Deltas <- smax/n # # f <- matrix(0, nrow=m+1, ncol=n+1) # A <- matrix(0, nrow=n-1, ncol=n-1) # # a1 <- (-r - sigma^2)*deltat/2 b1 <- 1 + (sigma^2 + r)*deltat c1 <- (r - sigma^2)*deltat/2 A[1,1] <- b1 A[1,2] <- a1 #a_1 #b_1 #c_1 for (j in 2:(N-2)){ a <- (-r*j - (sigma*j)^2)*deltat/2 b <- 1 + ((sigma*j)^2 + r)*deltat c <- (r*j - (sigma*j)^2)*deltat/2 A[j,j+1] <- a A[j,j] <- b A[j,j-1] <- c anm1 <- (-r*(n-1) - (sigma*(n-1))^2)*deltat/2 #a_n - 1 bnm1 <- 1 + ((sigma*(n-1))^2+r)*deltat #b_n - 1 cnm1 <- (r*(n-1) - (sigma*(n-1))^2)*deltat/2 #c_n - 1 A[N-1,N-1] <- bnm1 A[N-1,N-2] <- cnm1 # for (j in 0:N){ f[m+1, j+1] <- max(k - j*deltas, 0) # #

67 3 67 x <- numeric(n-1) # y <- numeric(n-1) # payoff <- pmax(k - seq(0, N*Deltas, Deltas), 0) # for (i in (M-1):0){ bl <- exp(-r*(t - i*deltat))*k bu <- 0 # f[i+1,1] <- bl f[i+1,n+1] <- bu y[1:(n-1)] <- f[i+2,2:n] y[1] <- y[1] - c1*bl y[n-1] <- y[n-1] - anm1*bu x <- solve(a,y) f[i+1,2:n] <- x # f[i+1,] <- pmax(payoff, f[i+1,]) return(f[1,]) # > x <- imaput(0.1,0.2,5/12,60,100,100,100) > x[63] [1] > plot(x,type="l") > (3.28) (3.27) Lamberton and Lapeyre P. Jaillet, D Lamberton, and B Lapeyre, Variational Inequalities and the Pricing of American Options, Acta Aplicandae Mathematicae, 21(1990) (3.27) Tavella and Randall 3.7

68 3 68 Option Price Stock Price 3.6 h A h(s(t ))1 (,A) (max 0<t<T S(t)) h(s(t ))1 (A, ) (min 0<t<T S(t)) h(s(t ))1 (A, ) (max 0<t<T S(t)) h(s(t ))1 (,A) (min 0<t<T S(t)) (a, b) 1 (a,b) (x) := { 1, x (a, b) 0, x / (a, b) A 2 ( ) (S(T ) K) + 1 (,A) max S(t) 0<t<T ( )] V 0 = E [(S(T ) K) + 1 (,A) max S(t) 0<t<T (3.29) (3.29) ( ) (S(T ), max 0<t<T S(t)) (3.29) Shreve

69 3 69 f (t, s) + rs f t s (t, s) σ2 s 2 2 f (t, s) rf(t, s) = 0, (t, s) [0, T ) [0, A) (3.30) s2 v(t, A) = 0 t [0, T ) (3.31) v(t, s) = (s K) +. (3.32) Black-Scholes v(t, A) = 0 (3.7) (3.31) s min = 0, s max = A s max Dirichlet Neumann Dirichlet { V (t, S) = 0 S 0 r(t t) V (t, S) = S Ke S 2 V (t, S) = S Neumann { V S = 0 S 0 V S = 1 S 3.9 (3.1) x = log(s) g(t, x) = f(t, e x )(= f(t, s)) g ( g (t, x) + r 1 ) g t 2 σ2 x (t, x) σ2 2 g x 2 (t, x) rg(t, x) = 0, g(t, x) = (ex K) + (3.33) s x x

70 X E[X] ( ( ) ) X 1, X 2,... E[X 1 ] = µ(< ), V [X 1 ] < 1 n n X i µ a.s. as n. ( ) i=1 ( ) X 1, X 2,... E[X 1 ] = µ(< ), V [X 1 ] = σ 2 (< ) [ ( n ) ] 1 P σ X i nµ x N(x) as n. n i=1 N X 1, X 2,..., X n X ( X 1, X 2,..., X n X ) E[X] = µ X 1, X 2,..., X n µ ( ) X n = 1 n n X i (4.1) i=1 X n µ a.s. as n n X n E[X] (4.1) 1 E[ X n ] = 1 n n E[X i ] = µ i=1 n X n X n X n X n X n e 2 [ ] e 2 = E[( X n µ) 2 1 n ] = V X i = 1 n n n 2 V [X i ] = 1 n 2 nσ2 X = σ2 X (4.2) n i=1 i=1

71 4 71 σ X = V [X 1 ] e X n µ 2 e = σx n 2 O(1/ n) σ X σ X s n = 1 n (X i n 1 X n ) 2 X n i=1 s n n n [ ] Xn µ P s n / n < x N(x) N( x) P (1 δ) 1 N(z δ/2 ) = δ/2 ( X n s n n z δ/2, X n + s n n z δ/2 ) [ Xn µ ] s n < x N(x) N( x) n 95% X n µ 0.01 δ = 0.05(5%) z δ/ s n n n f(x)dx (4.3) X U(0, 1) (4.3) E[f(X)] d (d 2) f(x)dx (4.4) [0,1] d d U([0, 1] d ) X E[f(X)] d O(1/ n) (4.3) ( ) [0, 1] m 0 = x 0 < x 1 < < x m = 1 w 0, w 1,..., w m m w i f(x i ) i=0

72 O(1/m 2 ) d [0, 1] d m m 1/d m 1/d i 1=0 i 2=0 m 1/d i d =0 w i1 w i2... w id f(x i1, x i2,..., x ik ) O(1/m 2/d ) d m * X X 1, X 2,..., X n E[X] 1/n n i=1 X i X X dx(t) = µ(x(t))dt + σ(x(t))dw (t) (4.5) (4.5) Euler- X((j + 1) t) = X(j t) + µ(x(j t)) t + σ(x(j t)) tz j+1 t = T/N (N [0, T ] ) X(T ) X N (T ) X(T ) N X N (T ) E[X N (T )] X N 1, X N 2,..., X N n n XN (T ) X N n = 1 n i=1 X N i E[X N (T )] E[X(T )] E[X N (T )] X(T ) E[X(T )] S = 1 T T 0 S(t)dt (4.6) ( S K) + S(t 1 ), S(t 2 ),..., S(t n ) S (4.6) 0 = t 0 < t 1 < < t m = T S m = 1 m m S(t j ) j=0 *15 d 4

73 4 73 E[( S m K) + ] E[( S K) + ] X µ = E[X] 2 e 2 = E[( X µ) 2 ] = E[( X E[ X]) 2 ] + (E[ X] µ) 2 = V [ X] + ( X ) 2. X s 2 /n (s X ) X E[ X] E[X] (1998) * 16 (Mersenne Twister) R runif() X U[0, 1] F 1 (X) F P (F 1 (X) < x) = P (X < F (x)) = F (x) F 1 F 1 R rnorm() Glasserman ( ) ds(t) = S(t) [rdt + σdw (t)] (4.7) E[e rt (S(T ) K) + ] (4.7) S(T ) = S 0 exp {(r 12 ) σ2 T + σw (T ) = S 0 exp {(r 12 ) σ2 T + σ T Z S(T ) Z *16

74 4 74 (1) i = 1, 2,..., n (1)-(4) (2) Z i (3) S i (T ) = S 0 exp { (r 1 2 σ2) T + σ T Z i S i (T ) (4) V i = e rt (S i (T ) K) + (5) V = (V 1 + V V n )/n R MCBScall <- function(r, sigma, T, K, S0){ #r: #sigma: #T: #K: #S0: N < # value <- numeric(n) S <- numeric(n) # for (i in 1:N){ S[i] <- S0 * exp((r * sigma^2) * T + sigma * sqrt(t) * rnorm(1)) value[i] <- exp(-r*t)*max(s[i] - K, 0) return(mean(value)) r = 0.1, S 0 = 62, K = 60, σ = 0.2, T = 5 > MCBScall(0.1,0.2,5/12,60,62) [1] > Black-Scholes MCBScall2 <- function(r, sigma, T, K, S0){ N < rnd <- rnorm(n)

75 4 75 S <- S0 * exp((r * sigma^2) * T + sigma * sqrt(t) * rnd) return(mean(exp(-r*t) * pmax(s - K, 0))) = t 0 < t 1 < < t m = T S = 1 m m S(t j ) j=1 ( S K) + E[e rt ( S K) + ] (S(t 1 ), S(t 2 ),..., S(t m )) (4.7) S(t j ) = S(t j 1 ) exp {(r 12 ) σ2 (t j t j 1 ) + σ t j t j 1 Z j {Z j (1) i = 1, 2,..., n (1)-(5) (2) Zj i (j = 1,..., m) (3) S i (t j ) = S i (t j 1 ) exp {( r 1 2 σ2) (t j t j 1 ) + σ t j t j 1 Zj i S i (t j ) (j = 1,..., m) (4) S i = (S i (t 1 ) + S i (t 2 ) + + S i (t m ))/m (5) V i = e rt ( S i K) + (6) V = (V 1 + V V n )/n AsianCall <- function(r, sigma, T, K, S0, M){ # M: [0,T] N < # Deltat <- T/M # value <- numeric(n) for (i in 1:N){ S <- numeric(m+1) S[1] <- S0 for (j in 1:M){

76 4 76 S[j+1] <- S[j]*exp((r-0.5*sigma^2)*Deltat + sigma*sqrt(deltat)*rnorm(1)) # SA <- mean(s) # value[i] <- exp(-r*t)*max(sa - K, 0) return(mean(value)) 4.2 Y E[Y ] σ Y n n, σ Y Y σ Y Y σ Y Glasserman Y E[Y ] X E[X] ( ) (X i, Y i ) (i = 1, 2,..., n) (X, Y ) b Y i (b) = Y i b(x i E[X i ]) Y i (b) (i = 1, 2,..., n) Ȳ (b) = 1 n n Y i (b) = 1 n i=1 n (Y i b(x i E[X])) = Ȳ b( X E[X]) i=1 Ȳ (b) E[Y ] E[Ȳ (b)] = E[Ȳ b( X E[X])] = E[Ȳ ] = E[Y ] Ȳ Y b Ȳ (b) Ȳ Y i(b) σ Y (b) σ 2 Y (b) = V [Y i (b)] = V [Y i b(x i E[X])] = V [Y i ] + V [b(x i E[X])] 2bCov(Y i, X i E[X]) = σ 2 Y + b 2 σ 2 X 2bσ X σ Y ρ XY (4.8) σx 2 = V [X], σ2 Y = V [Y ], ρ XY X, Y { b 2 σx 2 0 < b < σ Y 2σ 2bσ X σ Y ρ XY 0 X ρ XY, ρ XY > 0 σ Y 2σ X ρ XY < b < 0, ρ XY < 0 σ Y (b) < σ Y Ȳ Ȳ (b) b = b = σ Y Cov(X, Y) ρ XY = σ X V [X]

77 4 77 σ Y (b) V [Ȳ (b )] V [Ȳ ] V [Ȳ (b )] V [Ȳ ] = σ2 Y (1 ρ2 XY )/n σ 2 Y /n = 1 ρ 2 XY. X Y 1 1 X X Y b n i=1 b n = (X i X)(Y i Ȳ ) n i=1 (X i X) (4.9) 2 (4.9) (Y 1,..., Y n ) (X 1,..., X n ) Glasserman Black-Scholes E[e rt (S(T ) K) + ] S(T ) ( ) 1 n E[e rt S(T )] = S(0) E[S(T )] = e rt S(0) n [ e rt (S i (T ) K) + b (S i (T ) e rt S(0)) ] i=1 E[e rt (S(T ) K) + ] S i (T ) (i = 1, 2,..., n) S(T ) R CVEcall <- function(r, sigma, T, K, S0){ N < # S <- S0 * exp((r * sigma^2) * T + sigma * sqrt(t) * rnorm(n)) #S(T) Y <- exp(-r * T) * pmax(s - K,0) # b <- cov(s,y)/var(s) Yb <- Y - b * (S - exp(r * T) * S0) #Y return(mean(yb)) #mean() r = 0.1, σ = 0.2, T = 5, K = 60, S 0 = ( ) * 17 T = 1 n = 250, t i = i/m (i = 0, 1,..., M) *17 Jarrow and Protter A short history of stochastic integration adn mathematical finance. The early years, In The Herman Rubin Festschrift, IMS Lecture Notes 45, 2004 Samuelson Jarrow and Protter (2004) This is the paper that first coined the terms European and American options. According to a private commu-

78 4 78 Simple Monte Carlo Control Variate price price N N 4.1 ( ) ( ) N S A = 1 n S(t i ) (4.10) n i=0 (S(t i ) ) ( S A K) + ( ) ( ) e rt (S(T ) K) R BScall() AsianCallControlVariate <- function(r, sigma, K, S0){ N < # T <- 1 M <- 250 Deltat <- T/M # #[0,T] M # BScallprice <- BScall(r, sigma, T, K, S0) #BS nication with R.C. Merton, prior to writing the paper, P. Samuelson went to Wall Street to discuss options with industry professionals. His Wall Street contact explained that there were two types of options available, one more complex - that could be exercised any time prior to maturiy, and one more simple - that could be exercised only at the maturity date, and that only the more sophisticated European mind (as opposed to the American mind) could understand the former. In response, when Samuelson wrote the paper, he used these as prefixes and reversed the ordering.

79 4 79 Acall <- numeric(n) Ecall <- numeric(n) for (i in 1:N){ S <- numeric(m+1) S[1] <- S0 for (j in 1:M){ S[j+1] <- S[j]*exp((r-0.5*sigma^2)*Deltat + sigma*sqrt(deltat)*rnorm(1)) # SA <- mean(s) # Ecall[i] <- exp(-r*t)*max(s[m+1] - K,0) # Acall[i] <- exp(-r*t)*max(sa - K, 0) # b <- cov(acall,ecall)/var(ecall) Yb <- Acall - b * (Ecall - BScallprice) return(mean(yb)) (r = 0.1, σ = 0.2, T = 1, S 0 = 62, K = 60) >AsianCallControlVariate(0.1,0.2,60,62) [1] > S(T ) = S(0) exp{(r σ 2 /2)T + σw (T ) n Z 1, Z 2,..., Z n (4.11) S i = S(0) exp ) {(r σ2 T + σ T Z i 2 Z 1, Z 2,..., Z n (4.12)