I am facing ValueError: too many values to unpack (expected 2) while optimizing parameters of a system of ODEs using solve_ivp. In fact I get the same error when I tried to use solve_ivp instead of odeint in this SO answer, which you may use as a minimal working example since it has the same problem as far as I am concerned. The only changes I made to that code is swap positions of y, t in arguments for f and similarly while solving it using solve_ivp like so: x = solve_ivp(f, t, x0, args=(paras,)) instead of using odeint in g
Here's the full code for the sake of convenience:
# import libraries
import numpy as np
import matplotlib.pyplot as plt
from scipy.integrate import odeint, solve_ivp
from lmfit import minimize, Parameters, Parameter, report_fit
def f(t, y, paras):
"""
Your system of differential equations
"""
x1 = y[0]
x2 = y[1]
x3 = y[2]
try:
k0 = paras['k0'].value
k1 = paras['k1'].value
except KeyError:
k0, k1 = paras
# the model equations
f0 = -k0 * x1
f1 = k0 * x1 - k1 * x2
f2 = k1 * x2
return [f0, f1, f2]
def g(t, x0, paras):
"""
Solution to the ODE x'(t) = f(t,x,k) with initial condition x(0) = x0
"""
x = solve_ivp(f, t, x0, args=(paras,))
return x
def residual(paras, t, data):
"""
compute the residual between actual data and fitted data
"""
x0 = paras['x10'].value, paras['x20'].value, paras['x30'].value
model = g(t, x0, paras)
# you only have data for one of your variables
x2_model = model[:, 1]
return (x2_model - data).ravel()
# initial conditions
x10 = 5.
x20 = 0
x30 = 0
y0 = [x10, x20, x30]
# measured data
t_measured = np.linspace(0, 9, 10)
x2_measured = np.array([0.000, 0.416, 0.489, 0.595, 0.506, 0.493, 0.458, 0.394, 0.335, 0.309])
plt.figure()
plt.scatter(t_measured, x2_measured, marker='o', color='b', label='measured data', s=75)
# set parameters including bounds; you can also fix parameters (use vary=False)
params = Parameters()
params.add('x10', value=x10, vary=False)
params.add('x20', value=x20, vary=False)
params.add('x30', value=x30, vary=False)
params.add('k0', value=0.2, min=0.0001, max=2.)
params.add('k1', value=0.3, min=0.0001, max=2.)
# fit model
result = minimize(residual, params, args=(t_measured, x2_measured), method='leastsq') # leastsq nelder
# check results of the fit
data_fitted = g(np.linspace(0., 9., 100), y0, result.params)
# plot fitted data
plt.plot(np.linspace(0., 9., 100), data_fitted[:, 1], '-', linewidth=2, color='red', label='fitted data')
plt.legend()
plt.xlim([0, max(t_measured)])
plt.ylim([0, 1.1 * max(data_fitted[:, 1])])
# display fitted statistics
report_fit(result)
plt.show()
Here's the error traceback:
---------------------------------------------------------------------------
ValueError Traceback (most recent call last)
/home/swami/work/scrap/lmfit_example1.ipynb Cell 4 in <cell line: 67>()
64 params.add('k1', value=0.3, min=0.0001, max=2.)
66 # fit model
---> 67 result = minimize(residual, params, args=(t_measured, x2_measured), method='leastsq') # leastsq nelder
68 # check results of the fit
69 data_fitted = g(np.linspace(0., 9., 100), y0, result.params)
File ~/miniconda3/envs/dynamical/lib/python3.10/site-packages/lmfit/minimizer.py:2600, in minimize(fcn, params, method, args, kws, iter_cb, scale_covar, nan_policy, reduce_fcn, calc_covar, max_nfev, **fit_kws)
2460 """Perform the minimization of the objective function.
2461
2462 The minimize function takes an objective function to be minimized,
(...)
2594
2595 """
2596 fitter = Minimizer(fcn, params, fcn_args=args, fcn_kws=kws,
2597 iter_cb=iter_cb, scale_covar=scale_covar,
2598 nan_policy=nan_policy, reduce_fcn=reduce_fcn,
2599 calc_covar=calc_covar, max_nfev=max_nfev, **fit_kws)
-> 2600 return fitter.minimize(method=method)
File ~/miniconda3/envs/dynamical/lib/python3.10/site-packages/lmfit/minimizer.py:2369, in Minimizer.minimize(self, method, params, **kws)
2366 if (key.lower().startswith(user_method) or
2367 val.lower().startswith(user_method)):
2368 kwargs['method'] = val
-> 2369 return function(**kwargs)
File ~/miniconda3/envs/dynamical/lib/python3.10/site-packages/lmfit/minimizer.py:1693, in Minimizer.leastsq(self, params, max_nfev, **kws)
1691 result.call_kws = lskws
1692 try:
-> 1693 lsout = scipy_leastsq(self.__residual, variables, **lskws)
1694 except AbortFitException:
1695 pass
File ~/.local/lib/python3.10/site-packages/scipy/optimize/_minpack_py.py:410, in leastsq(func, x0, args, Dfun, full_output, col_deriv, ftol, xtol, gtol, maxfev, epsfcn, factor, diag)
408 if not isinstance(args, tuple):
409 args = (args,)
--> 410 shape, dtype = _check_func('leastsq', 'func', func, x0, args, n)
411 m = shape[0]
413 if n > m:
File ~/.local/lib/python3.10/site-packages/scipy/optimize/_minpack_py.py:24, in _check_func(checker, argname, thefunc, x0, args, numinputs, output_shape)
22 def _check_func(checker, argname, thefunc, x0, args, numinputs,
23 output_shape=None):
---> 24 res = atleast_1d(thefunc(*((x0[:numinputs],) + args)))
25 if (output_shape is not None) and (shape(res) != output_shape):
26 if (output_shape[0] != 1):
File ~/miniconda3/envs/dynamical/lib/python3.10/site-packages/lmfit/minimizer.py:586, in Minimizer.__residual(self, fvars, apply_bounds_transformation)
583 self.result.success = False
584 raise AbortFitException(f"fit aborted: too many function evaluations {self.max_nfev}")
--> 586 out = self.userfcn(params, *self.userargs, **self.userkws)
588 if callable(self.iter_cb):
589 abort = self.iter_cb(params, self.result.nfev, out,
590 *self.userargs, **self.userkws)
/home/swami/work/scrap/lmfit_example1.ipynb Cell 4 in residual(paras, t, data)
33 """
34 compute the residual between actual data and fitted data
35 """
37 x0 = paras['x10'].value, paras['x20'].value, paras['x30'].value
---> 38 model = g(t, x0, paras)
40 # you only have data for one of your variables
41 x2_model = model[:, 1]
/home/swami/work/scrap/lmfit_example1.ipynb Cell 4 in g(t, x0, paras)
23 def g(t, x0, paras):
24 """
25 Solution to the ODE x'(t) = f(t,x,k) with initial condition x(0) = x0
26 """
---> 27 x = solve_ivp(f, t, x0, args=(paras,))
28 return x
File ~/.local/lib/python3.10/site-packages/scipy/integrate/_ivp/ivp.py:512, in solve_ivp(fun, t_span, y0, method, t_eval, dense_output, events, vectorized, args, **options)
507 if method not in METHODS and not (
508 inspect.isclass(method) and issubclass(method, OdeSolver)):
509 raise ValueError("`method` must be one of {} or OdeSolver class."
510 .format(METHODS))
--> 512 t0, tf = map(float, t_span)
514 if args is not None:
515 # Wrap the user's fun (and jac, if given) in lambdas to hide the
516 # additional parameters. Pass in the original fun as a keyword
517 # argument to keep it in the scope of the lambda.
518 try:
ValueError: too many values to unpack (expected 2)
Any idea what the problem might be?
