Constrained Example Error

[bhclowers]

As of version 0.4.2 the example provided on the front of the github page for constraints appears to throw and error (at least on my system Python 3.5.2). Any ideas would be appreciated.

from symfit import parameters, Fit, Equality, GreaterThan
x, y = parameters('x, y')
model = 2 * x * y + 2 * x - x**2 - 2 * y**2
constraints = [
    Equality(x**3, y),
    GreaterThan(y, 1),
]

fit = Fit(model, constraints=constraints)
fit_result = fit.execute()

Error output:

---------------------------------------------------------------------------
TypeError                                 Traceback (most recent call last)
<ipython-input-1-5bb7e724f8a2> in <module>()
      9 
     10 fit = Fit(model, constraints=constraints)
---> 11 fit_result = fit.execute()

~/.local/lib/python3.5/site-packages/symfit/core/fit.py in execute(self, **minimize_options)
   1401         :return: FitResults instance
   1402         """
-> 1403         minimizer_ans = self.minimizer.execute(**minimize_options)
   1404         try: # to build covariance matrix
   1405             cov_matrix = minimizer_ans.covariance_matrix
~/.local/lib/python3.5/site-packages/symfit/core/minimizers.py in execute(self, **minimize_options)
    207             method='SLSQP',
    208             bounds=self.bounds,
--> 209             **minimize_options
    210         )
    211 
~/.local/lib/python3.5/site-packages/symfit/core/minimizers.py in execute(self, **minimize_options)
    194 
    195     def execute(self, **minimize_options):
--> 196         return super(ScipyGradientMinimize, self).execute(jacobian=self.wrapped_jacobian, **minimize_options)
    197 
    198 
~/.local/lib/python3.5/site-packages/symfit/core/support.py in wrapped_func(*args, **kwargs)
    280                     else:
    281                         bound_args.arguments[param.name] = param.default
--> 282             return func(*bound_args.args, **bound_args.kwargs)
    283         return wrapped_func
    284 
~/.local/lib/python3.5/site-packages/symfit/core/minimizers.py in execute(self, bounds, jacobian, **minimize_options)
    129             constraints=self.constraints,
    130             jac=jacobian,
--> 131             **minimize_options
    132         )
    133         # Build infodic
/usr/local/lib/python3.5/dist-packages/scipy/optimize/_minimize.py in minimize(fun, x0, args, method, jac, hess, hessp, bounds, constraints, tol, callback, options)
    609     elif meth == 'slsqp':
    610         return _minimize_slsqp(fun, x0, args, jac, bounds,
--> 611                                constraints, callback=callback, **options)
    612     elif meth == 'trust-constr':
    613         return _minimize_trustregion_constr(fun, x0, args, jac, hess, hessp,
/usr/local/lib/python3.5/dist-packages/scipy/optimize/slsqp.py in _minimize_slsqp(func, x0, args, jac, bounds, constraints, maxiter, ftol, iprint, disp, eps, callback, **unknown_options)
    311     # meq, mieq: number of equality and inequality constraints
    312     meq = sum(map(len, [atleast_1d(c['fun'](x, *c['args']))
--> 313               for c in cons['eq']]))
    314     mieq = sum(map(len, [atleast_1d(c['fun'](x, *c['args']))
    315                for c in cons['ineq']]))
/usr/local/lib/python3.5/dist-packages/scipy/optimize/slsqp.py in <listcomp>(.0)
    311     # meq, mieq: number of equality and inequality constraints
    312     meq = sum(map(len, [atleast_1d(c['fun'](x, *c['args']))
--> 313               for c in cons['eq']]))
    314     mieq = sum(map(len, [atleast_1d(c['fun'](x, *c['args']))
    315                for c in cons['ineq']]))
~/.local/lib/python3.5/site-packages/symfit/core/minimizers.py in <lambda>(p, x, c)
    178                 'type': types[constraint.constraint_type],
    179                 # Assume the lhs is the equation.
--> 180                 'fun': lambda p, x, c: c(*(list(x.values()) + list(p)))[0],
    181                 # Assume the lhs is the equation.
    182                 'jac': lambda p, x, c: [component(*(list(x.values()) + list(p)))
~/.local/lib/python3.5/site-packages/symfit/core/fit.py in __call__(self, *args, **kwargs)
    238         Ans = namedtuple('Ans', [var.name for var in self])
    239         # return Ans(*[component(**bound_arguments.arguments) for component in self.numerical_components])
--> 240         return Ans(*self.eval_components(**bound_arguments.arguments))
    241 
    242     @property
~/.local/lib/python3.5/site-packages/symfit/core/fit.py in eval_components(self, *args, **kwargs)
    504         :return: lambda functions of each of the components in model_dict, to be used in numerical calculation.
    505         """
