Multi-Risk Derivatives Valuation

A specialty of DX Analytics is the valuation of derivatives instruments defined on multiple risk factors and portfolios composed of such derivatives. This section of the documentation illustrates the usage of the dedicated multi-risk valuation classes.



In [1]:

    
from dx import *



In [2]:

    
import time
t0 = time.time()

There are the following multiple risk factor valuation classes available:

valuation_mcs_european_multi for the valuation of multi-risk derivatives with European exercise
valuation_mcs_american_multi for the valuation of multi-risk derivatives with American exercise

The handling of these classes is similar to building a portfolio of single-risk derivatives positions.

Market Environments

Market environments for the risk factors are the starting point.



In [3]:

    
r = constant_short_rate('r', 0.06)



In [4]:

    
me1 = market_environment('me1', dt.datetime(2015, 1, 1))
me2 = market_environment('me2', dt.datetime(2015, 1, 1))



In [5]:

    
me1.add_constant('initial_value', 36.)
me1.add_constant('volatility', 0.1)  # low volatility
me1.add_constant('currency', 'EUR')
me1.add_constant('model', 'gbm')



In [6]:

    
me2.add_environment(me1)
me2.add_constant('initial_value', 36.)
me2.add_constant('volatility', 0.5)  # high volatility

We assum a positive correlation between the two risk factors.



In [7]:

    
risk_factors = {'gbm1' : me1, 'gbm2' : me2}
correlations = [['gbm1', 'gbm2', 0.5]]

Valuation Environment

Similar to the instantiation of a derivatives_portfolio object, a valuation environment is needed (unifying certain parameters/assumptions for all relevant risk factors of a derivative).



In [8]:

    
val_env = market_environment('val_env', dt.datetime(2015, 1, 1))



In [9]:

    
val_env.add_constant('starting_date', val_env.pricing_date)
val_env.add_constant('final_date', dt.datetime(2015, 12, 31))
val_env.add_constant('frequency', 'W')
val_env.add_constant('paths', 5000)
val_env.add_curve('discount_curve', r)
val_env.add_constant('maturity', dt.datetime(2015, 12, 31))
val_env.add_constant('currency', 'EUR')

valuation_mcs_european_multi

As an example for a multi-risk derivative with European exercise consider a maximum call option. With multiple risk factors, payoff functions are defined by adding key (the name strings) to the maturity_value array object. As with the portfolio valuation class, the multi-risk factor valuation classes get passed market_environment objects only and not the risk factor model objects themsemselves.



In [10]:

    
# European maximum call option
payoff_func = "np.maximum(np.maximum(maturity_value['gbm1'], maturity_value['gbm2']) - 38, 0)"
vc = valuation_mcs_european_multi(
            name='European maximum call',  # name
            val_env=val_env,  # valuation environment
            risk_factors=risk_factors,  # the relevant risk factors
            correlations=correlations,  # correlations between risk factors
            payoff_func=payoff_func)  # payoff function



In [11]:

    
vc.risk_factors









    Out[11]:





{'gbm1': <dx.dx_frame.market_environment at 0x1098ca790>,
 'gbm2': <dx.dx_frame.market_environment at 0x1098cadd0>}

At instantiation, the respective risk factor model objects are instantiated as well.



In [12]:

    
vc.underlying_objects









    Out[12]:





{'gbm1': <dx.dx_models.geometric_brownian_motion at 0x1098cafd0>,
 'gbm2': <dx.dx_models.geometric_brownian_motion at 0x1098e2190>}

Correlations are stored as well and the resulting corrleation and Cholesky matrices are generated.



In [13]:

    
vc.correlations









    Out[13]:





[['gbm1', 'gbm2', 0.5]]



In [14]:

    
vc.correlation_matrix



In [15]:

    
vc.val_env.get_list('cholesky_matrix')









    Out[15]:





array([[ 1.       ,  0.       ],
       [ 0.5      ,  0.8660254]])

The payoff of a European option is a one-dimensional ndarray object.



In [16]:

    
np.shape(vc.generate_payoff())









    Out[16]:





(5000,)

Present value estimations are generated by a call of the present_value method.



In [17]:

    
vc.present_value()









    Out[17]:





7.599

The update method allows updating of certain parameters.



In [18]:

    
vc.update('gbm1', initial_value=50.)



