|
| 1 | +--- |
| 2 | +title: Kolmogorov Backward Equations |
| 3 | +author: Ashutosh Bharambe |
| 4 | +--- |
| 5 | + |
| 6 | +```julia |
| 7 | +using Flux, StochasticDiffEq |
| 8 | +using NeuralNetDiffEq |
| 9 | +using Plots |
| 10 | +using CuArrays |
| 11 | +using CUDAnative |
| 12 | +``` |
| 13 | +## Introduction on Backward Kolmogorov Equations |
| 14 | + |
| 15 | +The backward Kolmogorov Equation deals with a terminal condtion. |
| 16 | +The one dimensional backward kolmogorov equation that we are going to deal with is of the form : |
| 17 | + |
| 18 | +$$ |
| 19 | + \frac{\partial p}{\partial t} = -\mu(x)\frac{\partial p}{\partial x} - \frac{1}{2}{\sigma^2}(x)\frac{\partial^2 p}{\partial x^2} ,\hspace{0.5cm} p(T , x) = \varphi(x) |
| 20 | +$$ |
| 21 | +for all $ t \in{ [0 , T] } $ and for all $ x \in R^d $ |
| 22 | + |
| 23 | +#### The Black Scholes Model |
| 24 | + |
| 25 | +The Black-Scholes Model governs the price evolution of the European put or call option. In the below equation V is the price of some derivative , S is the Stock Price , r is the risk free interest |
| 26 | +rate and σ the volatility of the stock returns. The payoff at a time T is known to us. And this makes it a terminal PDE. In case of an European put option the PDE is: |
| 27 | +$$ |
| 28 | + \frac{\partial V}{\partial t} + rS\frac{\partial V}{\partial S} + \frac{1}{2}{\sigma^2}{S^2}\frac{\partial^2 V}{\partial S^2} -rV = 0 ,\hspace{0.5cm} V(T , S) = max\{\mathcal{K} - S , 0 \} |
| 29 | +$$ |
| 30 | +for all $ t \in{ [0 , T] } $ and for all $ S \in R^d $ |
| 31 | + |
| 32 | +In order to make the above equation in the form of the Backward - Kolmogorov PDE we should substitute |
| 33 | + |
| 34 | +$$ |
| 35 | + V(S , t) = e^{r(t-T)}p(S , t) |
| 36 | +$$ |
| 37 | +and thus we get |
| 38 | +$$ |
| 39 | + e^{r(t-T)}\frac{\partial p}{\partial t} + re^{r(t-T)}p(S , t) = -\mu(x)\frac{\partial p}{\partial x}e^{r(t-T)} - \frac{1}{2}{\sigma^2}(x)\frac{\partial^2 p}{\partial x^2}e^{r(t-T)} |
| 40 | + + re^{r(t-T)}p(S , t) |
| 41 | +$$ |
| 42 | +And the terminal condition |
| 43 | +$$ |
| 44 | + p(S , T) = max\{ \mathcal{K} - x , 0 \} |
| 45 | +$$ |
| 46 | +We will train our model and the model itself will be the solution of the equation |
| 47 | +## Defining the problem and the solver |
| 48 | +We should start defining the terminal condition for our equation: |
| 49 | +```julia |
| 50 | +function phi(xi) |
| 51 | + y = Float64[] |
| 52 | + K = 100 |
| 53 | + for x in eachcol(xi) |
| 54 | + val = max(K - maximum(x) , 0.00) |
| 55 | + y = push!(y , val) |
| 56 | + end |
| 57 | + y = reshape(y , 1 , size(y)[1] ) |
| 58 | + return y |
| 59 | +end |
| 60 | +``` |
| 61 | +Now we shall define the problem : |
| 62 | +We will define the σ and μ by comparing it to the orignal equation. The xspan is the span of initial stock prices. |
| 63 | +```julia |
| 64 | +d = 1 |
| 65 | +r = 0.02 |
| 66 | +sigma = 0.4 |
| 67 | +xspan = (80.00 , 115.0) |
| 68 | +tspan = (0.0 , 1.0) |
| 69 | +σ(du , u , p , t) = du .= sigma.*u |
| 70 | +μ(du , u , p , t) = du .= r.*u |
| 71 | +prob = KolmogorovPDEProblem(μ , σ , phi , xspan , tspan, d) |
| 72 | +``` |
| 73 | +Now once we have defined our problem it is necessary to define the parameters for the solver. |
| 74 | +```julia |
| 75 | +sdealg = EM() |
| 76 | +ensemblealg = EnsembleThreads() |
| 77 | +dt = 0.01 |
| 78 | +dx = 0.001 |
| 79 | +trajectories = 100000 |
| 80 | +``` |
| 81 | + |
| 82 | +Now lets define our model m and the optimiser |
| 83 | +```julia |
| 84 | +m = Chain(Dense(d, 8, leakyrelu),Dense(8, 16, leakyrelu),Dense(16 , 8 , leakyrelu) , Dense(8 , 1)) |
| 85 | +use_gpu = false |
| 86 | +if CUDAnative.functional() == true |
| 87 | + m = fmap(CuArrays.cu , m) |
| 88 | + use_gpu = true |
| 89 | +end |
| 90 | +opt = Flux.ADAM(0.01) |
| 91 | +``` |
| 92 | +And then finally call the solver |
| 93 | +```julia |
| 94 | +@time sol = solve(prob, NeuralNetDiffEq.NNKolmogorov(m, opt, sdealg, ensemblealg), verbose = true, dt = dt, |
| 95 | + dx = dx , trajectories = trajectories , abstol=1e-6, maxiters = 4200 , use_gpu = use_gpu) |
| 96 | +``` |
| 97 | +## Analyzing the solution |
| 98 | +Now let us find a Monte-Carlo Solution and plot the both: |
| 99 | +```julia |
| 100 | +monte_carlo_sol = [] |
| 101 | +x_out = collect(85:5.00:110.00) |
| 102 | +for x in x_out |
| 103 | + u₀= [x] |
| 104 | + g_val(du , u , p , t) = du .= 0.4.*u |
| 105 | + f_val(du , u , p , t) = du .= 0.02.*u |
| 106 | + dt = 0.01 |
| 107 | + tspan = (0.0,1.0) |
| 108 | + prob = SDEProblem(f_val,g_val,u₀,tspan) |
| 109 | + output_func(sol,i) = (sol[end],false) |
| 110 | + ensembleprob_val = EnsembleProblem(prob , output_func = output_func ) |
| 111 | + sim_val = solve(ensembleprob_val, EM(), EnsembleThreads() , dt=0.01, trajectories=100000,adaptive=false) |
| 112 | + s = reduce(hcat , sim_val.u) |
| 113 | + mean_phi = sum(phi(s))/length(phi(s)) |
| 114 | + global monte_carlo_sol = push!(monte_carlo_sol , mean_phi) |
| 115 | +end |
| 116 | + |
| 117 | +``` |
| 118 | + |
| 119 | +##Plotting the Solutions |
| 120 | +We should reshape the inputs and outputs to make it compatible with our model. This is the most important part. The algorithm gives a distributed function over all initial prices in the xspan. |
| 121 | +```julia |
| 122 | +x_model = reshape(x_out, 1 , size(x_out)[1]) |
| 123 | +if use_gpu == true |
| 124 | + m = fmap(cpu , m) |
| 125 | +end |
| 126 | +y_out = m(x_model) |
| 127 | +y_out = reshape(y_out , 6 , 1) |
| 128 | +``` |
| 129 | +And now finally we can plot the solutions |
| 130 | +```julia |
| 131 | +plot(x_out , y_out , lw = 3 , xaxis="Initial Stock Price", yaxis="Payoff" , label = "NNKolmogorov") |
| 132 | +plot!(x_out , monte_carlo_sol , lw = 3 , xaxis="Initial Stock Price", yaxis="Payoff" ,label = "Monte Carlo Solutions") |
| 133 | + |
| 134 | +``` |
0 commit comments