Theta-calculus for finance

Copyright 2004
Stefan Dirnstorfer, All rights reserved

License
This document is licensed under the terms of the GNU PUBLIC LICENSE, Version 2
See http://www.gnup.org/copyleft/gpl.html

Basic Operators

Shift operator

T[index_, power_, 
        op_] := ((op /. 
                index -> 
                  X00[t00]) /. DSolve[{X00'[
                        t00] == (power /. index -> X00[t00]), 
                    X00[0] == index}, X00[t00], 
                  t00][[1]][[1]]) /. t00 -> 1

T[x, a, f[x]]

T[x, r*x, f[x]]

Blur operator

B[index_, op_] := 
      Integrate[
        1/Sqrt[2*Pi]*
          Exp[-tt01^2/2]*(op /. 
              index -> 
                index + tt01), 
{tt01, -Infinity, 
          Infinity}]

B[index_, vola_, op_] := 
Simplify[
        T[index, -X2*vola, 
          B[X2, T[index, X2*vola, 
              op]]], 
{TimeConstraint -> 
            1}]

B[c, sigma,  V[c]]

B[x, sigma, B[x, sigma, a + b*x + c*x^2]]

B[x, sigma*x, V[x]]

Binomial blur

B2[index_,  vola_,  
        op_] := 
Expand[(T[index, vola, op]  + 
              T[index, -vola, op])/2]

B2[x, 1, a + b*x + c*x^2]

Financial operators

Interest rate curves

Y0 := {{0, 0.02}, {1, 0.034}, {5, 0.04},  {10, 0.06}}

List2Curve[list_] := 
      Function[M, 
          If[M ≤ 
              X2, 
((X2 - M)*Y1 + (M - X1)*
                    Y2)/(X2 - X1)
, 
            U1]]  /.  
{X1 -> list[[1]][[1]],  X2 -> list[[2]][[1]],  
          Y1 -> list[[1]][[2]],  
          Y2 -> list[[2]][[2]], 
U1 -> 
            If[Length[list] > 2, List2Curve[Rest[list]][
                M], list[[2]][[2]]]}

Z0 := List2Curve[Y0];

Plot[Z0[x], {x, 0, 10}];

Initial values

Z1[M_] := 
        M*(Z[M] + r1*1/100 + r2*1/100*(M - 2));

X0 := {if -> If,  Z -> Z0,  t ->  0,  
        c ->  0,  r1 -> 0,  r2 -> 0,  
        eur  ->  1}

Da[V_] := V  /.  X0;

Theta operator

Theta[dt_,  op_] :=  
T[
        eur,  -eur*(Z1[t + dt] - Z1[t]), 
T[t, 
          dt, 
op
]]
Theta[op_] := Theta[1, op]

Theta[5, t]
Da[%]

Theta[10, eur]
Da[%]
N[-Log[%]/10]

Interest rate extraction

forw[T1_, T2_] := 
      Theta[T1, -Log[Theta[T2 - T1, eur]/eur]/(T2 - T1)]

Da[forw[3, 5]]

-(Z0[3]*3 - 5*Z0[5])/(5 - 3)

Y[M_] := forw[0, M]

Da[Y[10]]
Da[Theta[3, Y[2]]]

Da[Y[1/5]]

Binomial tree Process

Theta2vola := 2/10

Theta2[V_] := 
B2[r1, Theta2vola, 
        Theta[V]]
Theta2[Dt_, V_] := 
B2[r1, Theta2vola*Sqrt[Dt], 
        Theta[Dt, V]]

Monte Carlo Process

Needs["Statistics`ContinuousDistributions`"]

randn := Random[NormalDistribution[0, 1]]

Table[randn, {i, 0, 9}]

Da[V_, n_] := Sum[Da[V  /.  N -> randn], {i, 1, n}]/n

Da[T[c, N, c], 1000]
Da[T[c, 1, T[c, c*N, c]], 1000]

Products

Zero Bond

P := Theta[10, eur]
Da[P]
Da[T[r1, 1, P]]
Da[D[P, r1]]
Da[D[P, r2]]
Da[D[P, r1, r1]]

Annuity

Op[V_] := Theta[T[c, a*eur, V]];

Op[c];
Da[%]

Op[Op[c]]
Da[%]

Nest[Op, c, 8];
Da[%]

L[l_, V_] := 
      If[l < 30,  
        Append[L[l + 1, Op[V]], {l, Da[V]}], {}]
L2 := 
      L[0, c] /. {eur -> 1,  a -> 1}
ListPlot[L2, PlotJoined -> True];

Zero Bond Option

u := Evaluate[Theta[10, eur]];
u
Da[u]

K := Da[u] - 1/100;
o := 
        Function[{r1, r2}, 
          eur*If[(E^((1/10) (-r1 - 8 r2 - 2 r2 t + 10 t Z[t] - 100 Z[10 + t] - 10 t Z[10 + t]))) > 
                K,  (E^((1/10) (-r1 - 8 r2 - 2 r2 t + 10 t Z[t] - 100 Z[10 + t] - 10 t Z[10 + t]))) - K,  
              0]];