In [19]:

    
vc.present_value()









    Out[19]:





16.671



In [20]:

    
vc.update('gbm2', volatility=0.6)



In [21]:

    
vc.present_value()









    Out[21]:





17.886

Let us reset the values to the original parameters.



In [22]:

    
vc.update('gbm1', initial_value=36., volatility=0.1)
vc.update('gbm2', initial_value=36., volatility=0.5)

When calculating Greeks the risk factor now has to be specified by providing its name.



In [23]:

    
vc.delta('gbm2', interval=0.5)









    Out[23]:





0.5679999999999996



In [24]:

    
vc.vega('gbm1')









    Out[24]:





6.099999999999994

Sensitivities Positive Correlation

Almos in complete analogy to the single-risk valuation classes, sensitivities can be estimated for the multi-risk valuation classes.

Sensitivities Risk Factor 1

Consider first the case from before with positive correlation between the two risk factors. The following estimates and plots the sensitivities for the first risk factor gbm1.



In [25]:

    
%%time
s_list = np.arange(28., 46.1, 2.)
pv = []; de = []; ve = []
for s in s_list:
    vc.update('gbm1', initial_value=s)
    pv.append(vc.present_value())
    de.append(vc.delta('gbm1', .5))
    ve.append(vc.vega('gbm1', 0.2))
vc.update('gbm1', initial_value=36.)









    



CPU times: user 434 ms, sys: 17.9 ms, total: 452 ms
Wall time: 457 ms



In [26]:

    
%matplotlib inline



In [27]:

    
plot_option_stats(s_list, pv, de, ve)

Sensitivities Risk Factor 2

Now the sensitivities for the second risk factor.



In [28]:

    
%%time
s_list = np.arange(28., 46.1, 2.)
pv = []; de = []; ve = []
for s in s_list:
    vc.update('gbm2', initial_value=s)
    pv.append(vc.present_value())
    de.append(vc.delta('gbm2', .5))
    ve.append(vc.vega('gbm2', 0.2))









    



CPU times: user 379 ms, sys: 18.4 ms, total: 398 ms
Wall time: 399 ms



In [29]:

    
plot_option_stats(s_list, pv, de, ve)

Sensitivities with Negative Correlation

The second case is for highly negatively correlated risk factors.



In [30]:

    
correlations = [['gbm1', 'gbm2', -0.9]]



In [31]:

    
# European maximum call option
payoff_func = "np.maximum(np.maximum(maturity_value['gbm1'], maturity_value['gbm2']) - 38, 0)"
vc = valuation_mcs_european_multi(
            name='European maximum call',
            val_env=val_env,
            risk_factors=risk_factors,
            correlations=correlations,
            payoff_func=payoff_func)

Sensitivities Risk Factor 1

Again, sensitivities for the first risk factor first.



In [32]:

    
%%time
s_list = np.arange(28., 46.1, 2.)
pv = []; de = []; ve = []
for s in s_list:
    vc.update('gbm1', initial_value=s)
    pv.append(vc.present_value())
    de.append(vc.delta('gbm1', .5))
    ve.append(vc.vega('gbm1', 0.2))
vc.update('gbm1', initial_value=36.)









    



CPU times: user 389 ms, sys: 17.8 ms, total: 407 ms
Wall time: 408 ms



In [33]:

    
plot_option_stats(s_list, pv, de, ve)

Sensitivities Risk Factor 2

Finally, the sensitivities for the second risk factor for this second scenario.



In [34]:

    
%%time
s_list = np.arange(28., 46.1, 2.)
pv = []; de = []; ve = []
for s in s_list:
    vc.update('gbm2', initial_value=s)
    pv.append(vc.present_value())
    de.append(vc.delta('gbm2', .5))
    ve.append(vc.vega('gbm2', 0.2))









    



CPU times: user 385 ms, sys: 29.9 ms, total: 415 ms
Wall time: 419 ms



In [35]:

    
plot_option_stats(s_list, pv, de, ve)

Surfaces for Positive Correlation Case

Let us return to the case of positive correlation between the two relevant risk factors.



In [36]:

    
correlations = [['gbm1', 'gbm2', 0.5]]



