Skip to content

Commit

Permalink
Fixed Bayes bugs (fitbenchmarking#1278)
Browse files Browse the repository at this point in the history
* fix bayes bugs

* update setting of param names

* small changes and tidy up

* fix failing test

* add more tests

* Remove print statement

Co-authored-by: Anthony <anthony.lim@stfc.ac.uk>

* e -> error_msg

Co-authored-by: Anthony <anthony.lim@stfc.ac.uk>

* review fixes

---------

Co-authored-by: Anthony <anthony.lim@stfc.ac.uk>
  • Loading branch information
jess-farmer and AnthonyLim23 authored Jun 5, 2024
1 parent 4e337ab commit 6f6cdf6
Showing 9 changed files with 226 additions and 103 deletions.
71 changes: 41 additions & 30 deletions fitbenchmarking/controllers/base_controller.py
Original file line number Diff line number Diff line change
@@ -6,6 +6,7 @@

import numpy
from scipy.optimize import curve_fit

from fitbenchmarking.utils.exceptions import (ControllerAttributeError,
IncompatibleHessianError,
IncompatibleJacobianError,
@@ -24,7 +25,7 @@ class Controller:

__metaclass__ = ABCMeta

VALID_FLAGS = [0, 1, 2, 3, 4, 5, 6, 7]
VALID_FLAGS = [0, 1, 2, 3, 4, 5, 6, 7, 8]

#: Within the controller class, you must
#: initialize a dictionary, ``algorithm_check``,
@@ -163,6 +164,7 @@ def flag(self):
| 5: `Solution doesn't respect parameter bounds`
| 6: `Solver has exceeded maximum allowed runtime`
| 7: `Validation of the provided options failed`
| 8: `Confidence in fit could not be calculated`
"""
return self._flag

@@ -248,36 +250,45 @@ def eval_confidence(self):
"""
Computes overall confidence in MCMC fit
"""
popt, pcov = curve_fit(self.problem.function,
xdata=self.data_x,
ydata=self.data_y,
p0=self.initial_params,
sigma=self.data_e,
maxfev=500)

perr = numpy.sqrt(numpy.diag(pcov))

self.params_pdfs['scipy_pfit'] = popt.tolist()
self.params_pdfs['scipy_perr'] = perr.tolist()

# calculate overall confidence within 2 sigma tolerance
par_conf = []
for i, name in enumerate(self.par_names):
tol = 2*perr[i]
hist, bin_edges = numpy.histogram(self.params_pdfs[name],
bins=100, density=True)
# check tol range is covered by hist range
tol_range = [popt[i]-tol, popt[i]+tol]
if tol_range[-1] < bin_edges[0] or tol_range[0] > bin_edges[-1]:
par_conf.append(0)
else:
width = numpy.diff(bin_edges)[0]
start_bin = numpy.argmin(abs(bin_edges-(popt[i]-tol)))
end_bin = numpy.argmin(abs(bin_edges-(popt[i]+tol)))
if start_bin == end_bin:
par_conf.append(hist[start_bin]*width)
self.params_pdfs['scipy_pfit'] = None
self.params_pdfs['scipy_perr'] = None
try:
popt, pcov = curve_fit(self.problem.function,
xdata=self.data_x,
ydata=self.data_y,
p0=self.initial_params,
sigma=self.data_e)

perr = numpy.sqrt(numpy.diag(pcov))

self.params_pdfs['scipy_pfit'] = popt.tolist()
self.params_pdfs['scipy_perr'] = perr.tolist()

# calculate overall confidence within 2 sigma tolerance
par_conf = []
for i, name in enumerate(self.par_names):
tol = 2*perr[i]
hist, bin_edges = numpy.histogram(
self.params_pdfs[name.replace('.', '_')],
bins=100, density=True
)
# check tol range is covered by hist range
tol_range = [popt[i]-tol, popt[i]+tol]
if tol_range[-1] < bin_edges[0] or \
tol_range[0] > bin_edges[-1]:
par_conf.append(0)
else:
par_conf.append(sum(hist[start_bin:end_bin]*width))
width = numpy.diff(bin_edges)[0]
start_bin = numpy.argmin(abs(bin_edges-(popt[i]-tol)))
end_bin = numpy.argmin(abs(bin_edges-(popt[i]+tol)))
if start_bin == end_bin:
par_conf.append(hist[start_bin]*width)
else:
par_conf.append(sum(hist[start_bin:end_bin]*width))
except RuntimeError as error_msg:
par_conf = 0
self.flag = 8
print("\n"+str(error_msg))

return numpy.prod(par_conf)

10 changes: 6 additions & 4 deletions fitbenchmarking/controllers/bumps_controller.py
Original file line number Diff line number Diff line change
@@ -181,12 +181,14 @@ def cleanup(self):
self.params_pdfs = {}
mcmc_draw = self._bumps_result.state.draw()
if self.fit_order != self._param_names:
for name in self._param_names:
ind = self.fit_order.index(name)
self.params_pdfs[name] = mcmc_draw.points[:, ind]
for name1, name2 in zip(self.problem.param_names,
self._param_names):
ind = self.fit_order.index(name2)
self.params_pdfs[name1] = mcmc_draw.points[:, ind].tolist()
else:
self.params_pdfs = {
self._param_names[i]: mcmc_draw.points[:, i].tolist()
self.problem.param_names[i]: mcmc_draw.points[:, i].tolist(
)
for i in range(len(self._param_names))
}

6 changes: 4 additions & 2 deletions fitbenchmarking/controllers/lmfit_controller.py
Original file line number Diff line number Diff line change
@@ -121,7 +121,7 @@ def lmfit_loglike(self, params):
"""
return self.cost_func.eval_loglike(
list(map(lambda name: params[name].value,
self.problem.param_names))
self._param_names))
)

