假设我有一个函数f1,它有一个自变量x和一个参数a。反过来,函数f1需要一个参数m,该参数 m 可以通过最小化第二个函数f2来获得,并且该函数需要变量x作为参数。
我用任何值调用f1 ( x , a )都没有问题,但是如果我尝试适应curve_fit
它,它会抛出 a ValueError: setting an array element with a sequence
。我猜它这样做是因为它分配的不是单个x值,而是整个数组。
有没有办法让curve_fit
工作?
玩具代码,检查问题:
import matplotlib.pyplot as plt
from scipy.optimize import minimize
from scipy.optimize import curve_fit
from random import uniform
def f1(x,a):
m = minimize(f2,0.,args=(x)).x[0]
return (x + a + m)**2
def f2(theta,a):
return (theta-a)**2
xs = np.linspace(-3.,3.,10)
ys = [round(f1(x,0.5)-uniform(-1,1),2) for x in xs]
plt.scatter(xs,ys)
popt,pcov = curve_fit(f1, xs, ys, [0.3])
plt.plot(xs,[f1(x,popt[0]) for x in xs])
完整追溯:
ValueErrorTraceback (most recent call last)
<ipython-input-87-5db669cee585> in <module>()
10 plt.scatter(xs,ys)
11 plt.plot(xs,[f1(x,0.5) for x in xs])
---> 12 popt,pcov = curve_fit(f1, xs, ys, [0.3])
13 plt.plot(xs,[f1(x,popt[0]) for x in xs])
/usr/lib64/python2.7/site-packages/scipy/optimize/minpack.pyc in curve_fit(f, xdata, ydata, p0, sigma, absolute_sigma, check_finite, bounds, method, jac, **kwargs)
749 # Remove full_output from kwargs, otherwise we're passing it in twice.
750 return_full = kwargs.pop('full_output', False)
--> 751 res = leastsq(func, p0, Dfun=jac, full_output=1, **kwargs)
752 popt, pcov, infodict, errmsg, ier = res
753 cost = np.sum(infodict['fvec'] ** 2)
/usr/lib64/python2.7/site-packages/scipy/optimize/minpack.pyc in leastsq(func, x0, args, Dfun, full_output, col_deriv, ftol, xtol, gtol, maxfev, epsfcn, factor, diag)
381 if not isinstance(args, tuple):
382 args = (args,)
--> 383 shape, dtype = _check_func('leastsq', 'func', func, x0, args, n)
384 m = shape[0]
385 if n > m:
/usr/lib64/python2.7/site-packages/scipy/optimize/minpack.pyc in _check_func(checker, argname, thefunc, x0, args, numinputs, output_shape)
25 def _check_func(checker, argname, thefunc, x0, args, numinputs,
26 output_shape=None):
---> 27 res = atleast_1d(thefunc(*((x0[:numinputs],) + args)))
28 if (output_shape is not None) and (shape(res) != output_shape):
29 if (output_shape[0] != 1):
/usr/lib64/python2.7/site-packages/scipy/optimize/minpack.pyc in func_wrapped(params)
461 if transform is None:
462 def func_wrapped(params):
--> 463 return func(xdata, *params) - ydata
464 elif transform.ndim == 1:
465 def func_wrapped(params):
<ipython-input-87-5db669cee585> in f1(x, a)
1 from scipy.optimize import minimize
2 def f1(x,a):
----> 3 m = minimize(f2,0.,args=(x)).x[0]
4 return (x + a + m)**2
5 def f2(theta,a):
/usr/lib64/python2.7/site-packages/scipy/optimize/_minimize.pyc in minimize(fun, x0, args, method, jac, hess, hessp, bounds, constraints, tol, callback, options)
595 return _minimize_cg(fun, x0, args, jac, callback, **options)
596 elif meth == 'bfgs':
--> 597 return _minimize_bfgs(fun, x0, args, jac, callback, **options)
598 elif meth == 'newton-cg':
599 return _minimize_newtoncg(fun, x0, args, jac, hess, hessp, callback,
/usr/lib64/python2.7/site-packages/scipy/optimize/optimize.pyc in _minimize_bfgs(fun, x0, args, jac, callback, gtol, norm, eps, maxiter, disp, return_all, **unknown_options)
961 else:
962 grad_calls, myfprime = wrap_function(fprime, args)
--> 963 gfk = myfprime(x0)
964 k = 0
965 N = len(x0)
/usr/lib64/python2.7/site-packages/scipy/optimize/optimize.pyc in function_wrapper(*wrapper_args)
291 def function_wrapper(*wrapper_args):
292 ncalls[0] += 1
--> 293 return function(*(wrapper_args + args))
294
295 return ncalls, function_wrapper
/usr/lib64/python2.7/site-packages/scipy/optimize/optimize.pyc in approx_fprime(xk, f, epsilon, *args)
721
722 """
--> 723 return _approx_fprime_helper(xk, f, epsilon, args=args)
724
725
/usr/lib64/python2.7/site-packages/scipy/optimize/optimize.pyc in _approx_fprime_helper(xk, f, epsilon, args, f0)
661 ei[k] = 1.0
662 d = epsilon * ei
--> 663 grad[k] = (f(*((xk + d,) + args)) - f0) / d[k]
664 ei[k] = 0.0
665 return grad
ValueError: setting an array element with a sequence.