In [37]:

    
# European maximum call option
payoff_func = "np.maximum(np.maximum(maturity_value['gbm1'], maturity_value['gbm2']) - 38, 0)"
vc = valuation_mcs_european_multi(
            name='European maximum call',
            val_env=val_env,
            risk_factors=risk_factors,
            correlations=correlations,
            payoff_func=payoff_func)

Value Surface

We are now interested in the value surface of the derivative instrument for both different initial values of the first and second risk factor.



In [38]:

    
asset_1 = np.arange(28., 46.1, 4.)  # range of initial values
asset_2 = asset_1
a_1, a_2 = np.meshgrid(asset_1, asset_2)
  # two-dimensional grids out of the value vectors
value = np.zeros_like(a_1)

The following estimates for all possible combinations of the initial values---given the assumptions from above---the present value of the European maximum call option.



In [39]:

    
%%time
for i in range(np.shape(value)[0]):
    for j in range(np.shape(value)[1]):
        vc.update('gbm1', initial_value=a_1[i, j])
        vc.update('gbm2', initial_value=a_2[i, j])
        value[i, j] = vc.present_value()









    



CPU times: user 450 ms, sys: 21.8 ms, total: 471 ms
Wall time: 472 ms

The resulting plot then looks as follows. Here, a helper plot function of DX Analytics is used.



In [40]:

    
plot_greeks_3d([a_1, a_2, value], ['gbm1', 'gbm2', 'present value'])

Delta Surfaces

Applying a very similar approach, a delta surface for all possible combinations of the intial values is as easily generated.



In [41]:

    
delta_1 = np.zeros_like(a_1)
delta_2 = np.zeros_like(a_1)



In [42]:

    
%%time
for i in range(np.shape(delta_1)[0]):
    for j in range(np.shape(delta_1)[1]):
        vc.update('gbm1', initial_value=a_1[i, j])
        vc.update('gbm2', initial_value=a_2[i, j])
        delta_1[i, j] = vc.delta('gbm1')
        delta_2[i, j] = vc.delta('gbm2')









    



CPU times: user 1.49 s, sys: 70.9 ms, total: 1.56 s
Wall time: 1.6 s

The plot for the delta surface of the first risk factor.



In [43]:

    
plot_greeks_3d([a_1, a_2, delta_1], ['gbm1', 'gbm2', 'delta gbm1'])

And the plot for the delta of the second risk factor.



In [44]:

    
plot_greeks_3d([a_1, a_2, delta_2], ['gbm1', 'gbm2', 'delta gbm2'])

Vega Surfaces

The same approach can of course be applied to generate vega surfaces.



In [45]:

    
vega_1 = np.zeros_like(a_1)
vega_2 = np.zeros_like(a_1)



In [46]:

    
for i in range(np.shape(vega_1)[0]):
    for j in range(np.shape(vega_1)[1]):
        vc.update('gbm1', initial_value=a_1[i, j])
        vc.update('gbm2', initial_value=a_2[i, j])
        vega_1[i, j] = vc.vega('gbm1')
        vega_2[i, j] = vc.vega('gbm2')

The surface for the first risk factor.



In [47]:

    
plot_greeks_3d([a_1, a_2, vega_1], ['gbm1', 'gbm2', 'vega gbm1'])

And the one for the second risk factor.



In [48]:

    
plot_greeks_3d([a_1, a_2, vega_2], ['gbm1', 'gbm2', 'vega gbm2'])

Finally, we reset the intial values and the volatilities for the two risk factors.



In [49]:

    
# restore initial values
vc.update('gbm1', initial_value=36., volatility=0.1)
vc.update('gbm2', initial_value=36., volatility=0.5)

valuation_mcs_american_multi

In general, the modeling and handling of the valuation classes for American exercise is not too different from those for European exercise. The major difference is in the definition of payoff function.

Present Values

This example models an American minimum put on the two risk factors from before.



In [50]:

    
# American put payoff
payoff_am = "np.maximum(34 - np.minimum(instrument_values['gbm1'], instrument_values['gbm2']), 0)"
# finer time grid and more paths
val_env.add_constant('frequency', 'B')
val_env.add_curve('time_grid', None)
  # delete existing time grid information
val_env.add_constant('paths', 5000)



In [51]:

    
# American put option on minimum of two assets
vca = valuation_mcs_american_multi(
            name='American minimum put',
            val_env=val_env,
            risk_factors=risk_factors,
            correlations=correlations,
            payoff_func=payoff_am)



