So my goal is to calculate option greeks with as few manual inputs as possible. I managed to get the IV for at the money options but then when I try further OTM strikes my results get completely screwed.
Im running this loop to estimate IV:
def find_vol(target_value, call_put, S, K, T, r):
MAX_ITERATIONS = 200
PRECISION = 0.01
sigma = 0.01
s = {}
for i in range(0, MAX_ITERATIONS):
price = bs_price(call_put, S, K, T, r, sigma)
vega = bs_vega(call_put, S, K, T, r, sigma)
price = price
diff = target_value - price # our root
print({'sigma':sigma,'vega':vega,'price':price, 'target_value':target_value,'diff':diff})
if abs(diff) < PRECISION:
return sigma
sigma = sigma + diff/vega
s[diff]=sigma
mindiff =min(list(s.keys()),key=abs)
sigma = s[mindiff]
# value wasn't found, return best guess so far
return sigma
My issue is that with certain options the loop will get stuck because Vega goes to 0 and then this line doesn't make sense anymore: sigma = sigma + diff/Vega
when the problem happens it looks like this :
The issue is probably due to the fact that I am using BS model on American options which is clearly wrong but I am unsure if my results will improve when using Binomial. Suggestions would be much appreciated!
Full working Code: '''python
import datetime
from scipy.stats import norm
from math import exp, log,sqrt
import quantsbin.derivativepricing as qbdp
def find_vol(target_value, call_put, S, K, T, r):
MAX_ITERATIONS = 200
PRECISION = 0.01
sigma = 0.5
s = {}
for i in range(0, MAX_ITERATIONS):
price = bs_price(call_put, S, K, T, r, sigma)
vega = bs_vega(call_put, S, K, T, r, sigma)
price = price
diff = target_value - price # our root
#print({'sigma':sigma,'vega':vega,'price':price, 'target_value':target_value,'diff':diff})
if abs(diff) < PRECISION:
return sigma
sigma = sigma + diff/vega # f(x) / f'(x)
s[diff]=sigma
mindiff =min(list(s.keys()),key=abs)
#print(mindiff)
sigma = s[mindiff]
#print(sigma)
# value wasn't found, return best guess so far
return sigma
n = norm.pdf
N = norm.cdf
def bs_price(cp_flag,S,K,T,r,v,q=0.0):
d1 = (log(S/K)+(r+v*v/2.)*T)/(v*sqrt(T))
d2 = d1-v*sqrt(T)
if cp_flag == 'c':
price = S*exp(-q*T)*N(d1)-K*exp(-r*T)*N(d2)
else:
price = K*exp(-r*T)*N(-d2)-S*exp(-q*T)*N(-d1)
return price
def bs_vega(cp_flag,S,K,T,r,v,q=0.0):
d1 = (log(S/K)+(r+v*v/2.)*T)/(v*sqrt(T))
return S * sqrt(T)*n(d1)
def calculate_iv_greeks(price_market,strike,expiration,date_calc,spot_underlying,rfr,kind,div_yield=0):
exp_y,exp_m,exp_d = expiration
calc_y, calc_m, calc_d = date_calc
expiration = datetime.date(exp_y, exp_m, exp_d)
date_calc = datetime.date(calc_y, calc_m, calc_d)
T = (expiration - date_calc).days / 365.
#T = (expiration - date_calc).days
r = ((rfr/100)/365)*(expiration - date_calc).days
iv = find_vol(price_market,kind,spot_underlying,strike,T,r)
if kind =='C': type = 'Call'
else: type = 'Put'
expiration = datetime.date(exp_y,exp_m,exp_d).strftime("%Y%m%d")
date_calc = datetime.date(calc_y, calc_m, calc_d).strftime("%Y%m%d")
equity_option1 = qbdp.EqOption(option_type=type, strike=strike, expiry_date=str(expiration))
if div_yield != 0:
eq1_engine = equity_option1.engine(model='Binomial', pricing_date=str(date_calc), spot0=spot_underlying, rf_rate=rfr/100,
volatility=iv, yield_div=div_yield)
print('Engine:Binomial')
else:
eq1_engine = equity_option1.engine(model='BSM', pricing_date=str(date_calc), spot0=spot_underlying, rf_rate=rfr/100,
volatility=iv)
print('Engine:BSM')
D,G,TH,V,R,P = eq1_engine.risk_parameters().values()
print(eq1_engine.risk_parameters())
return eq1_engine.valuation(),iv, D,TH
V_market = 5.85
K = 285
S = 289.85
cp = 'c'
greeks=calculate_iv_greeks(price_market=V_market,strike=K,expiration=. [2019,5,1],date_calc=[2019,4,17],
spot_underlying=S,rfr=2.41,kind=cp,div_yield=0)
print(greeks)
print('Broker values: ', {'IV':0.1217 ,'Delta':0.765 })
'''
