Newton_Raphson_Method,_Greeks.ipynb
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from scipy import stats
#call option의 BS formula
def bs_call(S,K,T,r,vol):
d1 = (np.log(S/K) + (r + 0.5 * vol ** 2) * T) / (vol * np.sqrt(T))
d2 = d1 - vol * np.sqrt(T)
return S * stats.norm.cdf(d1) - K * np.exp(-r * T) * stats.norm.cdf(d2)
def bs_put(S,K,T,r,vol):
d1 = (np.log(S/K) + (r + 0.5 * vol ** 2) * T) / (vol * np.sqrt(T))
d2 = d1 - vol * np.sqrt(T)
return S * stats.norm.cdf(-d1) + K * np.exp(-r * T) * stats.norm.cdf(d2)
def bs_vega(S,K,T,r,vol):
d1 = (np.log(S/K) + (r + 0.5 * vol ** 2) * T) / (vol * np.sqrt(T))
return S * stats.norm.pdf(d1) * np.sqrt(T) #cdf미분이 pdf이다
def find_vol(market_value, S, K, T, r, init_vol, option_type):
sigma = init_vol
max_iteration = 200 #newton방법을 최대 몇 번 반복할 것인지
precision = 0.00001 #오차가 0.00001보다 작으면 마치겠다
if option_type == 'call':
bs_price = bs_call(S, K, T, r, sigma)
vega = bs_vega(S, K, T, r, sigma)
diff = market_value - bs_price
if (abs(diff) < precision):
return sigma #우리가 찾는 implied volatility
sigma = sigma + diff / vega #Newton-Raphson solution = x0 - f(x) / f'(x)
return sigma
elif option_type == 'put':
bs_price = bs_put(S, K, T, r, sigma)
vega = bs_vega(S, K, T, r, sigma)
diff = market_value - bs_price
if (abs(diff) < precision):
return sigma #우리가 찾는 implied volatility
sigma = sigma + diff / vega #Newton-Raphson solution = x0 - f(x) / f'(x)
return sigma
S = 100
K = 100
T = 1
r = 0.01
vol = 0.25
init_vol = 0.35
market_value = bs_call(S, K, T, r, vol)
implied_vol = find_vol(market_value, S, K, T, r, init_vol, 'call') #이때 시그마는 임의로 정한 시그마
print('implied vol: %.2f%%' % (implied_vol * 100))
print('market price: %.2f' % market_value)
print('model price: %.2f' % bs_call(S, K, T, r, implied_vol))
implied vol: 24.96% market price: 10.40 model price: 10.39
import plotly.graph_objs as go
from plotly.offline import download_plotlyjs, init_notebook_mode, plot, iplot
def delta(option_type,S,K,r,T,sigma):
d1 = (np.log(S/K) + (r + 0.5 * (sigma ** 2)) * T) / (sigma * np.sqrt(T))
if option_type == 'call':
return stats.norm.cdf(d1)
elif option_type == 'put':
return stats.norm.cdf(d1) - 1
K,r,sigma = 100, 0.01, 0.25
n1,n2 = 100, 50
start_t,end_t,start_s,end_s = 0.000001, 1, 0.000001, 200 #격자
T = np.linspace(start_t,end_t,n1)
S = np.linspace(start_s,end_s,n1)
T,S = np.meshgrid(T,S)
call_delta = delta('call',S,K,r,T,sigma)
put_delta = delta('put',S,K,r,T,sigma)
trace = go.Surface(x=T,y=S,z=call_delta,colorscale='Jet',opacity=0.8,
contours_x=dict(show=True,color='black',start=start_t,end=end_t,size=(end_t-start_t)/n2,project_x=True),
contours_y=dict(show=True,color='black',start=start_s,end=end_s,size=(end_s-start_s)/n2,project_y=True)
)
data = [trace]
layout = go.Layout(title='Call Delta Surface',
scene={'xaxis':{'title':'maturity'},'yaxis':{'title':'spot price'},'zaxis':{'title':'delta'}},
width=600,height=600,autosize=True)
fig = go.Figure(data=data,layout=layout)
iplot(fig)
#직접 돌려서 그래프 확인하길
trace = go.Surface(x=T,y=S,z=put_delta,colorscale='Jet',opacity=0.8,
contours_x=dict(show=True,color='black',start=start_t,end=end_t,size=(end_t-start_t)/n2,project_x=True),
contours_y=dict(show=True,color='black',start=start_s,end=end_s,size=(end_s-start_s)/n2,project_y=True)
)
data = [trace]
layout = go.Layout(title='Put Delta Surface',
scene={'xaxis':{'title':'maturity'},'yaxis':{'title':'spot price'},'zaxis':{'title':'delta'}},
width=600,height=600,autosize=True)
fig = go.Figure(data=data,layout=layout)
iplot(fig)
def theta(option_type,S,K,r,T,sigma):
d1 = (np.log(S/K) + (r + 0.5 * (sigma ** 2)) * T) / (sigma * np.sqrt(T))
d2 = d1 - sigma*np.sqrt(T)
if option_type == 'call':
theta = (-sigma * S * stats.norm.pdf(d1)) / (2 * np.sqrt(T)) - r * K * np.exp(-r*T) * stats.norm.cdf(d2)
elif option_type == 'put':
theta = (-sigma * S * stats.norm.pdf(d1)) / (2 * np.sqrt(T)) + r * K * np.exp(-r*T) * stats.norm.cdf(-d2)
return theta
call_theta = theta('call',S,K,r,T,sigma)
put_theta = theta('put',S,K,r,T,sigma)
trace = go.Surface(x=T,y=S,z=call_theta,colorscale='Jet',opacity=0.8,
contours_x=dict(show=True,color='black',start=start_t,end=end_t,size=(end_t-start_t)/n2,project_x=True),
contours_y=dict(show=True,color='black',start=start_s,end=end_s,size=(end_s-start_s)/n2,project_y=True)
)
data = [trace]
layout = go.Layout(title='Call Theta Surface',
scene={'xaxis':{'title':'maturity'},'yaxis':{'title':'spot price'},'zaxis':{'title':'theta'}},
width=600,height=600,autosize=True)
fig = go.Figure(data=data,layout=layout)
iplot(fig)
trace = go.Surface(x=T,y=S,z=put_theta,colorscale='Jet',opacity=0.8,
contours_x=dict(show=True,color='black',start=start_t,end=end_t,size=(end_t-start_t)/n2,project_x=True),
contours_y=dict(show=True,color='black',start=start_s,end=end_s,size=(end_s-start_s)/n2,project_y=True)
)
data = [trace]
layout = go.Layout(title='Put Theta Surface',
scene={'xaxis':{'title':'maturity'},'yaxis':{'title':'spot price'},'zaxis':{'title':'theta'}},
width=600,height=600,autosize=True)
fig = go.Figure(data=data,layout=layout)
iplot(fig)
def gamma(S,K,r,T,sigma):
d1 = (np.log(S/K) + (r + 0.5 * (sigma ** 2)) * T) / (sigma * np.sqrt(T))
gamma = stats.norm.pdf(d1) / (sigma * S * np.sqrt(T))
return gamma
gamma = gamma(S,K,r,T,sigma)
trace = go.Surface(x=T,y=S,z=gamma,colorscale='Jet',opacity=0.8,
contours_x=dict(show=True,color='black',start=start_t,end=end_t,size=(end_t-start_t)/n2,project_x=True),
contours_y=dict(show=True,color='black',start=start_s,end=end_s,size=(end_s-start_s)/n2,project_y=True)
)
data = [trace]
layout = go.Layout(title='Gamma Surface',
scene={'xaxis':{'title':'maturity'},'yaxis':{'title':'spot price'},'zaxis':{'title':'gamma'}},
width=600,height=600,autosize=True)
fig = go.Figure(data=data,layout=layout)
iplot(fig)