Skip to content

Commit

Permalink
Merge branch 'master' into fitresults_cleanup
Browse files Browse the repository at this point in the history
  • Loading branch information
Martin Roelfs committed May 18, 2019
2 parents e8db5eb + 63a9ceb commit 84dd102
Show file tree
Hide file tree
Showing 12 changed files with 297 additions and 262 deletions.
Expand Up @@ -29,15 +29,13 @@ def setUpClass(cls):
k = Parameter('k', 900)
x0 = Parameter('x0', 1.5)

# You can NOT do this in one go. Blame Sympy. Not my fault.
cls.k = k
cls.x0 = x0

model = {y: distr(x, k, x0)}
x_data = np.linspace(0, 2.5, 50)
y_data = model[y](x=x_data, k=1000, x0=1)
cls.guess = interactive_guess.InteractiveGuess(model, x=x_data, y=y_data)
# plt.close(cls.fit.fig)

def test_number_of_sliders(self):
self.assertEqual(len(self.guess._sliders), 2)
Expand Down
246 changes: 140 additions & 106 deletions symfit/core/fit.py

Large diffs are not rendered by default.

8 changes: 5 additions & 3 deletions symfit/core/fit_results.py
Expand Up @@ -266,10 +266,12 @@ def r_squared(model, fit_result, data):
:param data: data with which the fit was performed.
"""
# First filter out the dependent vars
y_is = [data[var] for var in model if var in data]
x_is = [value for var, value in data.items() if var.name in model.__signature__.parameters]
y_is = [data[var] for var in model.dependent_vars if var in data]
x_is = [data[var] for var in model.independent_vars if var in data]
y_bars = [np.mean(y_i) if y_i is not None else None for y_i in y_is]
f_is = model(*x_is, **fit_result.params)
f_is = model(*x_is, **fit_result.params)._asdict()
# f_is also contains the evaluated interdependent_vars, skip those.
f_is = [f_is[var] for var in model.dependent_vars]
SS_res = np.sum([np.sum((y_i - f_i)**2) for y_i, f_i in zip(y_is, f_is) if y_i is not None])
SS_tot = np.sum([np.sum((y_i - y_bar)**2) for y_i, y_bar in zip(y_is, y_bars) if y_i is not None])
return 1 - SS_res/SS_tot
2 changes: 1 addition & 1 deletion symfit/core/minimizers.py
Expand Up @@ -72,7 +72,7 @@ def _baseobjective_from_callable(self, func, objective_type=MinimizeModel):
connectivity_mapping={y: set(self.parameters)}
)
return objective_type(model,
data={y: None for y in model.dependent_vars})
data={})

@abc.abstractmethod
def execute(self, **options):
Expand Down
54 changes: 46 additions & 8 deletions symfit/core/objectives.py
Expand Up @@ -18,6 +18,8 @@ def __init__(self, model, data):
"""
self.model = model
self.data = data
# Compares the model with the data to see if they are compatible.
self._sanity_checking()