--> 506         return [expr(*args, **kwargs) for expr in self.numerical_components]
    507         # return [sympy_to_py(expr, self.independent_vars, self.params)(*args, **kwargs) for expr in self.values()]
    508 
~/.local/lib/python3.5/site-packages/symfit/core/fit.py in <listcomp>(.0)
    504         :return: lambda functions of each of the components in model_dict, to be used in numerical calculation.
    505         """
--> 506         return [expr(*args, **kwargs) for expr in self.numerical_components]
    507         # return [sympy_to_py(expr, self.independent_vars, self.params)(*args, **kwargs) for expr in self.values()]
    508 
/usr/local/lib/python3.5/dist-packages/sympy/utilities/lambdify.py in wrapper(*argsx, **kwargsx)
    442                 newargs = [asarray(i) if isinstance(i, integer_types + (float,
    443                     complex)) else i for i in argsx]
--> 444                 return funcarg(*newargs, **kwargsx)
    445             return wrapper
    446         func = array_wrap(func)
TypeError: <lambda>() got an unexpected keyword argument 'y_1'

[tBuLi]

Thanks for reporting it, I submitted a patch. I actually think this embarrassing bug was already there for longer, because the tests made for this example explicitly gave a name to the dependent variable, and therefore this bug was never triggered. The tests have been modified to make sure this doesn't happen again.

[bhclowers]

So after cloning the most recent update and running: python3 setup.py install --user I get the following result which still seems a bit odd. Did I do something wrong here?

from symfit import parameters, Fit, Equality, GreaterThan
x, y = parameters('x, y')
model = 2 * x * y + 2 * x - x**2 - 2 * y**2
constraints = [
    Equality(x**3, y),
    GreaterThan(y, 1),
]

fit = Fit(model, constraints=constraints)
fit_result = fit.execute()

****Output 1:

/home/user/.local/lib/python3.5/site-packages/symfit/core/fit.py:1654: RuntimeWarning: invalid value encountered in double_scalars
  return 1 - SS_res/SS_tot

****Output 2:

>>fit_result.params
OrderedDict([('x', nan), ('y', nan)])

[tBuLi]

Actually the readme has not been updated properly it appears, try negating the model since Fit always minimizes and here we need to maximize:

fit = Fit(- model, constraints=constraints)

This was forgotten when changing the API in 0.4.x, so thanks again for pointing it out!

[bhclowers]

After banging my head against this for a while, I still appear stuck. I can get the simple example to complete but when I try to extend the concept I am getting some odd behavior. In short, I have the following code and after executing the fit the error thrown is saying that that parameter (theta1) is not present. Is this a syntax/keyword argument issue? If so can you please nudge me in the right direction?

import numpy as np
from symfit import parameters, Fit, GreaterThan

phi1, phi2, theta1, theta2 = parameters('phi1, phi2, theta1, theta2')
x, y = variables('x, y')

model_dict = {y: (1+x*theta1+theta2*x**2)/(1+phi1*x*theta1+phi2*theta2*x**2)}
constraints = [GreaterThan(theta1, theta2)]

xdata = np.array([0., 0.000376, 0.000752, 0.0015  , 0.00301 , 0.00601 , 0.00902 ])
ydata=np.array([1., 1.07968041, 1.08990638, 1.12151629, 1.13068452,1.15484109, 1.19883952])

phi1.value = 0.845251484373516
phi1.fixed = True

phi2.value = 0.7105427053026403
phi2.fixed = True

fit = Fit(model_dict, x=xdata, y=ydata, constraints=constraints)# Feed in the observed x and y data

fit_result = fit.execute()

***Error thrown:

TypeError: missing a required argument: 'theta1'

[tBuLi]

Damn, you are good at finding bugs. This one seems to be caused by fixing the values of the phi's, if I remove that the fit runs.

So after some searching it seems that the bug is caused by some of the changes made in PR #151, so I'm also going to mention @pckroon to also get his eyes on this. And I see that this PR was made in order to fix issue #150 which you also made, so thank you for your patience ;).

The problem lies in this: the Minimizer.initial_guesses property returns the initial guesses for the non-fixed parameters. The objective has been partialed to expect only those arguments, the fixed ones have been provided by partial:

self.objective = partial(objective, **{p.name: p.value for p in self._fixed_params})

However, the constraints have never been adapted to this, so they are still len(self._fixed_params) short. Therefore, the solution is to also partial the constraints. However, scipy_constraints is a staticmethod so it is agnostic about which parameters are fixed. Some thought is needed.

[pckroon]

Why is it that every time someone uses my code they find bugs? :')

I think the architectural issue is fairly simple to solve here. ConstrainedMinimizer wraps the symfit style constraints to get rid of the fixed parameters, and the Scipy wrapper remains unchanged (similar to what we do with the jacobian). I'm fairly swamped with work at the moment though, so I'll leave the implementation to @tBuLi.

[tBuLi]

In principle that's a nice idea. But it turns out not to be quite that simple, because of the MRO. I will watch Super Considered Super again, and then try again :P.

(BaseMinimizer is needed to know which parameters are fixed, and then we could do a simple partial after. But ScipyMinimize's init is called before BaseMinimizer, so the most obvious solution doesn't work.)

[bhclowers]

@tBuLi Wish I could help on the other issue, however, if you were able to get it to run without the phi parameters being fixed, could you please send me the code you ran and output values? I am having trouble making that happen on my end. Even with those values not being fixed, I might be able to make this work for our particular purposes. However, it would be great to be able to used constrained parameters along with fixed ones.

[bhclowers]

I'm not sure this is related but in the latest version that I cloned even the basic examples seem to be throwing an error associated with an AttributeError:

from symfit import parameters, variables, Fit, Model

xdata = [1.0, 2.0, 3.0, 4.0, 5.0]
ydata = [2.3, 3.3, 4.1, 5.5, 6.7]
yerr = [0.1, 0.1, 0.1, 0.1, 0.1]

a, b = parameters('a, b')
x, y = variables('x, y')
model = Model({y: a * x + b})

fit = Fit(model, x=xdata, y=ydata, sigma_y=yerr)
fit_result = fit.execute()

Error:

/anaconda2/envs/py36/lib/python3.6/site-packages/sympy/printing/codeprinter.py in _print_Variable(self, expr)
    342 
    343     def _print_Variable(self, expr):
--> 344         return self._print(expr.symbol)
    345 
    346     def _print_Statement(self, expr):

AttributeError: 'Variable' object has no attribute 'symbol'

[tBuLi]

To your first comment: I literally ran your example except that I removed the lines saying they had to be fixed:

import numpy as np
from symfit import parameters, Fit, GreaterThan, variables

phi1, phi2, theta1, theta2 = parameters('phi1, phi2, theta1, theta2')
x, y = variables('x, y')

model_dict = {y: (1+x*theta1+theta2*x**2)/(1+phi1*x*theta1+phi2*theta2*x**2)}
constraints = [GreaterThan(theta1, theta2)]

xdata = np.array([0., 0.000376, 0.000752, 0.0015  , 0.00301 , 0.00601 , 0.00902 ])
ydata=np.array([1., 1.07968041, 1.08990638, 1.12151629, 1.13068452,1.15484109, 1.19883952])

phi1.value = 0.845251484373516
phi2.value = 0.7105427053026403

fit = Fit(model_dict, x=xdata, y=ydata, constraints=constraints)# Feed in the observed x and y data

fit_result = fit.execute()
print(fit_result)

Parameter Value        Standard Deviation
phi1      8.549886e-01 1.977125e-01
phi2      1.600773e+02 3.472350e+07
theta1    2.074183e+03 4.091964e+03
theta2    -1.557679e+01 3.396675e+06
Fitting status message: Optimization terminated successfully.
Number of iterations:   117
Regression Coefficient: 0.950091684637884

I leave it up to you to plot it to see if it makes sense, but judging from the Regression Coefficient it should be pretty good.

To your second comment: that code example is not complete enough because it doesn't contain a print statement, whereas the traceback indicates that it's a printing issue. So I think the code snippet you posted is not exactly what you actually ran? And like you already said, it's probably not a related issue, so please start a separate issue for that.

Let me know if fitting without fixing works on your end too!

P.s., as a quick and dirty workaround in the meantime you could also set very tight bounds on the phi's, so that they are practically fixed.

[bhclowers]

@tBuLi with the reversion of sympy I was able to fix the values of phi1 and phi2 in the above example and get the constraints to work out. Cheers.

[tBuLi]

Fixed with PR #168. Thank you for reporting this problem.