def lmfit_jacobians(self, params):
@@ -184,7 +184,9 @@ def cleanup(self):
self.flag = 2

if self.minimizer == 'emcee':
self.params_pdfs = self.lmfit_out.flatchain.to_dict(orient='list')
params_pdf_dict = self.lmfit_out.flatchain.to_dict(orient='list')
self.params_pdfs = dict(zip(self.problem.param_names,
list(params_pdf_dict.values())))

self.final_params = list(map(lambda params: params.value,
self.lmfit_out.params.values()))
119 changes: 99 additions & 20 deletions fitbenchmarking/controllers/tests/test_controllers.py
Original file line number Diff line number Diff line change
@@ -6,21 +6,22 @@
import os
import platform
from unittest import TestCase
from unittest.mock import patch

import numpy as np
from pytest import test_type as TEST_TYPE # pylint: disable=no-name-in-module
from pytest import mark
from pytest import test_type as TEST_TYPE # pylint: disable=no-name-in-module

from conftest import run_for_test_types
from fitbenchmarking import test_files
from fitbenchmarking.controllers.base_controller import Controller
from fitbenchmarking.cost_func.weighted_nlls_cost_func import \
WeightedNLLSCostFunc
from fitbenchmarking.cost_func.loglike_nlls_cost_func import \
LoglikeNLLSCostFunc
from fitbenchmarking.cost_func.weighted_nlls_cost_func import \
WeightedNLLSCostFunc
from fitbenchmarking.hessian.scipy_hessian import Scipy as ScipyHessian
from fitbenchmarking.jacobian.default_jacobian import Default
from fitbenchmarking.jacobian.scipy_jacobian import Scipy
from fitbenchmarking.hessian.scipy_hessian import Scipy as ScipyHessian
from fitbenchmarking.parsing.parser_factory import parse_problem_file
from fitbenchmarking.utils import exceptions
from fitbenchmarking.utils.options import Options
@@ -30,41 +31,40 @@
from fitbenchmarking.controllers.controller_factory import \
ControllerFactory
from fitbenchmarking.controllers.dfo_controller import DFOController
from fitbenchmarking.controllers.lmfit_controller import LmfitController
from fitbenchmarking.controllers.minuit_controller import MinuitController
from fitbenchmarking.controllers.nlopt_controller import (NLoptController,
nlopt)
from fitbenchmarking.controllers.scipy_controller import ScipyController
from fitbenchmarking.controllers.scipy_go_controller import \
ScipyGOController
from fitbenchmarking.controllers.scipy_ls_controller import \
ScipyLSController
from fitbenchmarking.controllers.nlopt_controller import \
NLoptController, nlopt
from fitbenchmarking.controllers.lmfit_controller import LmfitController
if platform.system() != "Windows":
from fitbenchmarking.controllers.paramonte_controller import \
ParamonteController

if TEST_TYPE == 'all':
from fitbenchmarking.controllers.ceres_controller import CeresController
from fitbenchmarking.controllers.gofit_controller import GOFitController
from fitbenchmarking.controllers.gradient_free_controller import \
GradientFreeController
from fitbenchmarking.controllers.gsl_controller import GSLController
from fitbenchmarking.controllers.levmar_controller import LevmarController
from fitbenchmarking.controllers.mantid_controller import MantidController
from fitbenchmarking.controllers.ralfit_controller import RALFitController
from fitbenchmarking.controllers.gradient_free_controller import\
GradientFreeController
from fitbenchmarking.controllers.gofit_controller import GOFitController
from fitbenchmarking.controllers.ceres_controller import CeresController
from fitbenchmarking.controllers.theseus_controller import\
from fitbenchmarking.controllers.theseus_controller import \
TheseusController