@cached_property
def dependent_data(self):
Expand Down Expand Up @@ -91,7 +93,8 @@ def _shape_of_dependent_data(self, model_output, param_level=0):
"""
shaped_result = []
n_params = len(self.model.params)
for dep_data, component in zip(self.dependent_data.values(), model_output):
for dep_var, component in zip(self.model.dependent_vars, model_output):
dep_data = self.dependent_data.get(dep_var, None)
if dep_data is not None:
if dep_data.shape == component.shape:
shaped_result.append(component)
Expand Down Expand Up @@ -147,6 +150,16 @@ def __eq__(self, other):
return False
return True

def _sanity_checking(self):
"""
Check if the model and the provided data are compatible. Raises a
TypeError when this is not the case.
"""
# Simply checking for existence of these dicts will raise an error if
# they cannot be built (because these are properties).
self.dependent_data
self.independent_data
self.sigma_data


@add_metaclass(abc.ABCMeta)
Expand Down Expand Up @@ -224,7 +237,8 @@ def __call__(self, ordered_parameters=[], **parameters):

# zip together the dependent vars and evaluated component
for y, ans in zip(self.model.dependent_vars, evaluated_func):
if self.dependent_data[y] is not None:
dep_data = self.dependent_data.get(y, None)
if dep_data is not None:
result.append(((self.dependent_data[y] - ans) / self.sigma_data[self.model.sigmas[y]]) ** 2)
if flatten_components: # Flattens *within* a component
result[-1] = result[-1].flatten()
Expand All @@ -242,7 +256,8 @@ def eval_jacobian(self, ordered_parameters=[], **parameters):
result = len(self.model.params) * [0.0]
for ans, y, row in zip(evaluated_func, self.model.dependent_vars,
evaluated_jac):
if self.dependent_data[y] is not None:
dep_data = self.dependent_data.get(y, None)
if dep_data is not None:
for index, component in enumerate(row):
result[index] += component * (
(self.dependent_data[y] - ans) / self.sigma_data[self.model.sigmas[y]] ** 2
Expand Down Expand Up @@ -284,7 +299,7 @@ def __call__(self, ordered_parameters=[], **parameters):

chi2 = [0 for _ in evaluated_func]
for index, (dep_var, dep_var_value) in enumerate(zip(self.model.dependent_vars, evaluated_func)):
dep_data = self.dependent_data[dep_var]
dep_data = self.dependent_data.get(dep_var, None)
if dep_data is not None:
sigma = self.sigma_data[self.model.sigmas[dep_var]]
chi2[index] += np.sum(
Expand Down Expand Up @@ -312,7 +327,7 @@ def eval_jacobian(self, ordered_parameters=[], **parameters):
result = 0
for var, f, jac_comp in zip(self.model.dependent_vars, evaluated_func,
evaluated_jac):
y = self.dependent_data[var]
y = self.dependent_data.get(var, None)
sigma_var = self.model.sigmas[var]
if y is not None:
sigma = self.sigma_data[sigma_var]
Expand Down Expand Up @@ -344,7 +359,7 @@ def eval_hessian(self, ordered_parameters=[], **parameters):
for var, f, jac_comp, hess_comp in zip(self.model.dependent_vars,
evaluated_func, evaluated_jac,
evaluated_hess):
y = self.dependent_data[var]
y = self.dependent_data.get(var, None)
sigma_var = self.model.sigmas[var]
if y is not None:
sigma = self.sigma_data[sigma_var]
Expand Down Expand Up @@ -391,7 +406,30 @@ def eval_hessian(self, ordered_parameters=[], **parameters):
) for comp in result]


class LogLikelihood(HessianObjective):
class BaseIndependentObjective(BaseObjective):
"""
Some objective functions dependent only on independent variables, not
dependent and sigma variables. In this case, sanity checking is greatly
simplified.
"""
@cached_property
def dependent_data(self):
"""
:return: Empty OrderedDict.
:rtype: collections.OrderedDict
"""
return OrderedDict()

@cached_property
def sigma_data(self):
"""
:return: Empty OrderedDict.
:rtype: collections.OrderedDict
"""
return OrderedDict()


class LogLikelihood(HessianObjective, BaseIndependentObjective):
"""
Error function to be minimized by a minimizer in order to *maximize*
the log-likelihood.
Expand Down Expand Up @@ -468,7 +506,7 @@ def eval_hessian(self, ordered_parameters=[], **parameters):
return np.atleast_2d(np.squeeze(np.array(result)))


