• C++ Programming for Financial Engineering
    Highly recommended by thousands of MFE students. Covers essential C++ topics with applications to financial engineering. Learn more Join!
    Python for Finance with Intro to Data Science
    Gain practical understanding of Python to read, understand, and write professional Python code for your first day on the job. Learn more Join!
    An Intuition-Based Options Primer for FE
    Ideal for entry level positions interviews and graduate studies, specializing in options trading arbitrage and options valuation models. Learn more Join!

Volatility Surface Modeling for Equity markets


i am trying to construct volatility surface projection for stocks via the Vanna-Volga Implied Volatility method.

The Vol Smile i am getting is way wrong if i go by 3 inputs for the given vol curve:

Input TYPESigmaStrike
25 Delta PUT0.246814650
25 Delta Call0.208814950

Screenshot 2021-03-23 at 4.29.40 PM.png

Actual Vol Smile: 1.95 Days to expiry, Atm Price is 14850
Screenshot 2021-03-23 at 6.21.13 PM.png

Generated Vol Smile via Vanna-Volga Implied Volatility

I believe since the Volsmile is shifted to the right, this is causing problems, what adjustments should I make to get it done right

Here is the code i am using:

import numpy as np
import scipy.stats as sct
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
from matplotlib import cm

def volsurface(T,DP,DC,ATM,st,cmp,delta):
    # T is is time in days
    # DP,DC,ATM are Sigma 25 Delta Put,25 Delta Call, And Sigma of ATM
    # st is strike array
    # CMP is ATM
     #delta is =0.25
    T = np.array(T)
    sigma25DeltaPut = np.array(DP)
    sigmaATM = np.array(ATM)
    sigma25DeltaCall = np.array(DC)

    #discount factor Applicable in currencies
    r = np.array([0.00])
    q = np.array([0.00])
    mu = r - q
    Dd = np.exp(-np.multiply(r,T))
    Df = np.exp(-np.multiply(q,T))

    S0 = cmp

    F = S0*np.exp(np.multiply(mu,T))
    alpha = -sct.norm.ppf(delta*np.reciprocal(Df))
    K25DPut = F*np.exp(-(alpha*sigma25DeltaPut*np.sqrt(T)) +0.5*(sigma25DeltaPut**2)*T)
    KATM = F*np.exp(0.5*(sigmaATM**2)*T)
    K25DCall = F*np.exp((alpha*sigma25DeltaCall*np.sqrt(T))+0.5*(sigma25DeltaCall**2)*T)
    def d1(F,K,sigma,t):
        dplus = (np.log(F/K)+0.5*(sigma**2)*t)/(sigma*np.sqrt(t))
        return dplus

    def d2(F,K,sigma,t):
        dminus = d1(F,K,sigma,t) - sigma*np.sqrt(t)
        return dminus

    def VannavolgaImpliedVol(F,K,t,K1,K2,K3,sigma1,sigma2,sigma3):
        #sigma1,2,3 are 25d put,atm 25d call
        y1 = np.array(np.log(K2/K)*np.log(K3/K))/np.array(np.log(K2/K1)*np.log(K3/K1))
        y2 = np.array(np.log(K/K1)*np.log(K3/K))/np.array(np.log(K2/K1)*np.log(K3/K2))
        y3 = np.array(np.log(K/K1)*np.log(K/K2))/np.array(np.log(K3/K1)*np.log(K3/K2))
        P = (y1*sigma1 + y2*sigma2 + y3*sigma3) - sigma2
        Q = y1 * d1(F,K1,sigma1,t) * d2(F,K1,sigma1,t) * ((sigma1-sigma2)**2) \
            + y2 * d1(F,K2,sigma2,t) * d2(F,K2,sigma2,t) * ((sigma2-sigma2)**2) \
            + y3 * d1(F,K3,sigma3,t) * d2(F,K3,sigma3,t) * ((sigma3-sigma2)**2)
        d1d2 = d1(F,K,sigma2,t) * d2(F,K,sigma2,t)
        sigma = sigma2 + (-sigma2 + np.sqrt(sigma2**2 + d1d2 * (2*sigma2*P+Q)))/(d1d2)
        return sigma
    VVImpliedVol = VannavolgaImpliedVol(F, st.reshape(-1,1), T, K25DPut, KATM, K25DCall, sigma25DeltaPut, sigmaATM,sigma25DeltaCall)
    return VVImpliedVol
Last edited: