用sobol序列准随机标准法数在Python中进行Monte模拟，得到了错误的值
用sobol序列准随机标准法数在Python中进行Monte模拟，得到了错误的值
提问于 2019-04-22 03:18:06
我试着用准随机标准正规数进行蒙特卡罗模拟。我知道我们可以使用Sobol序列生成一致数，然后用概率积分变换将它们转换成标准的法线数。我的代码给出了模拟资产路径的不切实际的值：
import sobol_seq
import numpy as np
from scipy.stats import norm
def i4_sobol_generate_std_normal(dim_num, n, skip=1):
"""
Generates multivariate standard normal quasi-random variables.
Parameters:
Input, integer dim_num, the spatial dimension.
Input, integer n, the number of points to generate.
Input, integer SKIP, the number of initial points to skip.
Output, real np array of shape (n, dim_num).
"""
sobols = sobol_seq.i4_sobol_generate(dim_num, n, skip)
normals = norm.ppf(sobols)
return normals
def GBM(Ttm, TradingDaysInAYear, NoOfPaths, UnderlyingPrice, RiskFreeRate, Volatility):
dt = float(Ttm) / TradingDaysInAYear
paths = np.zeros((TradingDaysInAYear + 1, NoOfPaths), np.float64)
paths[0] = UnderlyingPrice
for t in range(1, TradingDaysInAYear + 1):
rand = i4_sobol_generate_std_normal(1, NoOfPaths)
lRand = []
for i in range(len(rand)):
a = rand[i][0]
lRand.append(a)
rand = np.array(lRand)
paths[t] = paths[t - 1] * np.exp((RiskFreeRate - 0.5 * Volatility ** 2) * dt + Volatility * np.sqrt(dt) * rand)
return paths
GBM(1, 252, 8, 100., 0.05, 0.5)
array([[1.00000000e+02, 1.00000000e+02, 1.00000000e+02, ...,
1.00000000e+02, 1.00000000e+02, 1.00000000e+02],
[9.99702425e+01, 1.02116774e+02, 9.78688323e+01, ...,
1.00978615e+02, 9.64128959e+01, 9.72154915e+01],
[9.99404939e+01, 1.04278354e+02, 9.57830834e+01, ...,
1.01966807e+02, 9.29544649e+01, 9.45085180e+01],
...,
[9.28295879e+01, 1.88049044e+04, 4.58249200e-01, ...,
1.14117599e+03, 1.08089096e-02, 8.58754653e-02],
[9.28019642e+01, 1.92029616e+04, 4.48483141e-01, ...,
1.15234371e+03, 1.04211828e-02, 8.34842557e-02],
[9.27743486e+01, 1.96094448e+04, 4.38925214e-01, ...,
1.16362072e+03, 1.00473641e-02, 8.11596295e-02]])不应该生成像8.11596295e-02这样的值，因此我认为代码中有错误。如果我使用的是从numpy库rand = np.random.standard_normal(NoOfPaths)中提取的标准法线，那么这个价格与黑斯科尔斯的价格相匹配。因此，我认为问题在于随机数发生器。价值8.11596295e-02指的是路径中的价格，价格不太可能从100 (初始价格)降到8.11596295e-02。
回答 2
Stack Overflow用户
回答已采纳
发布于 2019-04-22 17:28:23
sobol_seq中似乎有一个bug。Anaconda，python3.7，64位，Windows 10 x64，通过pip安装sobol_seq
pip install sobol_seq
# Name Version Build Channel
sobol-seq 0.1.2 pypi_0 pypi简单代码
print(sobol_seq.i4_sobol_generate(1, 1, 0))
print(sobol_seq.i4_sobol_generate(1, 1, 1))
print(sobol_seq.i4_sobol_generate(1, 1, 2))
print(sobol_seq.i4_sobol_generate(1, 1, 3))产出量
[[0.5]]
[[0.5]]
[[0.5]]
[[0.5]]来自src/sobol/sobol.html的代码，sobol_lib.py的行为是合理的(除了第一点)。
嗯，封闭的代码看起来可能会工作，把种子和采样数组放在一起。虽然很慢..。
import sobol_seq
import numpy as np
import matplotlib.pyplot as plt
from scipy.stats import norm
def i4_sobol_generate_std_normal(dim_num, seed, size=None):
"""
Generates multivariate standard normal quasi-random variables.
Parameters:
Input, integer dim_num, the spatial dimension.
Input, integer n, the number of points to generate.
Input, integer seed, initial seed
Output, real np array of shape (n, dim_num).
"""
if size is None:
q, seed = sobol_seq.i4_sobol(dim_num, seed)
normals = norm.ppf(q)
return (normals, seed)
if isinstance(size, int) or isinstance(size, np.int32) or isinstance(size, np.int64) or isinstance(size, np.int16):
rc = np.empty((dim_num, size))
for k in range(size):
q, seed = sobol_seq.i4_sobol(dim_num, seed)
rc[:,k] = norm.ppf(q)
return (rc, seed)
else:
raise ValueError("Size type is not recognized")
return None
seed = 1
x, seed = i4_sobol_generate_std_normal(1, seed)
print(x)
x, seed = i4_sobol_generate_std_normal(1, seed)
print(x)
seed = 1
x, seed = i4_sobol_generate_std_normal(1, seed, size=10)
print(x)
x, seed = i4_sobol_generate_std_normal(1, seed, size=1000)
print(x)
hist, bins = np.histogram(x, bins=20, range=(-2.5, 2.5), density=True)
plt.bar(bins[:-1], hist, width = 0.22, align='edge')
plt.show()这是这张照片
Stack Overflow用户
发布于 2021-07-21 09:35:06
作为将来的参考，我们在SciPy 1.7：scipy.stats.qmc.Sobol中添加了Sobol序列
页面原文内容由Stack Overflow提供。腾讯云小微IT领域专用引擎提供翻译支持
原文链接：
https://stackoverflow.com/questions/55788739
复制相关文章
相似问题
腾讯云开发者