class MinimizeModel(HessianObjective):
class MinimizeModel(HessianObjective, BaseIndependentObjective):
"""
Objective to use when the model itself is the quantity that should be
minimized. This is only supported for scalar models.
Expand Down
25 changes: 16 additions & 9 deletions tests/test_constrained.py
Expand Up @@ -17,7 +17,7 @@
SLSQP, MINPACK, TrustConstr, ScipyConstrainedMinimize, COBYLA
)
from symfit.core.support import key2str
from symfit.core.objectives import MinimizeModel, LogLikelihood
from symfit.core.objectives import MinimizeModel, LogLikelihood, LeastSquares
from symfit.core.fit import ModelError
from symfit import (
Symbol, MatrixSymbol, Inverse, CallableModel, sqrt, Sum, Idx, symbols
Expand Down Expand Up @@ -704,19 +704,23 @@ def test_data_for_constraint(self):

# Allowed
fit = Fit(model, x=xdata, y=ydata, Y=2, constraints=[constraint])
self.assertIsInstance(fit.objective, LeastSquares)
self.assertIsInstance(fit.minimizer.constraints[0], MinimizeModel)
fit = Fit(model, x=xdata, y=ydata)
self.assertIsInstance(fit.objective, LeastSquares)
fit = Fit(model, x=xdata, objective=LogLikelihood)
self.assertIsInstance(fit.objective, LogLikelihood)

# Not allowed
with self.assertRaises(TypeError):
fit = Fit(model, x=xdata, y=ydata, Y=2)
with self.assertRaises(TypeError):
fit = Fit(model, x=xdata, y=ydata, Y=2, Z=3, constraints=[constraint])
# TODO: Uncomment these next lines when #214 has been fixed.
# with self.assertRaises(TypeError):
# fit = Fit(model, x=xdata)
# with self.assertRaises(TypeError):
# fit = Fit(model, x=xdata, y=ydata, objective=LogLikelihood)
# Since #214 has been fixed, these properly raise an error.
with self.assertRaises(TypeError):
fit = Fit(model, x=xdata)
with self.assertRaises(TypeError):
fit = Fit(model, x=xdata, y=ydata, objective=LogLikelihood)


def test_constrained_dependent_on_model(self):
Expand Down Expand Up @@ -875,12 +879,14 @@ def test_constrained_dependent_on_matrixmodel(self):
self.assertEqual(fit.constraints[0].interdependent_vars, [B, y])
self.assertEqual(fit.constraints[0].params, [A, mu, sig])
self.assertEqual(fit.constraints[0].constraint_type, Eq)
self.assertIsInstance(fit.objective, LeastSquares)
self.assertIsInstance(fit.minimizer.constraints[0], MinimizeModel)

self.assertEqual(set(k for k, v in fit.data.items() if v is not None),
self.assertEqual({k for k, v in fit.data.items() if v is not None},
{x, y, dx, M, I, fit.model.sigmas[y]})
# These belong to internal variables
self.assertEqual(set(k for k, v in fit.data.items() if v is None),
{B, constraint.sigmas[Y], Y})
self.assertEqual({k for k, v in fit.data.items() if v is None},
{constraint.sigmas[Y], Y})

constr_result = fit.execute()
# The constraint should not be met for the unconstrained fit
Expand Down Expand Up @@ -993,6 +999,7 @@ def test_constrainedminimizers(self):
results = []
for minimizer in minimizers:
fit = Fit(- model, minimizer=minimizer)
self.assertIsInstance(fit.objective, MinimizeModel)
fit_result = fit.execute(tol=1e-15)
results.append(fit_result)

Expand Down
2 changes: 1 addition & 1 deletion tests/test_finite_difference.py
Expand Up @@ -197,7 +197,7 @@ def test_harmonic_oscillator_errors(self):
noise = np.random.normal(1, 0.05, size=t_data.shape)
x_data = ode_model(t=t_data, k=100).x * noise

ode_fit = sf.Fit(ode_model, t=t_data, x=x_data)
ode_fit = sf.Fit(ode_model, t=t_data, x=x_data, v=None)
ode_result = ode_fit.execute()

phi = 0
Expand Down

0 comments on commit 84dd102

Please sign in to comment.