if TEST_TYPE == 'matlab':
from fitbenchmarking.controllers.horace_controller import HoraceController
from fitbenchmarking.controllers.matlab_controller import MatlabController
from fitbenchmarking.controllers.matlab_opt_controller import\
from fitbenchmarking.controllers.matlab_curve_controller import \
MatlabCurveController
from fitbenchmarking.controllers.matlab_opt_controller import \
MatlabOptController
from fitbenchmarking.controllers.matlab_stats_controller import\
from fitbenchmarking.controllers.matlab_stats_controller import \
MatlabStatsController
from fitbenchmarking.controllers.matlab_curve_controller import\
MatlabCurveController
from fitbenchmarking.controllers.horace_controller import\
HoraceController


# pylint: disable=attribute-defined-outside-init, protected-access
@@ -237,6 +237,22 @@ def test_eval_conf(self):

self.assertAlmostEqual(controller.eval_confidence(), 0.192, 6)

@patch("fitbenchmarking.controllers.base_controller.curve_fit")
def test_eval_conf_failed_fit(self, mock):
"""
BaseSoftwareController: Test eval_confidence function handles
RuntimeError correctly
"""
controller = DummyController(self.cost_func)
controller.params_pdfs = {'A1': [4, 4, 4, 4, 4],
'A2': [3, 3.7, 3, 3, 3],
'A3': [2, 2, 2, 2.4, 2.5],
'A4': [0.5, 0.7, 1, 1, 1.2]}
mock.side_effect = RuntimeError
acc = controller.eval_confidence()
self.assertEqual(acc, 0)
self.assertEqual(controller.flag, 8)

def test_check_flag_attr(self):
"""
BaseSoftwareController: Test check_attributes function for _flag
@@ -1244,7 +1260,7 @@ def test_gradient_free(self):
self.shared_tests.check_diverged(controller)


@run_for_test_types(TEST_TYPE, 'default', 'all')
@run_for_test_types(TEST_TYPE, 'all')
@mark.skipif(
platform.system() == "Windows",
reason="Paramonte doesn't automatically detect MPI"
@@ -1256,7 +1272,8 @@ class BayesianControllerTests(TestCase):
"""

def setUp(self):
self.cost_func = make_cost_func(cost_func_type='loglike_nlls')
self.cost_func = make_cost_func(
'cubic-fba-test-go.txt', cost_func_type='loglike_nlls')
self.problem = self.cost_func.problem
self.shared_tests = ControllerSharedTesting()

@@ -1268,6 +1285,38 @@ def test_paramonte(self):
controller.minimizer = 'paraDram_sampler'
self.shared_tests.controller_run_test(controller)

assert len(controller.params_pdfs) == len(controller.final_params) + 1

def test_bumps(self):
"""
BumpsController: Test for output shape
"""
controller = BumpsController(self.cost_func)
controller.minimizer = 'dream'
self.shared_tests.controller_run_test(controller)

assert len(controller.params_pdfs) == len(controller.final_params)

def test_lmfit(self):
"""
LmfitController: Test for output shape
"""
controller = LmfitController(self.cost_func)
controller.minimizer = 'emcee'
self.shared_tests.controller_run_test(controller)

assert len(controller.params_pdfs) == len(controller.final_params)

def test_mantid(self):
"""
MantidController: Test for output shape
"""
controller = MantidController(self.cost_func)
controller.minimizer = 'FABADA'
self.shared_tests.controller_run_test(controller)

assert len(controller.params_pdfs) == len(controller.final_params)


@run_for_test_types(TEST_TYPE, 'all')
@mark.skipif(
@@ -1312,6 +1361,36 @@ def test_paramonte(self):

self.check_bounds(controller)

def test_bumps(self):
"""
ParamonteController: Test that parameter bounds are
respected for bounded problems
"""
controller = BumpsController(self.cost_func)
controller.minimizer = 'dream'

self.check_bounds(controller)

def test_lmfit(self):
"""
ParamonteController: Test that parameter bounds are
respected for bounded problems
"""
controller = LmfitController(self.cost_func)
controller.minimizer = 'emcee'

self.check_bounds(controller)

def test_mantid(self):
"""
ParamonteController: Test that parameter bounds are
respected for bounded problems
"""
controller = MantidController(self.cost_func)
controller.minimizer = 'FABADA'

self.check_bounds(controller)


@run_for_test_types(TEST_TYPE, 'default', 'all')
class FactoryTests(TestCase):
8 changes: 3 additions & 5 deletions fitbenchmarking/core/fitting_benchmarking.py
Original file line number Diff line number Diff line change
@@ -7,9 +7,9 @@
"""

import os
import platform
import timeit
import warnings
import platform
import time

import numpy as np
@@ -529,10 +529,8 @@ def __perform_fit(self, controller):
y=controller.data_y,
e=controller.data_e)
else:
if controller.eval_confidence != 0:
accuracy = 1/controller.eval_confidence()
else:
accuracy = np.inf
conf = controller.eval_confidence()
accuracy = 1 / conf if conf != 0 else np.inf

accuracy_check = any(np.isnan(n) for n in accuracy) \
if controller.problem.multifit else np.isnan(accuracy)
Loading

0 comments on commit 6f6cdf6

Please sign in to comment.