Da[o[r1, r2]]

p := Nest[Theta2, o2[r1, r2], 20]
p /. o2 -> o;
N[Da[%]]

p := Nest[Theta2, o2[r1, r2], 20]
pg := 
        Evaluate[
          Da[T[r1, b, p] /. o2 -> o]];
Plot[{pg, 
Da[T[r1, b, Theta[3.55*eur]]] - 
          Da[Theta[3.55*eur]] + 0.026}, {b, -0.5, 1.5}]

Da[D[p /. o2 -> o, r1]]
Da[D[Theta[a*eur], r1]]
Solve[% == %%, a]

Da[D[p /. o2 -> o, r2]]

Da[D[Theta[a*eur + Theta[b*eur]], r1]]
Da[D[Theta[a*eur + Theta[b*eur]], r2]]
Solve[{%% == -0.034,  % == -0.275}, {a, 
        b}]

Plot[{Da[
            T[r1, b, 
              Theta[-28.4511*eur + 
                  Theta[16.5868*eur]]]]
 - 
          Da[Theta[-28.4511*eur + Theta[16.5868*eur]]] + 
          0.0263
, pg}, {b, -0.5, 1.5}]

pg2 := Evaluate[Da[T[r2, b, p] /. o2 -> o]];

Plot[{Da[
            T[r2, b, 
              Theta[-28.4511*eur + 
                  Theta[16.5868*eur]]]]
 - 
          Da[Theta[-28.4511*eur + Theta[16.5868*eur]]] + 
          0.0263, 
Da[
            T[r2, b, Theta[3.55*eur]]]
 - 
          Da[Theta[3.55*eur]] + 
          0.0263, 
pg2
}, {b, -0.5, 
        0.5}]

Bond deal

Clear[q]; Clear[p]

q[V_] := 
        T[c, -0.5*eur, Nest[Theta2, T[c, 1*eur, V], 10]];

p := q[f[eur]]

pvar := Evaluate[Da[q[If[c/eur > var, 0, 1]]]]

Plot[pvar
, {var, 0.0, 0.2}];

index linked cash flow

Clear[ptheta];  Clear[p]

pdt := 5;

ptheta[V_] := 
B2[index,  
        Sqrt[1/pdt]*2/1000, 
T[index,  
          index*1/5*Y[1/pdt]/pdt, 
Theta2[1/pdt, 
            V]]]
p[V_] := Nest[ptheta,  V,  pdt] /. index -> 1

payoff := eur*index;

Expected values

p[index];
Da[%]

p[payoff];
Da[%]

Da[p[eur]]

Distribtion

Clear[o]; Clear[o2]
o[V_] := If[V > var, 0, 1]

pvar := 
      Evaluate[Da[p[o2[payoff]] /. o2 -> o] /. eur -> 1]

Plot[pvar
, {var, 0.964, 0.98}];

m2 := Evaluate[Da[p[payoff^2]]]
m2

m1 := Evaluate[Da[p[payoff]]]
m1

vola := Sqrt[Abs[m1^2 - m2]]
vola

QQPlot

pdist[x_] :=  .5 +  .5*Erf[(x - m1)/vola]

ParametricPlot[{{pvar, pdist[var]}, {pdist[var], pdist[var]}}, {var, 
        0.8, 1.0}]

Sensitivity

Da[D[p[payoff], r1]]

p1 := Da[p[payoff]]
p2 := Da[Theta[eur]]
p1
p2

Da[Theta[a eur]]
ps := Evaluate[Solve[% == p1, a]]
ps

Da[D[Theta[(a /. ps[[1]])eur],r1]]

Solve[Da[D[Theta[b*eur], r1]]==-0.00778 eur,b]

modif := 
      Da[-100*D[p[payoff], r1]/p[payoff]]
modif

po := 
      Evaluate[Da[T[r1, b, p[payoff]]] /.  
          eur -> 1]
bo := 
      Evaluate[Da[T[r1, b, Theta[0.8 eur] + 0.2 eur]] /.  
          eur -> 1]

Plot[{po, bo},  {b, -10, 10}]

Hazard rate model

Clear[ptheta];  Clear[p]

lambda := (0.1*t)^(1)

ptheta[V_, pay_] :=  lambda*pay + (1 - lambda)*Theta[V]

p[V_] := Nest[Function[u, ptheta[u, V]], V, 12]

Da[p[t]]

Da[p[eur]]
Da[D[p[eur], r1]]

ptheta2[V_] := 
      B2[index,  2/1000,  
        T[index,  index*0.2*Y[1],  ptheta[V,  index*eur]]]

p2 := 
      Evaluate[Nest[ptheta2, index*eur, 10] /.  index -> 1]

Da[p2]
Da[T[r1, 1, p2]] - %

stop := 
      If[And[t < t0 ,  t > t0 - 1], 1, 0]
pstop := Nest[Function[u, ptheta[u, stop]], stop, 12]

Evaluate[Da[pstop]];
Plot[%,  {t0, -1,  12}]