Black Scholes Greeks

The Black Scholes Greeks

The Option Greek formulae express the change in the option price with respect to a parameter change taking as fixed all the other inputs. (Haug explores multiple parameter changes at once.) One significant use of Greek measures is to calibrate risk exposure. A market-making financial institution with a portfolio of options, for instance, would want a snap shot of its exposure to asset price, interest rates, dividend fluctuations. It would try to establish impacts of volatility and time decay. In the formulae below, the Greeks merely evaluate change to only one input at a time. In reality, we might expect a conflagration of changes in interest rates and stock prices etc.

VBA code for Black Scholes Greeks

No discussion of Black Scholes is complete without estimating the Greeks. For a more exhaustive list of Greek formulae that incorporate dividends please go to wikipedia. To understand the intuition behind Greeks and parameter sensitivities please view this playlist. Parameters are represented by the notation S, K, T, r, q and Sigma. See VBA code below. See Martin Haugh notes for formulae and clear exposition of same. Parameter Sensitivity can be estimated more precisely using the Black Scholes Greeks. The Greeks are important because we rely upon these to assist in hedging. The VBA code below can be used to estimate the effect of altering parameter inputs on the value of options. These are important for hedging and grasping the importance of variables changing. I have adapted here the coding style set by the late Simon Benninga. For a deep dive into estimating the Greeks using both analytical and numerically estimation please check out this playlist graphing included. Below we implement Analytical Greeks in VBA. Also see https://github.com/jrvarma/MSExcel-Black-Scholes

'S = Asset Price

'K = Option Exercise

'T = Option Maturity

'r = Risk free rate

'q = Dividend yield

'Sigma = Asset price return volatility

'Calculating d1

Function d1One(S, K, T, r, q, Sigma)

d1One = (Log(S / K) + (r - q) * T) / (Sigma * Sqr(T)) + 0.5 * Sigma * Sqr(T)

End Function

'Calculating d2

Function d2Two(S, K, T, r, q, Sigma)

d2Two = d1One(S, K, T, r, q, Sigma) - Sigma * Sqr(T)

End Function

'The normal probability density, N'(x), is defined as pdf

Function pdf(x)

pdf = Exp(-x ^ 2 / 2) / (Sqr(2 * Application.Pi()))

End Function

'The Black Schole Merton Call Function

' Application.NormSDist() calls up the built in excel function to estimate Normal Cumulative Probability

Function BSMCall(S, K, T, r, q, Sigma)

BSMCall = S * Exp(-q * T) * Application.NormSDist(d1One(S, K, T, r, q, Sigma)) - K * Exp(-T * r) * _

Application.NormSDist(d2Two(S, K, T, r, q, Sigma))

End Function

'The put is estimated using Put Call Parity

Function BSMPut(S, K, T, r, q, Sigma)

BSMPut = BSMCall(S, K, T, r, q, Sigma) + K * Exp(-r * T) - S * Exp(-q * T)

End Function

'The Delta for the Call

Function DeltaC(S, K, T, r, q, Sigma)

DeltaC = Exp(-q * T) * Application.NormSDist(d1One(S, K, T, r, q, Sigma))

End Function

'The Delta for the Put

Function DeltaP(S, K, T, r, q, Sigma)

DeltaP = -Exp(-q * T) * Application.NormSDist(-d1One(S, K, T, r, q, Sigma))

End Function

'The Gamma for the Call or Put

Function GammaCP(S, K, T, r, q, Sigma)

GammaCP = ((pdf(d1One(S, K, T, r, q, Sigma))) * Exp(-q * T)) / (S * Sigma * Sqr(T))

End Function

'The Vega for the Call or Put

Function Vega(S, K, T, r, q, Sigma)

Vega = S * Sqr(T) * pdf(d1One(S, K, T, r, q, Sigma)) * Exp(-q * T)

End Function

'The Theta for the Call

Function ThetaC(S, K, T, r, q, Sigma)

ThetaC = -S * pdf(d1One(S, K, T, r, q, Sigma)) * Sigma * Exp(-q * T) / (2 * Sqr(T)) + q * S * Application.NormSDist(d1One(S, K, T, r, q, Sigma)) * Exp(-q * T) - r * K * Exp(-r * T) * _

