如何将参数传递给其他函数(通常通过 scipy)?
How do I pass through arguments to other functions (generally and via scipy)?
我正在尝试通过 scipy minimize a function that outputs chi-square 并找到最适合高斯叠加的 mu、sigma、normc。
from math import exp
from math import pi
from scipy.integrate import quad
from scipy.optimize import minimize
from scipy.stats import chisquare
import numpy as np
# guess intitial values for minimized chi-square
mu, sigma = np.mean(mydata), np.std(mydata) # mydata is my data points
normc = 1/(sigma * (2*pi)**(1/2))
gauss = lambda x: normc * exp( (-1) * (x - mu)**2 / ( 2 * (sigma **2) ) ) # Gaussian Distribution
# assume I have pre-defined bin-boundaries as a list called binbound
def expvalperbin(binbound,mu,sigma,normc):
# calculates expectation value per bin
ans = []
for index in range(len(binbound)):
if index != len(binbound)-1:
ans.append( quad( gauss, binbound[index], binbound[index+1])[0] )
return ans
expvalguess = expvalperbin(binbound,mu,sig,normc)
obsval = countperbin(binbound,mydata)
arglist = [mu,sig,norm]
def chisquareopt(obslist,explist):
return chisquare(obslist,explist)[0]
chisquareguess = chisquareopt((obsval,expvalguess), expvalguess, args=arglist)
result = minimize( chisquareopt(obsval,expvalguess), chisquareguess )
print(result)
运行 此代码为我提供了此错误:
TypeError: chisquareopt() got an unexpected keyword argument 'args'
我有几个问题:
1) 如何编写函数以允许将参数传递到我的函数 chisquareopt?
2) 我如何判断 scipy 是否会优化给出最小卡方的参数 [mu, sigma, normc]?我如何从优化中找到这些参数?
3) 很难知道我是否在这里取得进展。我在正确的轨道上吗?
编辑:如果相关,我有一个输入 [mu, sigma, normc] 并输出子列表列表的函数,每个子列表包含 [mu, sigma, normc] 的可能组合(其中外部列表涵盖指定范围内所有可能的参数组合。
通常这些 scipy
函数将 args
元组值原封不动地传递给您的代码。我应该仔细检查代码,但是
minimize(myfunc, x0, args=(y,z))
def myfunc(x, y, z):
<do something>
minimize
获取变量 x
的当前值(标量或数组,取决于 x0
的样子)和 args
参数,并构造
args = tuple(x) + args
myfunc(*args)
换句话说,它将 args
元组与迭代变量连接起来,并将其传递给您的函数。因此任何中间函数定义都需要使用该模式。
为了说明,定义一个接受通用参数元组的函数。
In [665]: from scipy.optimize import minimize
In [666]: def myfunc(*args):
...: print(args)
...: return np.abs(args[0])**2
...:
In [667]: myfunc(1,2,3)
(1, 2, 3)
Out[667]: 1
In [668]: myfunc(2,2,3)
(2, 2, 3)
Out[668]: 4
In [669]: minimize(myfunc, 10, args=(2,3))
(array([ 10.]), 2, 3)
(array([ 10.00000001]), 2, 3)
(array([ 10.]), 2, 3)
(array([ 8.99]), 2, 3)
....
(array([-0.00000003]), 2, 3)
Out[669]:
fun: 1.7161984122524196e-15
hess_inv: array([[ 0.50000001]])
jac: array([-0.00000007])
message: 'Optimization terminated successfully.'
nfev: 15
nit: 4
njev: 5
status: 0
success: True
x: array([-0.00000004])
(删除了关于正在最小化哪些参数的混淆的讨论。请参阅其他答案或我的编辑历史记录)
我已经稍微简化了您的问题,以便您对问题 2) 有一个想法。
首先,我已将您的直方图 obslist
和数据点数 N
硬编码为全局变量(这稍微简化了函数签名)。其次,我在 expvalperbin
中对 bin 边界进行了硬编码,假设 9 个 bin 具有固定宽度 5
并且第一个 bin 从 30
开始(因此直方图的范围从 30 到 75)。
第三,我使用 optimize.fmin
(Nelder-Mead) 而不是 optimize.minimize
。使用 fmin
而不是 minimize
的原因是通过 args=(x,y)
传递附加参数似乎在附加参数保持在固定值的意义上不起作用第一次调用。这不是你想要的:你想同时优化 mu
和 sigma
。
鉴于这些简化,我们有以下(肯定非常非 pythonic)脚本:
from math import exp
from math import pi
from scipy.integrate import quad
from scipy.optimize import fmin
from scipy.stats import chisquare
obslist = [12, 51, 144, 268, 264, 166, 75, 18, 2] # histogram, 1000 observations
N = 1000 # no. of data points
def gauss(x, mu, sigma):
return 1/(sigma * (2*pi)**(1/2)) * exp( (-1) * (x - mu)**2 / ( 2 * (sigma **2) ) )
def expvalperbin(mu, sigma):
e = []
# hard-coded bin boundaries
for i in range(30, 75, 5):
e.append(quad(gauss, i, i + 5, args=(mu, sigma))[0] * N)
return e
def chisquareopt(args):
# args[0] = mu
# args[1] = sigma
return chisquare(obslist, expvalperbin(args[0], args[1]))[0]
# initial guesses
initial_mu = 35.5
initial_sigma = 14
result = fmin(chisquareopt, [initial_mu, initial_sigma])
print(result)
Optimization terminated successfully.
Current function value: 2.010966
Iterations: 49
Function evaluations: 95
[ 50.57590239 7.01857529]
顺便说一句,obslist
直方图是来自 N(50.5, 7.0)
正态分布的 1000 点随机样本。请记住,这是我的第一个 Python 代码行,所以请不要根据风格来评判我。我只是想给你一个关于问题的一般结构的想法。
我正在尝试通过 scipy minimize a function that outputs chi-square 并找到最适合高斯叠加的 mu、sigma、normc。
from math import exp
from math import pi
from scipy.integrate import quad
from scipy.optimize import minimize
from scipy.stats import chisquare
import numpy as np
# guess intitial values for minimized chi-square
mu, sigma = np.mean(mydata), np.std(mydata) # mydata is my data points
normc = 1/(sigma * (2*pi)**(1/2))
gauss = lambda x: normc * exp( (-1) * (x - mu)**2 / ( 2 * (sigma **2) ) ) # Gaussian Distribution
# assume I have pre-defined bin-boundaries as a list called binbound
def expvalperbin(binbound,mu,sigma,normc):
# calculates expectation value per bin
ans = []
for index in range(len(binbound)):
if index != len(binbound)-1:
ans.append( quad( gauss, binbound[index], binbound[index+1])[0] )
return ans
expvalguess = expvalperbin(binbound,mu,sig,normc)
obsval = countperbin(binbound,mydata)
arglist = [mu,sig,norm]
def chisquareopt(obslist,explist):
return chisquare(obslist,explist)[0]
chisquareguess = chisquareopt((obsval,expvalguess), expvalguess, args=arglist)
result = minimize( chisquareopt(obsval,expvalguess), chisquareguess )
print(result)
运行 此代码为我提供了此错误:
TypeError: chisquareopt() got an unexpected keyword argument 'args'
我有几个问题:
1) 如何编写函数以允许将参数传递到我的函数 chisquareopt?
2) 我如何判断 scipy 是否会优化给出最小卡方的参数 [mu, sigma, normc]?我如何从优化中找到这些参数?
3) 很难知道我是否在这里取得进展。我在正确的轨道上吗?
编辑:如果相关,我有一个输入 [mu, sigma, normc] 并输出子列表列表的函数,每个子列表包含 [mu, sigma, normc] 的可能组合(其中外部列表涵盖指定范围内所有可能的参数组合。
通常这些 scipy
函数将 args
元组值原封不动地传递给您的代码。我应该仔细检查代码,但是
minimize(myfunc, x0, args=(y,z))
def myfunc(x, y, z):
<do something>
minimize
获取变量 x
的当前值(标量或数组,取决于 x0
的样子)和 args
参数,并构造
args = tuple(x) + args
myfunc(*args)
换句话说,它将 args
元组与迭代变量连接起来,并将其传递给您的函数。因此任何中间函数定义都需要使用该模式。
为了说明,定义一个接受通用参数元组的函数。
In [665]: from scipy.optimize import minimize
In [666]: def myfunc(*args):
...: print(args)
...: return np.abs(args[0])**2
...:
In [667]: myfunc(1,2,3)
(1, 2, 3)
Out[667]: 1
In [668]: myfunc(2,2,3)
(2, 2, 3)
Out[668]: 4
In [669]: minimize(myfunc, 10, args=(2,3))
(array([ 10.]), 2, 3)
(array([ 10.00000001]), 2, 3)
(array([ 10.]), 2, 3)
(array([ 8.99]), 2, 3)
....
(array([-0.00000003]), 2, 3)
Out[669]:
fun: 1.7161984122524196e-15
hess_inv: array([[ 0.50000001]])
jac: array([-0.00000007])
message: 'Optimization terminated successfully.'
nfev: 15
nit: 4
njev: 5
status: 0
success: True
x: array([-0.00000004])
(删除了关于正在最小化哪些参数的混淆的讨论。请参阅其他答案或我的编辑历史记录)
我已经稍微简化了您的问题,以便您对问题 2) 有一个想法。
首先,我已将您的直方图 obslist
和数据点数 N
硬编码为全局变量(这稍微简化了函数签名)。其次,我在 expvalperbin
中对 bin 边界进行了硬编码,假设 9 个 bin 具有固定宽度 5
并且第一个 bin 从 30
开始(因此直方图的范围从 30 到 75)。
第三,我使用 optimize.fmin
(Nelder-Mead) 而不是 optimize.minimize
。使用 fmin
而不是 minimize
的原因是通过 args=(x,y)
传递附加参数似乎在附加参数保持在固定值的意义上不起作用第一次调用。这不是你想要的:你想同时优化 mu
和 sigma
。
鉴于这些简化,我们有以下(肯定非常非 pythonic)脚本:
from math import exp
from math import pi
from scipy.integrate import quad
from scipy.optimize import fmin
from scipy.stats import chisquare
obslist = [12, 51, 144, 268, 264, 166, 75, 18, 2] # histogram, 1000 observations
N = 1000 # no. of data points
def gauss(x, mu, sigma):
return 1/(sigma * (2*pi)**(1/2)) * exp( (-1) * (x - mu)**2 / ( 2 * (sigma **2) ) )
def expvalperbin(mu, sigma):
e = []
# hard-coded bin boundaries
for i in range(30, 75, 5):
e.append(quad(gauss, i, i + 5, args=(mu, sigma))[0] * N)
return e
def chisquareopt(args):
# args[0] = mu
# args[1] = sigma
return chisquare(obslist, expvalperbin(args[0], args[1]))[0]
# initial guesses
initial_mu = 35.5
initial_sigma = 14
result = fmin(chisquareopt, [initial_mu, initial_sigma])
print(result)
Optimization terminated successfully.
Current function value: 2.010966
Iterations: 49
Function evaluations: 95
[ 50.57590239 7.01857529]
顺便说一句,obslist
直方图是来自 N(50.5, 7.0)
正态分布的 1000 点随机样本。请记住,这是我的第一个 Python 代码行,所以请不要根据风格来评判我。我只是想给你一个关于问题的一般结构的想法。