可以使用字典作为参数引用在PyMinuit中进行最小化吗
有没有可能通过传递一个参数字典来执行 PyMinuit 的函数最小化?
比如,通常使用 PyMinuit 的方式是这样的:
def f(x, a, b): return a + b*x
def chi2(a,b):
c2 = 0.
for x, y, yerr in data:
c2 += (f(x, a, b) - y)**2 / yerr**2
return c2
m = minuit.Minuit(chi2)
m.migrad()
根据这个问题,我了解到 PyMinuit 使用了一种叫做“反射”的技术来确定参数 x 和 y(不过我不太明白这是什么意思)。理想情况下,我希望能够做类似这样的事情:
p = dict()
p['x'] = 0.
p['y'] = 0.
def f(x,a,b): return a + b*x
def chi2():
c2 = 0.
for x, y, yerr in data:
c2 += (f(x, a, b) - y)**2 / yerr**2
return c2
m = minuit.Minuit(chi2,**p)
m.migrad()
甚至可以这样:
p = <dictionary of parameters + initial values>
model = <list containing strings representing functions e.g. 'a*b+a**2*x'>
data = x, y, yerr, model
def chi2():
c2 = 0.
for x, y, yerr, model in data:
c2 += (eval(model,{"__builtins__":None},p) - y)**2 / yerr**2
return c2
m = minuit.Minuit(chi2)
m.migrad()
我在谷歌小组的问题页面上看到一个类似问题的解决方法,他们通过整数输入生成了“假代码”和“假函数”(点击链接查看)。我尝试用我的字典 p 做了类似的事情:
class fake_code:
def __init__(self,p):
self.co_argcount = len(p)
self.co_varnames = tuple(p.keys())
print tuple(p.keys())
class fake_function:
def __init__(self,p):
self.func_code = fake_code(p)
def __call__(self,*args):
c2 = 0.
print args
for x, y, yerr in data:
c2 += (f(x, a, b) - y)**2 / yerr**2
return c2
但由于某种原因,所有参数都被标记为“固定”,我似乎无法将它们“解锁”。
我觉得这样做应该是可行的,但我对 Python 了解不够,不知道这是否是最佳方法,甚至是否应该尝试。如果有人能对此提供一些见解,我将非常感激。:)
3 个回答
1
可能回答得有点晚。你可以试试这个iminuit。我写这个是因为在其他工具中缺少这个特定的功能。
http://iminuit.github.com/iminuit/
这里有个例子,教你怎么写一个通用的成本函数:
不过,虽然写一个chi^2/似然函数很简单,但在probfit中已经为你写好了。
1
下面的内容大部分没有经过测试,我通常会避免这样做,但为了更好地向你解释我在评论中提到的简化方法,我决定做个例外。这种方法是基于这里第一个示例的。
import minuit
def minuit_call(func, **kwargs):
CALL_TEMPLATE = "minuit.Minuit({0.__name__}, {1})"
arg_str = ', '.join('{}={}'.format(k, v) for k,v in kwargs.iteritems())
return eval(CALL_TEMPLATE.format(func, arg_str))
def f(x, y):
return ((x-2) / 3)**2 + y**2 + y**4
m = minuit_call(f, x=0, y=0)
m.migrad()
正如你所看到的,使用的模板相当简单,创建它并不需要手动将要进行最小化的函数体中的任何代码翻译成格式化字符串。
1
好的,我不太喜欢自己回答自己的问题,但我想我找到了一个用 exec
的解决办法。如果在一个模板中定义 chi2
函数,并通过一个叫 make_chi_squared
的函数在运行时构建它,那就可以实现。下面是我想到的解决方案。
import minuit
import numpy
chi_squared_template = """
def chi_squared(%(params)s):
li = [%(params)s]
for i,para in enumerate(li):
p[l[i]] = para
return (((f(data_x, p) - data_y) / errors) ** 2).sum()
"""
l = ['a1','a2','a3','a4']
p = dict()
p['a1'] = 1.
p['a2'] = 1.
p['a3'] = 1.
p['a4'] = 1.
def make_chi_squared(f, data_x, data_y, errors):
params = ", ".join(l)
exec chi_squared_template % {"params": params}
return chi_squared
def f(x,p):
return eval('a1 + a2*x + a3*x**2 + a4*x**3',
{"__builtins__":locals()},
p)
data_x = numpy.arange(50)
errors = numpy.random.randn(50) * 0.3
data_y = data_x**3 + errors
chi_squared = make_chi_squared(f, data_x, data_y, errors)
m = minuit.Minuit(chi_squared)
m.printMode = 1
m.migrad()
print m.values
p = m.values
print p
这有点乱,我不确定这是不是处理这种问题的最佳方法,但它确实有效!