optSensByMertonNI

Option price and sensitivities by Merton76 model using numerical integration

Syntax

PriceSens = optSensByMertonNI(Rate,AssetPrice,Settle,Maturity,OptSpec,Strike,Sigma,MeanJ,JumpVol,JumpFreq)

PriceSens = optSensByMertonNI(___,Name,Value)

Description

PriceSens = optSensByMertonNI(Rate,AssetPrice,Settle,Maturity,OptSpec,Strike,Sigma,MeanJ,JumpVol,JumpFreq) computes vanilla European option price and sensitivities by the Merton76 model, using numerical integration.

Note

Alternatively, you can use the Vanilla object to calculate price or sensitivities for vanilla options. For more information, see Get Started with Workflows Using Object-Based Framework for Pricing Financial Instruments.

PriceSens = optSensByMertonNI(___,Name,Value) adds optional name-value pair arguments.

Examples

Workflow for Plotting an Option Sensitivity Surface Using the Merton76 Model

optSensByMertonNI uses numerical integration to compute option sensitivities and then plot option sensitivity surfaces.

Define Option Variables and Merton76 Model Parameters

AssetPrice = 80;
Rate = 0.03;
DividendYield = 0.02;
OptSpec = 'call';

Sigma = 0.16;
MeanJ = 0.02;
JumpVol = 0.08;
JumpFreq = 2;

Compute the Option Sensitivity for a Single Strike

Settle = datetime(2017,6,29);
Maturity = datemnth(Settle, 6);
Strike = 80; 

Delta = optSensByMertonNI(Rate, AssetPrice, Settle, Maturity, OptSpec, Strike, ...
    Sigma, MeanJ, JumpVol, JumpFreq, 'DividendYield', DividendYield, ...
    'OutSpec', "delta")

Delta = 
0.5283

Compute the Option Sensitivities for a Vector of Strikes

The Strike input can be a vector.

Settle = datetime(2017,6,29);
Maturity = datemnth(Settle, 6);
Strike = (76:2:84)';

Delta = optSensByMertonNI(Rate, AssetPrice, Settle, Maturity, OptSpec, Strike, ...
    Sigma, MeanJ, JumpVol, JumpFreq, 'DividendYield', DividendYield, ...
    'OutSpec', "delta")

Delta = 5×1

    0.6727
    0.6013
    0.5283
    0.4565
    0.3883

Compute the Option Sensitivities for a Vector of Strikes and a Vector of Dates of the Same Lengths

Use the Strike input to specify the strikes. Also, the Maturity input can be a vector, but it must match the length of the Strike vector if the ExpandOutput name-value pair argument is not set to "true".

Settle = datetime(2017,6,29);
Maturity = datemnth(Settle, [12 18 24 30 36]); % Five maturities
Strike = [76 78 80 82 84]'; % Five strikes

Delta = optSensByMertonNI(Rate, AssetPrice, Settle, Maturity, OptSpec, Strike, ...
    Sigma, MeanJ, JumpVol, JumpFreq, 'DividendYield', DividendYield, ...
    'OutSpec', "delta")

Delta = 5×1

    0.6419
    0.5907
    0.5565
    0.5311
    0.5110

    % Five values in vector output

Expand the Output for a Surface

Set the ExpandOutput name-value pair argument to "true" to expand the output into a NStrikes-by-NMaturities matrix. In this case, it is a square matrix.

Delta = optSensByMertonNI(Rate, AssetPrice, Settle, Maturity, OptSpec, Strike, ...
    Sigma, MeanJ, JumpVol, JumpFreq, 'DividendYield', DividendYield, ...
    'OutSpec', "delta", 'ExpandOutput', true) % (5 x 5) matrix output

Delta = 5×5

    0.6419    0.6305    0.6245    0.6204    0.6173
    0.5922    0.5907    0.5905    0.5905    0.5905
    0.5422    0.5507    0.5565    0.5607    0.5637
    0.4927    0.5112    0.5229    0.5311    0.5372
    0.4447    0.4725    0.4898    0.5020    0.5110

Compute the Option Sensitivities for a Vector of Strikes and a Vector of Dates of Different Lengths

When ExpandOutput is "true", NStrikes do not have to match NMaturities. That is, the output NStrikes-by-NMaturities matrix can be rectangular.

