请优化下面的代码使其能够通过输入一组行权价来绘制波动率微笑曲线 import numpy as np from scipy.stats import norm from scipy.optimize import minimize import matplotlib.pyplot as plt def bs_option_price(S, K, r, q, sigma, T, option_type): d1 = (np.log(S/K) + (r - q + sigma**2/2) * T) / (sigma * np.sqrt(T)) d2 = d1 - sigma * np.sqrt(T) if option_type == 'call': Nd1 = norm.cdf(d1) Nd2 = norm.cdf(d2) option_price = S * np.exp(-q * T) * Nd1 - K * np.exp(-r * T) * Nd2 elif option_type == 'put': Nd1 = norm.cdf(-d1) Nd2 = norm.cdf(-d2) option_price = K * np.exp(-r * T) * (1 - Nd2) - S * np.exp(-q * T) * (1 - Nd1) else: raise ValueError('Invalid option type') return option_price def implied_volatility(S, K, r, q, T, option_price, option_type): obj_fun = lambda sigma: (bs_option_price(S, K, r, q, sigma, T, option_type) - option_price)**2 res = minimize(obj_fun, x0=0.2) return res.x[0] def smile_curve(S, r, q, T, option_type, strike_range, option_prices): vols = [] for K, option_price in zip(strike_range, option_prices): vol = implied_volatility(S, K, r, q, T, option_price, option_type) vols.append(vol) plt.plot(strike_range, vols) plt.xlabel('Strike') plt.ylabel('Implied Volatility') plt.title(f'{option_type.capitalize()} Implied Volatility Smile') plt.show() S = 100 r = 0.05 q = 0.02 T = 0.25 option_type = 'call' strike_range = np.linspace(80, 120, 41) option_prices = [13.05, 10.40, 7.93, 5.75, 4.00, 2.66, 1.68, 1.02, 0.58, 0.31, 0.15, 0.07, 0.03, 0.01, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.01, 0.03, 0.07, 0.14, 0.25, 0.42, 0.67, 1.00, 1.44, 2.02, 2.74, 3.60, 4.60, 5.73, 7.00, 8.39, 9.92, 11.57, 13.34, 15.24] smile_curve(S, r, q, T, option_type, strike_range, option_prices)
时间: 2023-11-23 11:06:51 浏览: 47
可以通过向量化计算和使用更高效的求解器来优化代码。下面是优化后的代码:
```
import numpy as np
from scipy.stats import norm
from scipy.optimize import root_scalar
import matplotlib.pyplot as plt
def bs_option_price(S, K, r, q, sigma, T, option_type):
d1 = (np.log(S/K) + (r - q + sigma**2/2) * T) / (sigma * np.sqrt(T))
d2 = d1 - sigma * np.sqrt(T)
if option_type == 'call':
Nd1 = norm.cdf(d1)
Nd2 = norm.cdf(d2)
option_price = S * np.exp(-q * T) * Nd1 - K * np.exp(-r * T) * Nd2
elif option_type == 'put':
Nd1 = norm.cdf(-d1)
Nd2 = norm.cdf(-d2)
option_price = K * np.exp(-r * T) * (1 - Nd2) - S * np.exp(-q * T) * (1 - Nd1)
else:
raise ValueError('Invalid option type')
return option_price
def implied_volatility(S, K, r, q, T, option_price, option_type):
obj_fun = lambda sigma: (bs_option_price(S, K, r, q, sigma, T, option_type) - option_price)**2
res = root_scalar(obj_fun, bracket=[0.01, 2], method='brentq')
return res.root
def smile_curve(S, r, q, T, option_type, strike_range, option_prices):
vols = np.vectorize(implied_volatility)(S, strike_range, r, q, T, option_prices, option_type)
plt.plot(strike_range, vols)
plt.xlabel('Strike')
plt.ylabel('Implied Volatility')
plt.title(f'{option_type.capitalize()} Implied Volatility Smile')
plt.show()
S = 100
r = 0.05
q = 0.02
T = 0.25
option_type = 'call'
strike_range = np.linspace(80, 120, 41)
option_prices = [13.05, 10.40, 7.93, 5.75, 4.00, 2.66, 1.68, 1.02, 0.58, 0.31, 0.15, 0.07, 0.03, 0.01, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.01, 0.03, 0.07, 0.14, 0.25, 0.42, 0.67, 1.00, 1.44, 2.02, 2.74, 3.60, 4.60, 5.73, 7.00, 8.39, 9.92, 11.57, 13.34, 15.24]
smile_curve(S, r, q, T, option_type, strike_range, option_prices)
```
这个代码中,我们使用了 `np.vectorize` 函数对 `implied_volatility` 函数进行向量化计算,从而避免了使用循环。同时,我们使用了 `root_scalar` 函数来代替 `minimize` 函数,因为 `root_scalar` 函数通常比 `minimize` 函数更高效。我们还增加了一个 `bracket` 参数来指定求解器的搜索范围,从而加快了求解的速度。