In [52]:

    
vca.present_value()









    Out[52]:





4.601



In [53]:

    
for key, obj in vca.instrument_values.items():
    print np.shape(vca.instrument_values[key])









    



(261, 5000)
(261, 5000)

The present value surface is generated in the same way as before for the European option on the two risk factors. The computational burden is of course much higher for the American option, which are valued by the use of the Least-Squares Monte Carlo approach (LSM) according to Longstaff-Schwartz (2001).



In [54]:

    
asset_1 = np.arange(28., 44.1, 4.)
asset_2 = asset_1
a_1, a_2 = np.meshgrid(asset_1, asset_2)
value = np.zeros_like(a_1)



In [55]:

    
%%time
for i in range(np.shape(value)[0]):
    for j in range(np.shape(value)[1]):
        vca.update('gbm1', initial_value=a_1[i, j])
        vca.update('gbm2', initial_value=a_2[i, j])
        value[i, j] = vca.present_value()









    



CPU times: user 21.2 s, sys: 1.35 s, total: 22.5 s
Wall time: 22.3 s



In [56]:

    
plot_greeks_3d([a_1, a_2, value], ['gbm1', 'gbm2', 'present value'])

Delta Surfaces

The same exercise as before for the two delta surfaces.



In [57]:

    
delta_1 = np.zeros_like(a_1)
delta_2 = np.zeros_like(a_1)



In [58]:

    
%%time
for i in range(np.shape(delta_1)[0]):
    for j in range(np.shape(delta_1)[1]):
        vca.update('gbm1', initial_value=a_1[i, j])
        vca.update('gbm2', initial_value=a_2[i, j])
        delta_1[i, j] = vca.delta('gbm1')
        delta_2[i, j] = vca.delta('gbm2')









    



CPU times: user 1min 14s, sys: 4.67 s, total: 1min 19s
Wall time: 1min 18s



In [59]:

    
plot_greeks_3d([a_1, a_2, delta_1], ['gbm1', 'gbm2', 'delta gbm1'])



In [60]:

    
plot_greeks_3d([a_1, a_2, delta_2], ['gbm1', 'gbm2', 'delta gbm2'])

Vega Surfaces

And finally for the vega surfaces.



In [61]:

    
vega_1 = np.zeros_like(a_1)
vega_2 = np.zeros_like(a_1)



In [62]:

    
%%time
for i in range(np.shape(vega_1)[0]):
    for j in range(np.shape(vega_1)[1]):
        vca.update('gbm1', initial_value=a_1[i, j])
        vca.update('gbm2', initial_value=a_2[i, j])
        vega_1[i, j] = vca.vega('gbm1')
        vega_2[i, j] = vca.vega('gbm2')









    



CPU times: user 1min 19s, sys: 4.81 s, total: 1min 23s
Wall time: 1min 24s



In [63]:

    
plot_greeks_3d([a_1, a_2, vega_1], ['gbm1', 'gbm2', 'vega gbm1'])



In [64]:

    
plot_greeks_3d([a_1, a_2, vega_2], ['gbm1', 'gbm2', 'vega gbm2'])

More than Two Risk Factors

The principles of working with multi-risk valuation classes can be illustrated quite well in the two risk factor case. However, there is---in theory---no limitation on the number of risk factors used for derivatives modeling.

Four Asset Basket Option

Consider a maximum basket option on four different risk factors. We add a jump diffusion as well as a stochastic volatility model to the mix



In [65]:

    
me3 = market_environment('me3', dt.datetime(2015, 1, 1))
me4 = market_environment('me4', dt.datetime(2015, 1, 1))



In [66]:

    
me3.add_environment(me1)
me4.add_environment(me1)



In [67]:

    
# for jump-diffusion
me3.add_constant('lambda', 0.5)
me3.add_constant('mu', -0.6)
me3.add_constant('delta', 0.1)
me3.add_constant('model', 'jd')



In [68]:

    
# for stoch volatility model
me4.add_constant('kappa', 2.0)
me4.add_constant('theta', 0.3)
me4.add_constant('vol_vol', 0.2)
me4.add_constant('rho', -0.75)
me4.add_constant('model', 'sv')



In [69]:

    
val_env.add_constant('paths', 5000)
val_env.add_constant('frequency', 'W')
val_env.add_curve('time_grid', None)