Settle = datetime(2017,6,29);
Maturity = datemnth(Settle, 12*(0.5:0.5:3)'); % Six maturities
Strike = (76:2:84)'; % Five strikes

Delta = optSensByMertonNI(Rate, AssetPrice, Settle, Maturity, OptSpec, Strike, ...
    Sigma, MeanJ, JumpVol, JumpFreq, 'DividendYield', DividendYield, ...
    'OutSpec', "delta", 'ExpandOutput', true)  % (5 x 6) matrix output

Delta = 5×6

    0.6727    0.6419    0.6305    0.6245    0.6204    0.6173
    0.6013    0.5922    0.5907    0.5905    0.5905    0.5905
    0.5283    0.5422    0.5507    0.5565    0.5607    0.5637
    0.4565    0.4927    0.5112    0.5229    0.5311    0.5372
    0.3883    0.4447    0.4725    0.4898    0.5020    0.5110

Compute the Option Sensitivities for a Vector of Strikes and a Vector of Asset Prices

When ExpandOutput is "true", the output can also be a NStrikes-by-NAssetPrices rectangular matrix by accepting a vector of asset prices.

Settle = datetime(2017,6,29);
Maturity = datemnth(Settle, 12); % Single maturity
ManyAssetPrices = [70 75 80 85]; % Four asset prices
Strike = (76:2:84)'; % Five strikes

Delta = optSensByMertonNI(Rate, ManyAssetPrices, Settle, Maturity, OptSpec, Strike, ...
    Sigma, MeanJ, JumpVol, JumpFreq, 'DividendYield', DividendYield, ...
    'OutSpec', "delta", 'ExpandOutput', true) % (5 x 4) matrix output

Delta = 5×4

    0.3796    0.5157    0.6419    0.7472
    0.3315    0.4637    0.5922    0.7043
    0.2874    0.4137    0.5422    0.6592
    0.2474    0.3664    0.4927    0.6128
    0.2117    0.3224    0.4447    0.5657

Plot Option Sensitivity Surfaces

The Strike and Maturity inputs can be vectors. Set ExpandOutput to "true" to output the surfaces as NStrikes-by-NMaturities matrices.

Settle = datetime(2017,6,29);
Maturity = datemnth(Settle, 12*[1/12 0.25 (0.5:0.5:3)]');
Times = yearfrac(Settle, Maturity);
Strike = (2:2:200)';

[Delta, Gamma, Rho, Theta, Vega] = optSensByMertonNI(...
    Rate, AssetPrice,Settle, Maturity, OptSpec, Strike, ...
    Sigma, MeanJ, JumpVol, JumpFreq, 'DividendYield', DividendYield, ...
    'OutSpec', ["delta", "gamma", "rho", "theta", "vega"], ...
    'ExpandOutput', true);

[X,Y] = meshgrid(Times,Strike);

figure;
surf(X,Y,Delta);
title('Delta');
xlabel('Years to Option Expiry');
ylabel('Strike');
view(-112,34);
xlim([0 Times(end)]);

figure;
surf(X,Y,Gamma)
title('Gamma')
xlabel('Years to Option Expiry')
ylabel('Strike')
view(-112,34);
xlim([0 Times(end)]);

figure;
surf(X,Y,Rho)
title('Rho')
xlabel('Years to Option Expiry')
ylabel('Strike')
view(-112,34);
xlim([0 Times(end)]);

figure;
surf(X,Y,Theta)
title('Theta')
xlabel('Years to Option Expiry')
ylabel('Strike')
view(-112,34);
xlim([0 Times(end)]);

figure;
surf(X,Y,Vega)
title('Vega')
xlabel('Years to Option Expiry')
ylabel('Strike')
view(-112,34);
xlim([0 Times(end)]);

Input Arguments

Rate — Continuously compounded risk-free interest rate

decimal

Continuously compounded risk-free interest rate, specified as a scalar decimal value.

Data Types: double

AssetPrice — Current underlying asset price

numeric

Current underlying asset price, specified as numeric value using a scalar or a NINST-by-1 or NColumns-by-1 vector.

For more information on the proper dimensions for AssetPrice, see the name-value pair argument ExpandOutput.

Data Types: double

Settle — Option settlement date

datetime array | string array | date character vector

Option settlement date, specified as a NINST-by-1 or NColumns-by-1 vector using a datetime array, string array, or date character vectors. The Settle date must be before the Maturity date.

To support existing code, optSensByMertonNI also accepts serial date numbers as inputs, but they are not recommended.

For more information on the proper dimensions for Settle, see the name-value pair argument ExpandOutput.

Maturity — Option maturity date

datetime array | string array | date character vector

Option maturity date, specified as a NINST-by-1 or NColumns-by-1 vector using a datetime array, string array, or date character vectors.

To support existing code, optSensByMertonNI also accepts serial date numbers as inputs, but they are not recommended.

For more information on the proper dimensions for Maturity, see the name-value pair argument ExpandOutput.

OptSpec — Definition of option

cell array of character vector with values 'call' or 'put' | string array with values "call" or "put"

Definition of the option, specified as a NINST-by-1 or NColumns-by-1 vector using a cell array of character vectors or string arrays with values 'call' or 'put'.

For more information on the proper dimensions for OptSpec, see the name-value pair argument ExpandOutput.

Data Types: cell | string

Strike — Option strike price value

numeric

Option strike price value, specified as a NINST-by-1, NRows-by-1, NRows-by-NColumns vector of strike prices.

For more information on the proper dimensions for Strike, see the name-value pair argument ExpandOutput.

Data Types: double

Sigma — Volatility of underlying asset

numeric

Volatility of the underlying asset, specified as a scalar numeric value.

Data Types: double

MeanJ — Mean of the random percentage jump size

decimal

Mean of the random percentage jump size (J), specified as a scalar decimal value where log(1+J) is normally distributed with mean (log(1+MeanJ)-0.5*JumpVol^2) and the standard deviation JumpVol.

Data Types: double

JumpVol — Standard deviation of log(1+J)

decimal

Standard deviation of log(1+J) where J is the random percentage jump size, specified as a scalar decimal value.

Data Types: double

JumpFreq — Annual frequency of Poisson jump process

numeric

Annual frequency of Poisson jump process, specified as a scalar numeric value.

Data Types: double

Name-Value Arguments

Specify optional pairs of arguments as Name1=Value1,...,NameN=ValueN, where Name is the argument name and Value is the corresponding value. Name-value arguments must appear after other arguments, but the order of the pairs does not matter.

Before R2021a, use commas to separate each name and value, and enclose Name in quotes.

Example: Price = optSensByMertonNI(Rate,AssetPrice,Settle,Maturity,OptSpec,Strike,Sigma,MeanJ,JumpVol,JumpFreq,'Basis',7)

Basis — Day-count basis of instrument

0 (default) | numeric values: 0,1, 2, 3, 4, 6, 7, 8, 9, 10, 11, 12, 13

Day-count of the instrument, specified as the comma-separated pair consisting of 'Basis' and a scalar using a supported value:

• 0 = actual/actual

• 1 = 30/360 (SIA)

• 2 = actual/360

• 3 = actual/365

• 4 = 30/360 (PSA)

• 5 = 30/360 (ISDA)

• 6 = 30/360 (European)

• 7 = actual/365 (Japanese)

• 8 = actual/actual (ICMA)

• 9 = actual/360 (ICMA)

• 10 = actual/365 (ICMA)

• 11 = 30/360E (ICMA)

• 12 = actual/365 (ISDA)

• 13 = BUS/252

For more information, see Basis.

Data Types: double

DividendYield — Continuously compounded underlying asset yield

0 (default) | numeric

Continuously compounded underlying asset yield, specified as the comma-separated pair consisting of 'DividendYield' and a scalar numeric value.

Data Types: double

OutSpec — Define outputs

["price"] (default) | string array with values "price", "delta", "gamma", "vega", "rho", and "theta" | cell array of character vectors with values 'price', 'delta', 'gamma', 'vega', 'rho', and 'theta'

Define outputs, specified as the comma-separated pair consisting of 'OutSpec' and a NOUT- by-1 or a 1-by-NOUT string array or cell array of character vectors with supported values.

Note

"vega" is the sensitivity with respect the initial volatility sqrt(V0).

Example: OutSpec = ["price","delta","gamma","vega","rho","theta"]

Data Types: string | cell

AbsTol — Absolute error tolerance for numerical integration

1e-10 (default) | numeric

Absolute error tolerance for numerical integration, specified as the comma-separated pair consisting of 'AbsTol' and a scalar numeric value.

Data Types: double

RelTol — Relative error tolerance for numerical integration

1e-6 (default) | numeric

Relative error tolerance for numerical integration, specified as the comma-separated pair consisting of 'RelTol' and a scalar numeric value.

Data Types: double

IntegrationRange — Numerical integration range used to approximate continuous integral over [0 Inf]

[1e-9 Inf] (default) | vector

Numerical integration range used to approximate the continuous integral over [0 Inf], specified as the comma-separated pair consisting of 'IntegrationRange' and a 1-by-2 vector representing [LowerLimit UpperLimit].

Data Types: double

Framework — Framework for computing option prices and sensitivities using numerical integration of models

"heston1993" (default) | string with values "heston1993" or "lewis2001" | character vector with values 'heston1993' or 'lewis2001'

Framework for computing option prices and sensitivities using numerical integration of models, specified as the comma-separated pair consisting of 'Framework' and a scalar string or character vector with the following values:

• "heston1993" or 'heston1993' — Method used in Heston (1993)

• "lewis2001" or 'lewis2001' — Method used in Lewis (2001)

Data Types: char | string

ExpandOutput — Flag to expand the outputs

false (outputs are NINST-by-1 vectors) (default) | logical with value of true or false

Flag to expand the outputs, specified as the comma-separated pair consisting of 'ExpandOutput' and a logical:

• true — If true, the outputs are NRows-by- NColumns matrices. NRows is the number of strikes for each column and it is determined by the Strike input. For example, Strike can be a NRows-by-1 vector, or a NRows-by-NColumns matrix. NColumns is determined by the sizes of AssetPrice, Settle, Maturity, and OptSpec, which must all be either scalar or NColumns-by-1 vectors.

• false — If false, the outputs are NINST-by-1 vectors. Also, the inputs Strike, AssetPrice, Settle, Maturity, and OptSpec must all be either scalar or NINST-by-1 vectors.

Data Types: logical

Output Arguments

PriceSens — Option prices or sensitivities

numeric

Option prices or sensitivities, returned as a NINST-by-1, or NRows-by-NColumns, depending on ExpandOutput. The name-value pair argument OutSpec determines the types and order of the outputs.

More About

Vanilla Option

A vanilla option is a category of options that includes only the most standard components.

A vanilla option has an expiration date and straightforward strike price. American-style options and European-style options are both categorized as vanilla options.

The payoff for a vanilla option is as follows:

• For a call: $$\max (St-K,0)$$

• For a put: $$\max (K-St,0)$$

where:

St is the price of the underlying asset at time t.

K is the strike price.

For more information, see Vanilla Option.

Merton Jump Diffusion Model

The Merton jump diffusion model (Merton (1976)) is a different extension of the Black-Scholes model, where sudden asset price movements (both up and down) are modeled by adding the jump diffusion parameters with the Poisson process.

The stochastic differential equation is:

$$\begin{array}{l}d{S}_t=(r-q-{\lambda }_p{\mu }_j)S_{t}dt+\sigma S_{t}d{W}_t+J{S}_{t}d{P}_t\\ \text{prob(}d{P}_t=1)={\lambda }_{p}dt\end{array}$$

where

r is the continuous risk-free rate.

q is the continuous dividend yield.

W_t is the Wiener process.

J is the random percentage jump size conditional on the jump occurring, where ln(1+J) is normally distributed with mean $$\ln (1+{\mu }_J)-\frac{{\delta }^2}{2}$$ and the standard deviation δ, and (1+J) has a lognormal distribution:

$$\frac{1}{(1+J)\delta \sqrt{2\pi }}\exp \left\{{\frac{-\left[\ln (1+J)-\left(\ln (1+{\mu }_J\right)-\frac{{\delta }^2}{2}\right]}{2{\delta }^2}}^2\right\}$$

μ_J is the mean of J for (μ_J > -1).

δ is the standard deviation of ln(1+J) for (δ≥ 0).

ƛ_p is the annual frequency (intensity) of Poisson process P_tfor (ƛ_p ≥ 0).

σ is the volatility of the asset price for (σ > 0).

The characteristic function $$f_{Merton{76}_j}(\varphi )$$ for j = 1 (asset prices measure) and j = 2 (risk-neutral measure) is:

$$\begin{array}{l}{f}_{Merton{76}_j}=f_{B{S}_j}\exp \left({\lambda }_p\tau {\left(1+{\mu }_j\right)}^{m_j+\frac{1}{2}}\left[{\left(1+{\mu }_j\right)}^{i\varphi }{e}^{{\delta }^2\left(m_{j}i\varphi +\frac{{(i\varphi )}^2}{2}\right)}-1\right]-{\lambda }_p\tau {\mu }_{j}i\varphi \right)\\ \text{where for }j=1,2:\\ f_{B{S}_1}(\varphi )=\frac{f_{B{S}_2}\left(\varphi -i\right)}{f_{B{S}_2}\left(-i\right)}\\ f_{BS2}(\varphi )=\exp \left(i\varphi \left[\ln S_t+\left(r-q-\frac{{\sigma }^2}{2}\right)\tau \right]-\frac{{\varphi }^2{\sigma }^2}{2}\tau \right)\\ m_1=\frac{1}{2},m_2=-\frac{1}{2}\end{array}$$

where

ϕ is the characteristic function variable.

τ is the time to maturity (τ = T- t).

i is the unit imaginary number ( i² = -1).

Numerical Integration Method Under Heston (1993) Framework

Numerical integration is used to evaluate the continuous integral for the inverse Fourier transform.

The numerical integration method under the Heston (1993) framework is based on the following expressions:

$$\begin{array}{l}Call(K)=S_t{e}^{-q\tau }{P}_1-K{e}^{-r\tau }{P}_2\\ Put(K)=Call(K)+K{e}^{-r\tau }-S_t{e}^{-q\tau }\\ P_j=\frac{1}{2}+\frac{1}{\pi }{\displaystyle \underset{0}{\overset{\infty }{\int }}\mathrm{Re}\left[\frac{e^{-i\varphi \ln (K)}{f}_j(\varphi )}{i\varphi }\right]d\varphi }\end{array}$$

where

r is the continuous risk-free rate.

q is the continuous dividend yield.

S_t is the asset price at time t.

K is the strike.

τ is time to maturity (τ = T-t).

Call(K) is the call price at strike K.

Put(K) is the put price at strike K.

i is a unit imaginary number (i²= -1).

ϕ is the characteristic function variable.

f_j(ϕ) is the characteristic function for P_j(j = 1,2).

P₁ is the probability of S_t > K under the asset price measure for the model.

P₂ is the probability of S_t > K under the risk-neutral measure for the model.

Where j = 1,2 so that f₁(ϕ) and f₂(ϕ) are the characteristic functions for probabilities P₁ and P₂, respectively.

This framework is chosen with the default value “Heston1993” for the Framework name-value pair argument.

Numerical Integration Method Under Lewis (2001) Framework

Numerical integration is used to evaluate the continuous integral for the inverse Fourier transform.

The numerical integration method under the Lewis (2001) framework is based on the following expressions:

$$\begin{array}{l}Call(k)=S_t{e}^{-q\tau }-\frac{\sqrt{K}{e}^{-\tau t}}{\pi }{\displaystyle \underset{0}{\overset{\infty }{\int }}\mathrm{Re}\left[K^{-iu}{f}_2\left(\varphi =\left(u-\frac{i}{2}\right)\right)\frac{1}{u^2+\frac{1}{4}}\right]}du\\ Put(K)=Call(K)=K{e}^{-\tau t}-S_t{e}^{-q\tau }\end{array}$$

where

r is the continuous risk-free rate.

q is the continuous dividend yield.

S_t is the asset price at time t.

K is the strike.

τ is time to maturity (τ = T-t).

Call(K) is the call price at strike K.

Put(K) is the put price at strike K.

i is a unit imaginary number (i²= -1).

ϕ is the characteristic function variable.

u is the characteristic function variable for integration, where $$\varphi =\left(u-\frac{i}{2}\right)$$.

f₂(ϕ) is the characteristic function for P₂.

P₂ is the probability of S_t > K under the risk-neutral measure for the model.

This framework is chosen with the value “Lewis2001” for the Framework name-value pair argument.

References

[1] Bates, D. S. “Jumps and Stochastic Volatility: Exchange Rate Processes Implicit in Deutsche Mark Options.” The Review of Financial Studies. Vol 9. No. 1. 1996.

[2] Cont, R. and P. Tankov. Financial Modeling with Jump Processes. Chapman & Hall/CRC Press, 2004.

[3] Heston, S. L. “A Closed-Form Solution for Options with Stochastic Volatility with Applications to Bond and Currency Options.” The Review of Financial Studies. Vol 6. No. 2. 1993.

[4] Lewis, A. L. “A Simple Option Formula for General Jump-Diffusion and Other Exponential Levy Processes.” Envision Financial Systems and OptionCity.net, 2001.

[5] Merton, R. “Option Pricing When Underlying Stock Returns are Discontinuous.” Journal of Financial Economics. Vol 3. 1976.

Version History

Introduced in R2018a

R2022b: Serial date numbers not recommended

Although optSensByMertonNI supports serial date numbers, datetime values are recommended instead. The datetime data type provides flexible date and time formats, storage out to nanosecond precision, and properties to account for time zones and daylight saving time.

To convert serial date numbers or text to datetime values, use the datetime function. For example:

t = datetime(738427.656845093,"ConvertFrom","datenum");
y = year(t)

y =

        2021

There are no plans to remove support for serial date number inputs.