Application.NormSDist(d2Two(S, K, T, r, q, Sigma))

End Function

'The Theta for the Put

Function ThetaP(S, K, T, r, q, Sigma)

ThetaP = -S * pdf(d1One(S, K, T, r, q, Sigma)) * Sigma * Exp(-q * T) / (2 * Sqr(T)) - q * S * Application.NormSDist(-d1One(S, K, T, r, q, Sigma)) * _

Exp(-q * T) + r * K * Exp(-r * T) * Application.NormSDist(-d2Two(S, K, T, r, q, Sigma))

End Function

'The Rho for the Call

Function RhoC(S, K, T, r, q, Sigma)

RhoC = K * T * Exp(-r * T) * Application.NormSDist(d2Two(S, K, T, r, q, Sigma))

End Function

'The Rho for the Put

Function RhoP(S, K, T, r, q, Sigma)

RhoP = -K * T * Exp(-r * T) * Application.NormSDist(-d2Two(S, K, T, r, q, Sigma))

End Function

C++ code for Black Scholes Greeks

//Fabrice Rouah and Brian Byrne

#include <vector>

#include <stdio.h>

#include <math.h>

#include <iostream>

using namespace std;

// N(0,1) density

double f(double x) {

double pi = 4.0*atan(1.0);

return exp(-x*x*0.5) / sqrt(2 * pi);

}

// Boole's Rule

double Boole(double StartPoint, double EndPoint, int n) {

vector<double> X(n + 1, 0.0);

vector<double> Y(n + 1, 0.0);

double delta_x = (EndPoint - StartPoint) / double(n);

for (int i = 0; i <= n; i++) {

X[i] = StartPoint + i*delta_x;

Y[i] = f(X[i]);

}

double sum = 0;

for (int t = 0; t <= (n - 1) / 4; t++) {

int ind = 4 * t;

sum += (1 / 45.0)*(14 * Y[ind] + 64 * Y[ind + 1] + 24 * Y[ind + 2] + 64 * Y[ind + 3] + 14 * Y[ind + 4])*delta_x;

}

return sum;

}

// N(0,1) cdf by Boole's Rule

double N(double x) {

return Boole(-10.0, x, 240);

}

// Black-Scholes Call Price

double BSPrice(double S, double K, double T, double r, double q, double v, char OpType)

{

double d = (log(S / K) + T*(r - q + 0.5*v*v)) / (v*sqrt(T));

double call = S* exp(-q*T)*N(d) - exp(-r*T)*K*N(d - v*sqrt(T));

if (OpType == 'C')

return call;

else

return call - S* exp(-q*T) + K*exp(-r*T);

}

// Black-Scholes Delta

double BSDelta(double S, double K, double T, double r, double q, double v, char OpType)

{

double d = (log(S / K) + T*(r - q + 0.5*v*v)) / (v*sqrt(T));

if (OpType == 'C')

return exp(-q*T)*N(d);

else

return exp(-q*T)*(N(d) - 1);

}

// Black-Scholes Gamma

double BSGamma(double S, double K, double T, double r, double q, double v)

{

double d = (log(S / K) + T*(r - q + 0.5*v*v)) / (v*sqrt(T));

return f(d) * exp(-q*T) / S / v / sqrt(T);

}

// Black-Scholes Vega

double BSVega(double S, double K, double T, double r, double q, double v)

{

double d = (log(S / K) + T*(r - q + 0.5*v*v)) / (v*sqrt(T));

return S*f(d) * exp(-q*T) *sqrt(T);

}

// Black-Scholes Rho

double BSRho(double S, double K, double T, double r, double q, double v, char OpType)

{

double d = (log(S / K) + T*(r - q + 0.5*v*v)) / (v*sqrt(T));

if (OpType == 'C')

return T*K*exp(-r*T)*N(d - v*sqrt(T));

else

return -T*K*exp(-r*T)*N(v*sqrt(T) - d);

}

// Black-Scholes Theta

double BSTheta(double S, double K, double T, double r, double q, double v, char OpType)

{

double d = (log(S / K) + T*(r - q + 0.5*v*v)) / (v*sqrt(T));

if (OpType == 'C')

return (-S*f(d) * exp(-q*T) *v / 2 / sqrt(T)) + q*exp(-q*T)*S*N(d) - r*K*exp(-r*T)*N(d - v*sqrt(T));

else

return (-S*f(d) * exp(-q*T) *v / 2 / sqrt(T)) - q*exp(-q*T)*S*N(-d) + r*K*exp(-r*T)*N(v*sqrt(T) - d);

}

