diff --git a/fitbenchmarking/controllers/base_controller.py b/fitbenchmarking/controllers/base_controller.py index 84c084c08..179af6395 100644 --- a/fitbenchmarking/controllers/base_controller.py +++ b/fitbenchmarking/controllers/base_controller.py @@ -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) diff --git a/fitbenchmarking/controllers/bumps_controller.py b/fitbenchmarking/controllers/bumps_controller.py index 352731d35..0ff114f2d 100644 --- a/fitbenchmarking/controllers/bumps_controller.py +++ b/fitbenchmarking/controllers/bumps_controller.py @@ -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)) } diff --git a/fitbenchmarking/controllers/lmfit_controller.py b/fitbenchmarking/controllers/lmfit_controller.py index 0522fb2e9..a3d73c4b2 100644 --- a/fitbenchmarking/controllers/lmfit_controller.py +++ b/fitbenchmarking/controllers/lmfit_controller.py @@ -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())) diff --git a/fitbenchmarking/controllers/tests/test_controllers.py b/fitbenchmarking/controllers/tests/test_controllers.py index 81deeb94c..180c5f3db 100644 --- a/fitbenchmarking/controllers/tests/test_controllers.py +++ b/fitbenchmarking/controllers/tests/test_controllers.py @@ -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): diff --git a/fitbenchmarking/core/fitting_benchmarking.py b/fitbenchmarking/core/fitting_benchmarking.py index 6c19cc3b5..f12782a6e 100644 --- a/fitbenchmarking/core/fitting_benchmarking.py +++ b/fitbenchmarking/core/fitting_benchmarking.py @@ -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) diff --git a/fitbenchmarking/core/tests/test_fitting_benchmarking.py b/fitbenchmarking/core/tests/test_fitting_benchmarking.py index 99450dcb6..f1212f2c8 100644 --- a/fitbenchmarking/core/tests/test_fitting_benchmarking.py +++ b/fitbenchmarking/core/tests/test_fitting_benchmarking.py @@ -3,35 +3,35 @@ """ Tests for fitbenchmarking.core.fitting_benchmarking.Fit """ +import inspect +import json import os import unittest -from unittest.mock import patch -import json -import inspect import warnings from pathlib import Path +from unittest.mock import patch + import numpy as np from fitbenchmarking import test_files from fitbenchmarking.controllers.scipy_controller import ScipyController -from fitbenchmarking.jacobian.analytic_jacobian import Analytic from fitbenchmarking.core.fitting_benchmarking import Fit -from fitbenchmarking.cost_func.weighted_nlls_cost_func import ( - WeightedNLLSCostFunc) from fitbenchmarking.cost_func.nlls_cost_func import NLLSCostFunc +from fitbenchmarking.cost_func.weighted_nlls_cost_func import \ + WeightedNLLSCostFunc +from fitbenchmarking.jacobian.analytic_jacobian import Analytic from fitbenchmarking.parsing.parser_factory import parse_problem_file from fitbenchmarking.utils.checkpoint import Checkpoint -from fitbenchmarking.utils.options import Options -from fitbenchmarking.utils.fitbm_result import FittingResult from fitbenchmarking.utils.exceptions import (FitBenchmarkException, IncompatibleCostFunctionError, - UnsupportedMinimizerError, - UnknownMinimizerError, IncompatibleMinimizerError, + MaxRuntimeError, NoHessianError, NoJacobianError, - NoHessianError, - ValidationException, - MaxRuntimeError) + UnknownMinimizerError, + UnsupportedMinimizerError, + ValidationException) +from fitbenchmarking.utils.fitbm_result import FittingResult +from fitbenchmarking.utils.options import Options FITTING_DIR = "fitbenchmarking.core.fitting_benchmarking" TEST_FILES_DIR = os.path.dirname(inspect.getfile(test_files)) @@ -151,6 +151,7 @@ class PerformFitTests(unittest.TestCase): Verifies the output of the __perform_fit method in the Fit class when run with different options. """ + def setUp(self): """ Initializes the fit class for the tests @@ -167,20 +168,20 @@ def test_perform_fit_method(self): """ testcases = [{ - 'file': "ENSO.dat", - 'results': [111.70773805099354, - 107.53453144913736] - }, - { - 'file': "Gauss3.dat", - 'results': [76.64279628070524, - 76.65043476327958] - }, - { - 'file': "Lanczos1.dat", - 'results': [0.0009937705466940194, - 0.06269418241377904] - }] + 'file': "ENSO.dat", + 'results': [111.70773805099354, + 107.53453144913736] + }, + { + 'file': "Gauss3.dat", + 'results': [76.64279628070524, + 76.65043476327958] + }, + { + 'file': "Lanczos1.dat", + 'results': [0.0009937705466940194, + 0.06269418241377904] + }] for case in testcases: @@ -204,8 +205,7 @@ def test_perform_fit_method(self): assert emissions != np.inf @patch("fitbenchmarking.controllers." + - "base_controller.Controller.eval_confidence", - return_value=2) + "base_controller.Controller.eval_confidence") def test_eval_confidence_branches(self, mock): """ The test checks the eval_confidence branches in @@ -219,11 +219,12 @@ def test_eval_confidence_branches(self, mock): fit = Fit(options=self.options, data_dir='test', checkpointer=self.cp) + mock.return_value = 2 accuracy, _, _ = fit._Fit__perform_fit(controller) assert accuracy == 0.5 assert mock.call_count == 1 - controller.eval_confidence = 0 + mock.return_value = 0 accuracy, _, _ = fit._Fit__perform_fit(controller) assert accuracy == np.inf @@ -760,6 +761,7 @@ class StartingValueTests(unittest.TestCase): Verifies the output of the __loop_over_starting_values method in the Fit class when run with different options. """ + def setUp(self): """ Initializes the fit class for the tests diff --git a/fitbenchmarking/results_processing/plots.py b/fitbenchmarking/results_processing/plots.py index da108f4af..4491f0fd3 100644 --- a/fitbenchmarking/results_processing/plots.py +++ b/fitbenchmarking/results_processing/plots.py @@ -208,18 +208,24 @@ def plot_posteriors(self, result): fig = make_subplots(rows=len(par_names), cols=1, subplot_titles=par_names) - scipy_fit = result.params_pdfs['scipy_pfit'] - scipy_err = result.params_pdfs['scipy_perr'] for i, name in enumerate(par_names): fig.append_trace(go.Histogram(x=result.params_pdfs[name], histnorm='probability density'), row=i+1, col=1) - fig.add_vline(x=result.params_pdfs['scipy_pfit'][i], - row=i+1, col=1, line_color='red') - fig.add_vline(x=scipy_fit[i]-2*scipy_err[i], - row=i+1, col=1, line_color='red', line_dash='dash') - fig.add_vline(x=scipy_fit[i]+2*scipy_err[i], - row=i+1, col=1, line_color='red', line_dash='dash') + + if result.params_pdfs['scipy_pfit'] is not None: + scipy_fit = result.params_pdfs['scipy_pfit'] + scipy_err = result.params_pdfs['scipy_perr'] + + for i, name in enumerate(par_names): + fig.add_vline(x=scipy_fit[i], + row=i+1, col=1, line_color='red') + fig.add_vline(x=scipy_fit[i]-2*scipy_err[i], + row=i+1, col=1, line_color='red', + line_dash='dash') + fig.add_vline(x=scipy_fit[i]+2*scipy_err[i], + row=i+1, col=1, line_color='red', + line_dash='dash') fig.update_layout(showlegend=False) @@ -305,7 +311,7 @@ def plot_summary(cls, categories, title, options, figures_dir): plotlyfig.update_layout( title=title - ) + ) if result.plot_scale in ["loglog", "logx"]: plotlyfig.update_xaxes(type="log") diff --git a/fitbenchmarking/results_processing/tables.py b/fitbenchmarking/results_processing/tables.py index de56a0338..b8321e721 100644 --- a/fitbenchmarking/results_processing/tables.py +++ b/fitbenchmarking/results_processing/tables.py @@ -23,7 +23,8 @@ 4: "Solver doesn't support bounded problems", 5: "Solution doesn't respect parameter bounds", 6: "Solver has exceeded maximum allowed runtime", - 7: "Validation of the provided options failed"} + 7: "Validation of the provided options failed", + 8: "Confidence in fit could not be calculated"} SORTED_TABLE_NAMES = ["compare", "acc", "runtime", "local_min", 'emissions'] diff --git a/fitbenchmarking/results_processing/tests/test_plots.py b/fitbenchmarking/results_processing/tests/test_plots.py index 4b129d7b0..edf393d86 100644 --- a/fitbenchmarking/results_processing/tests/test_plots.py +++ b/fitbenchmarking/results_processing/tests/test_plots.py @@ -6,17 +6,18 @@ import os import unittest from tempfile import TemporaryDirectory +from unittest import mock import pandas as pd import plotly.graph_objects as go from fitbenchmarking import test_files +from fitbenchmarking.core.results_output import (create_directories, + create_index_page) from fitbenchmarking.results_processing import plots from fitbenchmarking.utils.checkpoint import Checkpoint from fitbenchmarking.utils.exceptions import PlottingError from fitbenchmarking.utils.options import Options -from fitbenchmarking.core.results_output import \ - create_directories, create_index_page def load_mock_result(): @@ -157,6 +158,27 @@ def test_plotly_fit_create_files(self): self.assertTrue(os.path.exists(path)) self.assertEqual(find_error_bar_count(path), 2) + def test_plot_posteriors_create_files(self): + """ + Test that plot_posteriors creates a file + """ + + self.plot.result.param_names = ['a', 'b'] + + result = mock.Mock() + result.params_pdfs = {"scipy_pfit": [1.0, 1.0], "scipy_perr": [ + 0.1, 0.2], "a": [1.5, 1.0, 1.2], "b": [0.9, 1.6, 1.1]} + result.sanitised_min_name.return_value = 'm10_[s1]_jj0' + result.sanitised_name = 'cf1_prob_1' + + file_name = self.plot.plot_posteriors(result) + + self.assertEqual(file_name, + 'm10_[s1]_jj0_posterior_' + 'pdf_plot_for_cf1_prob_1.html') + path = os.path.join(self.figures_dir, file_name) + self.assertTrue(os.path.exists(path)) + def test_multivariate_plot(self): """ Test that the plotting fails gracefully for multivariate problems.