In this case, we need to specify three correlation values.



In [70]:

    
risk_factors = {'gbm1' : me1, 'gbm2' : me2, 'jd' : me3, 'sv' : me4}
correlations = [['gbm1', 'gbm2', 0.5], ['gbm2', 'jd', -0.5], ['gbm1', 'sv', 0.7]]

The payoff function in this case gets a bit more complex.



In [71]:

    
# European maximum call payoff
payoff_1 = "np.maximum(np.maximum(np.maximum(maturity_value['gbm1'], maturity_value['gbm2']),"
payoff_2 = " np.maximum(maturity_value['jd'], maturity_value['sv'])) - 40, 0)"
payoff = payoff_1 + payoff_2



In [72]:

    
payoff









    Out[72]:





"np.maximum(np.maximum(np.maximum(maturity_value['gbm1'], maturity_value['gbm2']), np.maximum(maturity_value['jd'], maturity_value['sv'])) - 40, 0)"

However, the instantiation of the valuation classe remains the same.



In [73]:

    
vc = valuation_mcs_european_multi(
            name='European maximum call',
            val_env=val_env,
            risk_factors=risk_factors,
            correlations=correlations,
            payoff_func=payoff)

Example Output and Calculations

The following just displays some example output and the results from certain calculations.



In [74]:

    
vc.risk_factors









    Out[74]:





{'gbm1': <dx.dx_frame.market_environment at 0x1098ca790>,
 'gbm2': <dx.dx_frame.market_environment at 0x1098cadd0>,
 'jd': <dx.dx_frame.market_environment at 0x10ac41d50>,
 'sv': <dx.dx_frame.market_environment at 0x10218eb90>}



In [75]:

    
vc.underlying_objects









    Out[75]:





{'gbm1': <dx.dx_models.geometric_brownian_motion at 0x10ac51e90>,
 'gbm2': <dx.dx_models.geometric_brownian_motion at 0x10c335050>,
 'jd': <dx.dx_models.jump_diffusion at 0x10acbb3d0>,
 'sv': <dx.dx_models.stochastic_volatility at 0x10c335610>}



In [76]:

    
vc.present_value()









    Out[76]:





13.171

The correlation and Cholesky matrices now are of shape 4x4.



In [77]:

    
vc.correlation_matrix



In [78]:

    
vc.val_env.get_list('cholesky_matrix')









    Out[78]:





array([[ 1.        ,  0.        ,  0.        ,  0.        ],
       [ 0.5       ,  0.8660254 ,  0.        ,  0.        ],
       [ 0.        , -0.57735027,  0.81649658,  0.        ],
       [ 0.7       , -0.40414519, -0.2857738 ,  0.51478151]])

Delta and vega estimates are generated in exactly the same fashion as in the two risk factor case.



In [79]:

    
vc.delta('jd', interval=0.1)









    Out[79]:





0.4499999999999993



In [80]:

    
vc.delta('sv')









    Out[80]:





0.37083333333333507



In [81]:

    
vc.vega('jd')









    Out[81]:





6.400000000000006



In [82]:

    
vc.vega('sv')









    Out[82]:





1.29999999999999

Delta for Jump Diffusion and Stochastic Vol Process

Of course, we cannot visualize Greek surfaces dependent on initial values for all four risk factors but still for two. In what follows we generate the delta surfaces with respect to the jump diffusion- and stochastic volatility-based risk factors.



In [83]:

    
delta_1 = np.zeros_like(a_1)
delta_2 = np.zeros_like(a_1)



In [84]:

    
%%time
for i in range(np.shape(delta_1)[0]):
    for j in range(np.shape(delta_1)[1]):
        vc.update('jd', initial_value=a_1[i, j])
        vc.update('sv', initial_value=a_2[i, j])
        delta_1[i, j] = vc.delta('jd')
        delta_2[i, j] = vc.delta('sv')









    



CPU times: user 6.68 s, sys: 921 ms, total: 7.6 s
Wall time: 7.62 s



In [85]:

    
plot_greeks_3d([a_1, a_2, delta_1], ['jump diffusion', 'stochastic vol', 'delta jd'])



In [86]:

    
plot_greeks_3d([a_1, a_2, delta_2], ['jump diffusion', 'stochastic vol', 'delta sv'])