int main()

{

double S = 100.0; // Stock Price

double K = 100.0; // Strike Price

double T = 1; // Years to maturity

double r = 0.05; // Risk free interest rate

double q = 0.0;

double v = 0.20; // Yearly volatility

char OpType = 'P'; // 'C'all or 'P'ut

cout << "Black Scholes Price " << BSPrice(S, K, T, r, q, v, OpType) << endl;

cout << "Black Scholes Delta " << BSDelta(S, K, T, r, q, v, OpType) << endl;

cout << "Black Scholes Gamma " << BSGamma(S, K, T, r, q, v) << endl;

cout << "Black Scholes Vega " << BSVega(S, K, T, r, q, v) << endl;

cout << "Black Scholes Rho " << BSRho(S, K, T, r, q, v, OpType) << endl;

cout << "Black Scholes Theta " << BSTheta(S, K, T, r, q, v, OpType)

<< " per year, or " << BSTheta(S, K, T, r, q, v, OpType) / 365.0 << " per day." << endl;

system("PAUSE");

};

Black Scholes Valuation, Greeks and Implied Volatility from Diethelm Weurtz

The statistical software package R is one of the most promising tools for rapid prototyping of financial applications. Here we will introduce programs or functions written in R code to estimate the Black Scholes model, their Greeks (or parameter sensitivities) and Implied Volatility. The code is of a very high standard and follows Espen Haug. Haug E.G., The Complete Guide to Option Pricing Formulas. The R code is written to achieve high efficiency and optimization.

R code for Black Scholes Greeks with Graphing

Robert McDonald, Professor of Finance at the Kellogg School of Management, Northwestern University and author of "Derivatives Markets" (Pearson, 2013) has put together the derivmkts package for R. Here I used the template for estimating Greeks as presented in Robert's derivmkts package and in his book's namesake.

library("derivmkts")

s=100; k=100; v=0.2; r=0.05; tt=1; d=0;

greeks(bscall(s, k, v, r, tt, d), complete=FALSE, long=FALSE, initcaps=TRUE)

greeks2(bscall, list(s=s, k=k, v=v, r=r, tt=tt, d=d))

greeks2(bscall, list(s=s, k=k, v=v, r=r, tt=tt, d=d))[c('Delta', 'Gamma'), ]

bsopt(s, k, v, r, tt, d)

bsopt(s, c(90, 100, 110), v, r, tt, d)

bsopt(s, c(90, 100, 110), v, r, tt, d)[['Call']][c('Delta', 'Gamma'), ]

## plot Greeks for calls and puts for 500 different stock prices

k <- 100; v <- 0.20; r <- 0.05; tt <- 1; d <- 0

S <- seq(10, 200, by=5)

Call <- greeks(bscall(S, k, v, r, tt, d))

Put <- greeks(bsput(S, k, v, r, tt, d))

y <- list(Call=Call, Put=Put)

par(mfrow=c(4, 4), mar=c(2, 2, 2, 2)) ## create a 4x4 plot

for (i in names(y)) {

for (j in rownames(y[[i]])) { ## loop over greeks

plot(S, y[[i]][j, ], main=paste(i, j), ylab=j, type='l')

}

}

Derivmkts R Package using Google Colab

Options greeks are mathematical derivatives of the option price with respect to inputs; see McDonald (2013), Chapters 12 and 13, for a discussion of the greeks for vanilla options. Greeks for vanilla and barrier options can be computed using the greeks function, which is a wrapper for any pricing function that returns the option price and which uses the default naming of inputs. In the Derivmkts package ,Prefessor McDonald has two alternative functions that return Greeks:

The bsopt function by default produces prices and Greeks for European calls and puts.

The greeks2 function takes as arguments the name of the pricing function and then inputs as a list

Python Code for Black Scholes Greeks

The Python code below was obtained from Davis Edwards github. We stripped back the full code to focus on estimating the Black Scholes Greeks and we follow the Generalized Black Scholes approach referenced in Espen Haug.

