output[0] = y0 ValueError: 使用序列设置数组元素

Question

我已经为此苦苦挣扎了几天。我正在尝试估计分段高斯函数的密度。谁能告诉我为什么我现在收到错误

TypeError:     output[0] = y0

ValueError: 使用序列设置数组元素。

它发生在这一行：

Zero_RG = integrate.romberg(gaussian(q,x,mu,sigma), Q1, Q2).`

这是脚本：

import numpy as np
import sympy as sp
from sympy import *
from scipy import integrate
from sympy import Integral, log, exp, sqrt, pi
import matplotlib.pyplot as plt 
from scipy.stats import norm, gaussian_kde
from quantecon import LAE

from sympy import symbols
var('Q1 Q2 x q sigma mu')
#q= symbols('q')

## == Define parameters == #
mu=80
sigma=20
b=0.2
Q=80
Q1=Q*(1-b)
Q2=Q*(1+b)
d = (sigma*np.sqrt(2*np.pi))
phi = norm()
n = 500

def p(x, y):
    x, y = np.array(x, dtype=float), np.array(y, dtype=float)
    Positive_RG = norm.pdf(x-y+Q1, mu, sigma)
    print('Positive_R = ', Positive_RG)
    Negative_RG = norm.pdf(x-y+Q2, mu, sigma) 
    print('Negative_RG = ', Negative_RG)
    gaussian = lambda q,x,mu,sigma: 1/(sigma*np.sqrt(2*np.pi))*np.exp(-(x+q-mu)**2 /(2*sigma**2))
    wrapped_gaussian = lambda q: gaussian(q, x, mu, sigma)
    Zero_RG = integrate.romberg(wrapped_gaussian, Q1, Q2)
    print('pdf',gaussian)
    #Zero_RG = scipy.integrate.quad(norm.pdf(x + q, mu, sigma))
#    Int_zerocase= lambda q: norm.pdf(x + q, u, sigma)
#    Zero_RG = scipy.integrate.quad(Int_zerocase, Q1, Q2)
#    print(Zero_RG)

    if y>0.0 and x -y>=-Q1:
        #print('printA', Positive_RG)
        return Positive_RG
    elif y<0.0 and x -y>=-Q2:
        #print('printC', Negative_RG)
        return Negative_RG
    elif y==0.0 and x >=-Q1:
        #print('printB', Zero_RG)
        return Zero_RG
    return 0.0 


Z = phi.rvs(n)
X = np.empty(n)
for t in range(n-1):
    X[t+1] = X[t] + Z[t]
    #X[t+1] = np.abs(X[t]) + Z[t]
psi_est = LAE(p, X)
k_est = gaussian_kde(X)

fig, ax = plt.subplots(figsize=(10,7))
ys = np.linspace(-200.0, 200.0, 200)
ax.plot(ys, psi_est(ys), 'g-', lw=2, alpha=0.6, label='look ahead estimate')
ax.plot(ys, k_est(ys), 'k-', lw=2, alpha=0.6, label='kernel based estimate')
ax.legend(loc='upper left')
plt.show()

Answer 1

romberg 的文档非常清楚，第一个参数是一个集成的单个变量的函数。

首先，一个小问题。优先使用np.exp 而不是np.e**。

在 Python 中，表达式

gaussian = lambda q,x,mu,sigma: 1/(sigma*np.sqrt(2*np.pi))*np.exp(-(x+q-mu)**2 /(2*sigma**2))

将gaussian 设置为一个有四个参数的函数。表达式gaussian(q, x, mu, sigma) 只是该函数的返回值。

有两种方法可以将所需的参数传递给romberg。最简单的方法是使用args 参数将三个附加参数作为元组传递：

Zero_RG = integrate.romberg(gaussian, Q1, Q2, args=(x,mu,sigma))

另一种方法是创建一个包装函数，它将为您传递参数：

wrapped_gaussian = lambda q: gaussian(q, x, mu, sigma)
Zero_RG = integrate.romberg(wrapped_gaussian, Q1, Q2)

我会推荐第一种方法，因为它使用现有机制，并且更短且更易于阅读。