Vega for Jump Diffusion and Stochastic Vol Process

Now the same exercise for the vega surfaces for the same two risk factors.



In [87]:

    
vega_1 = np.zeros_like(a_1)
vega_2 = np.zeros_like(a_1)



In [88]:

    
%%time
for i in range(np.shape(vega_1)[0]):
    for j in range(np.shape(vega_1)[1]):
        vc.update('jd', initial_value=a_1[i, j])
        vc.update('sv', initial_value=a_2[i, j])
        vega_1[i, j] = vc.vega('jd')
        vega_2[i, j] = vc.vega('sv')









    



CPU times: user 7.3 s, sys: 1.04 s, total: 8.34 s
Wall time: 8.42 s



In [89]:

    
plot_greeks_3d([a_1, a_2, vega_1], ['jump diffusion', 'stochastic vol', 'vega jd'])



In [90]:

    
plot_greeks_3d([a_1, a_2, vega_2], ['jump diffusion', 'stochastic vol', 'vega sv'])

American Exercise

As a final illustration consider the case of an American minimum put option on the four risk factors. This again is a step that leads to a much increased computational burden due to the necessity to apply the least-squares regression approach.



In [91]:

    
# payoff of American minimum put option
payoff_am_1 = "np.maximum(40 - np.minimum(np.minimum(instrument_values['gbm1'], instrument_values['gbm2']),"
payoff_am_2 = "np.minimum(instrument_values['jd'], instrument_values['sv'])), 0)"
payoff_am = payoff_am_1 + payoff_am_2



In [92]:

    
vca = valuation_mcs_american_multi(
            name='American minimum put',
            val_env=val_env,
            risk_factors=risk_factors,
            correlations=correlations,
            payoff_func=payoff_am)

However, another illustration that even such a complex instrument can be handled as elegantly as the most simple one (i.e. European option on single risk factor). Let us compare the present value estimates for both the European and American maximum basket options.



In [93]:

    
# restore initial values
vc.update('jd', initial_value=36., volatility=0.1)
vc.update('sv', initial_value=36., volatility=0.1)
%time vc.present_value()









    



CPU times: user 2.91 ms, sys: 0 ns, total: 2.91 ms
Wall time: 2.92 ms






    Out[93]:





13.171



In [94]:

    
%time vca.present_value()









    



CPU times: user 791 ms, sys: 154 ms, total: 945 ms
Wall time: 886 ms






    Out[94]:





14.444



In [95]:

    
%time vca.delta('gbm1')









    



CPU times: user 1.33 s, sys: 302 ms, total: 1.63 s
Wall time: 1.5 s






    Out[95]:





-0.0027777777777787055



In [96]:

    
%time vca.delta('gbm2')









    



CPU times: user 1.39 s, sys: 305 ms, total: 1.69 s
Wall time: 1.63 s






    Out[96]:





-0.21250000000000066



In [97]:

    
%time vca.vega('jd')









    



CPU times: user 1.35 s, sys: 244 ms, total: 1.59 s
Wall time: 1.49 s






    Out[97]:





3.0999999999998806



In [98]:

    
%time vca.vega('sv')









    



CPU times: user 1.35 s, sys: 224 ms, total: 1.58 s
Wall time: 1.44 s






    Out[98]:





1.3999999999999346



In [99]:

    
print "Duration for whole notebook %.2f in min" % ((time.time() - t0) / 60)









    



Duration for whole notebook 3.77 in min

Copyright, License & Disclaimer

DX Analytics (the "dx library") is licensed under the GNU Affero General Public License version 3 or later (see http://www.gnu.org/licenses/).

DX Analytics comes with no representations or warranties, to the extent permitted by applicable law.

http://www.pythonquants.com | analytics@pythonquants.com | http://twitter.com/dyjh

Python Quant Platform | http://quant-platform.com

Derivatives Analytics with Python (Wiley Finance) | http://derivatives-analytics-with-python.com

Python for Finance (O'Reilly) | http://shop.oreilly.com/product/0636920032441.do

	gbm1	gbm2	jd	sv
gbm1	1.0	0.5	0.0	0.7
gbm2	0.5	1.0	-0.5	0.0
jd	0.0	-0.5	1.0	0.0
sv	0.7	0.0	0.0	1.0