# This document demonstrates a Python implementation of some option models described in books written by Davis

# Edwards: "Energy Trading and Investing", "Risk Management in Trading", "Energy Investing Demystified".

# The gbs approach follows Espen Haug

# for backward compatability with Python 2.7

from __future__ import division

# import necessary libaries

import math

import numpy as np

from scipy.stats import norm

from scipy.stats import mvn

# The primary class for calculating Generalized Black Scholes option prices and deltas

# Inputs: option_type = "p" or "c", fs = price of underlying, x = strike, t = time to expiration, r = risk free rate

# b = cost of carry, v = implied volatility

# Outputs: value, delta, gamma, theta, vega, rho

def _gbs(option_type, fs, x, t, r, b, v):

# Create preliminary calculations

t__sqrt = math.sqrt(t)

d1 = (math.log(fs / x) + (b + (v * v) / 2) * t) / (v * t__sqrt)

d2 = d1 - v * t__sqrt

if option_type == "c":

# it's a call

# _debug(" Call Option")

value = fs * math.exp((b - r) * t) * norm.cdf(d1) - x * math.exp(-r * t) * norm.cdf(d2)

delta = math.exp((b - r) * t) * norm.cdf(d1)

gamma = math.exp((b - r) * t) * norm.pdf(d1) / (fs * v * t__sqrt)

theta = -(fs * v * math.exp((b - r) * t) * norm.pdf(d1)) / (2 * t__sqrt) - (b - r) * fs * math.exp(

(b - r) * t) * norm.cdf(d1) - r * x * math.exp(-r * t) * norm.cdf(d2)

vega = math.exp((b - r) * t) * fs * t__sqrt * norm.pdf(d1)

rho = x * t * math.exp(-r * t) * norm.cdf(d2)

else:

# it's a put

# _debug(" Put Option")

value = x * math.exp(-r * t) * norm.cdf(-d2) - (fs * math.exp((b - r) * t) * norm.cdf(-d1))

delta = -math.exp((b - r) * t) * norm.cdf(-d1)

gamma = math.exp((b - r) * t) * norm.pdf(d1) / (fs * v * t__sqrt)

theta = -(fs * v * math.exp((b - r) * t) * norm.pdf(d1)) / (2 * t__sqrt) + (b - r) * fs * math.exp(

(b - r) * t) * norm.cdf(-d1) + r * x * math.exp(-r * t) * norm.cdf(-d2)

vega = math.exp((b - r) * t) * fs * t__sqrt * norm.pdf(d1)

rho = -x * t * math.exp(-r * t) * norm.cdf(-d2)

return value, delta, gamma, theta, vega, rho

# This is the public interface for European Options

# Each call does a little bit of processing and then calls the calculations located in the _gbs module

# Inputs:

# option_type = "p" or "c"

# fs = price of underlying

# x = strike

# t = time to expiration

# v = implied volatility

# r = risk free rate

# q = dividend payment

# b = cost of carry

# Outputs:

# value = price of the option

# delta = first derivative of value with respect to price of underlying

# gamma = second derivative of value w.r.t price of underlying

# theta = first derivative of value w.r.t. time to expiration

# vega = first derivative of value w.r.t. implied volatility

# rho = first derivative of value w.r.t. risk free rates

# Black Scholes: stock Options (no dividend yield)

def black_scholes(option_type, fs, x, t, r, v):

b = r

return _gbs(option_type, fs, x, t, r, b, v)

# Merton Model: Stocks Index, stocks with a continuous dividend yields

def merton(option_type, fs, x, t, r, q, v):

b = r - q

return _gbs(option_type, fs, x, t, r, b, v)

def black_76(option_type, fs, x, t, r, v):

b = 0

return _gbs(option_type, fs, x, t, r, b, v)

# FX Options

def garman_kohlhagen(option_type, fs, x, t, r, rf, v):

b = r - rf

return _gbs(option_type, fs, x, t, r, b, v)

#black_scholes('c', fs=100, x=100, t=1, r=0.05, v=0.20)

# merton(option_type, fs, x, t, r, q, v):

# b = r - q

merton('c', fs=100, x=100, t=1, r=0.05, q = 0, v=0.20)

# c = black_scholes('c', fs=100, x=100, t=1, r=0.05, v=0.20)

# print( c )