用pythonPrice this barrier option by the Heston stochastic volatility model using Monte Carlo simulation by the time discretization, also plot one sample path which hits the barrier and one path that never hits the barrier.
时间: 2024-02-25 14:56:14 浏览: 20
下面是使用 Python 对 Heston 随机波动率模型进行蒙特卡罗模拟定价障碍期权的示例代码,同时画出一条触碰障碍的样本路径和一条未触碰障碍的样本路径:
```python
import numpy as np
import matplotlib.pyplot as plt
# Heston 随机波动率模型参数
V0 = 0.05
theta = 0.05
kappa = 1.5
sigma = 0.3
rho = -0.3
# 期权参数
S0 = 100
K = 110
r = 0.03
T = 1
H = 120
option_type = 'call'
# 蒙特卡罗模拟参数
N = 10000
dt = 1/250
# 计算障碍期权价格
def barrier_option_price():
# 初始化模拟路径
S = np.zeros((N, T * 250 + 1))
V = np.zeros((N, T * 250 + 1))
S[:, 0] = S0
V[:, 0] = V0
# 生成随机数
Zv = np.random.normal(size=(N, T * 250))
Zs = rho * Zv + np.sqrt(1 - rho**2) * np.random.normal(size=(N, T * 250))
# 模拟路径
for i in range(1, T * 250 + 1):
V[:, i] = np.maximum(V[:, i-1] + kappa * (theta - V[:, i-1]) * dt + sigma * np.sqrt(V[:, i-1] * dt) * Zv[:, i-1], 0)
S[:, i] = S[:, i-1] * np.exp((r - V[:, i-1]/2) * dt + np.sqrt(V[:, i-1] * dt) * Zs[:, i-1])
# 计算到期时是否触碰障碍
hit_barrier = np.any(S[:, 1:] > H, axis=1)
# 计算期望收益
if option_type == 'call':
payoff = np.maximum(S[:, -1] - K, 0)
elif option_type == 'put':
payoff = np.maximum(K - S[:, -1], 0)
else:
raise ValueError('Invalid option type')
# 计算障碍期权价格
price = np.mean(payoff * np.logical_not(hit_barrier))
# 返回结果
return price
# 计算障碍期权价格
price = barrier_option_price()
print('障碍期权价格:', price)
# 绘制一条触碰障碍的样本路径
S = np.zeros(T * 250 + 1)
V = np.zeros(T * 250 + 1)
S[0] = S0
V[0] = V0
for i in range(1, T * 250 + 1):
V[i] = np.maximum(V[i-1] + kappa * (theta - V[i-1]) * dt + sigma * np.sqrt(V[i-1] * dt) * np.random.normal(), 0)
S[i] = S[i-1] * np.exp((r - V[i-1]/2) * dt + np.sqrt(V[i-1] * dt) * np.random.normal())
if S[i] > H:
break
plt.plot(np.arange(T * 250 + 1), S)
plt.title('Sample Path Hits Barrier')
plt.xlabel('Time')
plt.ylabel('Stock Price')
plt.show()
# 绘制一条未触碰障碍的样本路径
S = np.zeros(T * 250 + 1)
V = np.zeros(T * 250 + 1)
S[0] = S0
V[0] = V0
for i in range(1, T * 250 + 1):
V[i] = np.maximum(V[i-1] + kappa * (theta - V[i-1]) * dt + sigma * np.sqrt(V[i-1] * dt) * np.random.normal(), 0)
S[i] = S[i-1] * np.exp((r - V[i-1]/2) * dt + np.sqrt(V[i-1] * dt) * np.random.normal())
plt.plot(np.arange(T * 250 + 1), S)
plt.title('Sample Path Does Not Hit Barrier')
plt.xlabel('Time')
plt.ylabel('Stock Price')
plt.show()
```
需要注意的是,这里使用的是欧式障碍期权,即只有到期时才判断是否触碰障碍,因此计算结果可能会存在一定的误差。如果需要更准确的定价,可以使用美式障碍期权,即在到期前任意时刻都可能触碰障碍。