Coding towards CFA (17) – Implied Volatility Calculation

Coding towards CFA (17) – Implied Volatility Calculation

Implied volatility (IV) is a key metric in options trading and risk management. It is derived from the market price of an option and reflects the market’s consensus view on the expected future volatility of the underlying asset.

There is a lot of information that can be interpreted from IV. For example, high IV often signals uncertainty, fear, or the expectation of significant market-moving events, which drives up option prices. Traders can profit by identifying discrepancies between implied volatility and the expected future volatility of the underlying asset.

IV can be calculated using an option pricing model, such as the Black-Scholes Model (BSM), by finding the implied volatility that makes the theoretical price of an option equal to its market price. This process involves solving an equation where the market price of the option is set equal to the option price calculated from the model, with implied volatility as the variable (sigma).

Market Option Price = bsm(S, K, T, r, Sigma)

In this way, the calculation of IV is essentially a root-finding problem to solve the equation. The Newton-Raphson method or the scipy brentq function can be used for root-finding. In this blog post, I will implement the Newton-Raphson method in Python and DolphinDB as examples.

First, we create the function for calculate option price and vega using BSM model.

To calculate the implied volatility, we first provide an initial guess for the IV and use it to calculate the option price and vega with the function we created earlier. We then compare the market option price with the calculated option price based on the guessed IV.

If the difference is smaller than the given threshold, the guessed IV is the right IV we are after. If not, we need to update the sigma guess with the Newton-Raphson formula, and repeat the process.

sigma = sigma – (calculated option price – market option price) / option vega

The full Python code and DolphinDB code can be found at the bottom of this blog post.

With the function for calculating the IV, we can further generate the volatility surface, which is a 3D visualisation of the market’s expectations for volatility at various strike prices over time.

Full Code – Python

import yfinance as yf
import numpy as np
from scipy.stats import norm
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D

def calculate_option_vega(s, k, T, r, sigma, option_type="call"):

    # Calculate d1 and d2
    d1 = (np.log(s / k) + (r + 0.5 * sigma ** 2) * T) / (sigma * np.sqrt(T))
    d2 = d1 - sigma * np.sqrt(T)

    # Calculate vega
    option_vega = s * norm.pdf(d1) * np.sqrt(T)

    # Calculate option price
    if option_type == "call":
        option_price = s * norm.cdf(d1) - k * np.exp(-r * T) * norm.cdf(d2)
    else:
        option_price = k * np.exp(-r * T) * norm.cdf(-d2) - s * norm.cdf(-d1)

    return option_vega, option_price

# Resolve implied volatility with Newton-Raphson method 
def calculate_implied_volatility(s, k, T, r, marekt_option_price, 
        initial_guess=0.2, option_type="call", tolerance=1e-6, iterations=100):
    
    # Set initial guess 
    sigma = initial_guess
    
    for i in range(iterations):
        print("sigma = "+str(sigma))

        # Calculate the price and vega with current sigma guess
        option_vega, option_price = calculate_option_vega(s, k, T, r, sigma)
        
        # Difference between the market option price and the 
        # option price calculated with the sigma guess
        price_diff = option_price - marekt_option_price
        
        # return implied volatility when the difference 
        # is smaller than the given threshold
        if abs(price_diff) < tolerance:
            return sigma
        
        # Calculate new sigma using the Newton-Raphson formula
        if option_vega <= 0:
            return None
        else:
            sigma = sigma - price_diff / option_vega
    
    return None

symbol = 'TSLA'
stock = yf.Ticker(symbol)
expiration_dates = stock.options[:6]
option_chains = {e:stock.option_chain(e) for e in expiration_dates}
TTMs = {e:-(np.datetime64('today', 'D') - np.datetime64(e, 'D')).astype(int) /365
            for e in expiration_dates}
s = stock.history(period="1d")['Close'].iloc[-1]
r = 0.042  

strikes = [435, 440, 445, 450, 455.0, 460, 465, 470, 475]
implied_vols = []

for e, option_chain in option_chains.items():

    calls = option_chain.calls
    prices = [calls[calls['strike'] == s]['lastPrice'].iloc[0] for s in strikes]
    ivs = []
    for i in range(len(strikes)):
        iv = calculate_implied_volatility(s, strikes[i], TTMs[e], r, prices[i])
        ivs.append(iv)
    implied_vols.append(ivs)

fig = plt.figure(figsize=(8, 5))
ax = fig.add_subplot(111, projection='3d')
X, Y = np.meshgrid( np.array(strikes), np.array([v for v in TTMs.values()]))

ax.plot_surface(X, Y, np.array(implied_vols), cmap='viridis', edgecolor='none')
ax.set_xlabel('Strike Price')
ax.set_ylabel('Time to Expiration')
ax.set_zlabel('Implied Volatility')
ax.set_title('Option Volatility Surface')
plt.show()

Full Code – DolphinDB

def calculate_option_vega(s, k, T, r, sigma, option_type="call"){

    // Calculate d1 and d2
    d1 = (log(s/k) + (r + 0.5 * pow(sigma,2)) * T) / (sigma * sqrt(T))
    d2 = d1 - sigma * sqrt(T)

    // Calculate vega
    option_vega = s * pdfNormal(0, 1, d1) * sqrt(T)

    // Calculate option price
    if (option_type == "call") {
        option_price = s * cdfNormal(0, 1, d1) - k * exp(-r * T) * cdfNormal(0, 1, d2)
    }
    else {
        option_price = k * exp(-r * T) * cdfNormal(0, 1, -d2) - s * cdfNormal(0, 1, -d1)
    }

    return option_vega, option_price
}

// Resolve implied volatility with Newton-Raphson method 
def calculate_implied_volatility(s, k, T, r, marekt_option_price, option_type="call", tolerance=1e-6, iterations=100){
    
    // Set initial guess 
    sigma = 0.2 
    
    for (i in 0:iterations){

        // Calculate the price and vega with current sigma guess
        option_vega, option_price = calculate_option_vega(s, k, T, r, sigma)
        
        // Difference between the market option price 
        // and the option price calculated with the sigma guess
        price_diff = option_price - marekt_option_price
        
        // return implied volatility when the difference is smaller than the given threshold
        if (abs(price_diff) < tolerance){
            return sigma
        }
        
        // Calculate new sigma using the Newton-Raphson formula
        if (option_vega <= 0) {
            return NULL
        } else {
            sigma = sigma - price_diff / option_vega
        }
    }
    
    return NULL
}

s, k, T, r, marekt_option_price = 100, 100, 0.5, 0.05, 4

implied_volatility = calculate_implied_volatility(s, k, r, T, marekt_option_price)
print("Implied Volatility: "+implied_volatility)

Leave a comment