2017-03-01 89 views
0

我一直在为此奋斗了几天。我试图估计分段高斯函数的密度。谁能告诉我为什么现在我得到错误输出[0] = y0 ValueError:设置一个序列的数组元素

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() 
+2

“integeate.romberg”函数的第一个参数必须是一个函数。你给它'gausian(...)',用参数调用lambda的结果。 – hpaulj

+0

@hpaulj,仍然无法正常工作。我已经忍受了2天的尝试:-) .... –

回答

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) 

我会推荐拳头的方法,因为它使用现有的机制,以及更短,更容易阅读。

+0

这两种方式的工作,但创建另一个错误。例如在第二种情况下声明:Zero_RG = integrate.romberg(wrapped_gaussian,Q1,Q2),output [0] = y0 ValueError:使用序列设置数组元素。 –

+0

第一个案子呢? –

+1

另外,将'np.e **(...)'改为'np.exp(...)'。我已经更新了答案。 –

相关问题