0

如果我使用泊松分布对非中心卡方分布进行采样,我无法更改大小,只能输入平均值“nc / 2”(我必须设置 size = 1 否则它也会返回相同的错误) :

n = np.random.poisson(nc / 2, 1)  # generates a random variable from the poisson distribution with
# mean: non-centrality parameter / 2
x[t] = c * mp.nsum(lambda i: np.random.standard_normal() ** 2, [0, v + 2 * n])

如果我尝试将大小增加到正在运行的模拟数量

n = np.random.poisson(nc / 2, simulations)

其中模拟 = 10000,我收到:

“ValueError:具有多个元素的数组的真值不明确。使用 a.any() 或 a.all()”

使用 1 次模拟运行代码会产生一个期望的结果,并且每次运行都会产生另一个随机路径。

在 10,000 次模拟下创建的图表,大小 = 1

但是,有必要让图形由模拟的每次迭代确定的路径组成。在不同的条件下,非中心卡方分布由代码确定:

x[t] = c * ((np.random.standard_normal(simulations) + nc ** 0.5) ** 2 + mp.nsum(
                lambda i: np.random.standard_normal(simulations) ** 2, [0, v - 1]))

这确实产生了预期的结果

由上面的代码行生成的图表

尽管无法更改泊松分布的大小,我如何获得 x[t] 的不同路径(即,10,000 个模拟中的每一个都没有相同的路径)

如果需要:

import numpy as np
import matplotlib.pyplot as plt
from scipy.stats import stats
import mpmath as mp

T = 1
beta = 1.5
x0 = 0.05
q = 0
mu = x0 - q
alpha = - (2 - beta) * mu
sigma0 = 0.1
sigma = (2 - beta) * sigma0
b = - (1 - beta) / (2 * mu) * sigma ** 2
simulations = 10000
M = 50
dt = T / M


def srd_sampled_nxc2():
    x = np.zeros((M + 1, simulations))
    x[0] = x0
    for t in range(1, M + 1):
        v = 4 * b * alpha / sigma ** 2
        c = (sigma ** 2 * (1 - np.exp(-alpha * dt))) / (4 * alpha)
        nc = np.exp(-alpha * dt) / c * x[t - 1]  # The non-centrality parameter lambda
        if v > 1:
            x[t] = c * ((np.random.standard_normal(simulations) + nc ** 0.5) ** 2 + mp.nsum(
                lambda i: np.random.standard_normal(simulations) ** 2, [0, v - 1]))
        else:
            n = np.random.poisson(nc / 2, 1)
            x[t] = c * mp.nsum(lambda i: np.random.standard_normal() ** 2, [0, v + 2 * n])
    return x


x1 = srd_sampled_nxc2()


plt.figure(figsize=(10, 6))
plt.plot(x1[:, :10], lw=1)
plt.xlabel('time')
plt.ylabel('index')
plt.show()
4

1 回答 1

0

我已经意识到大于 1 的变量 beta 会创建一个负 v 和一个非常大的 nc。由于 v 不能变为正数,因此无法创建分布,因此没有任何东西可以填充数组。我的印象是必须将 b 设为正数,从而解决负数 v 并允许程序运行。

于 2019-11-24T11:52:53.743 回答