From 46992b73750decbb548c4587bf5b6ce521672d7a Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Yannik=20Sch=C3=A4lte?= <31767307+yannikschaelte@users.noreply.github.com> Date: Mon, 2 Oct 2023 15:35:11 +0200 Subject: [PATCH 1/2] Documentation fixes (#1120) * document hierarchical optimization sources; add aesara+jax to objective docs * document origin of the mh + pt samplers * fix typo * fix typo * Update pypesto/sample/parallel_tempering.py Typo * streamline doc deps * update setup req versions * try cyipopt fix * test * tesT * test * fix test_visualize::test_optimization_scatter_with_x_None * test * fixup --------- Co-authored-by: Paul Jonas Jost <70631928+PaulJonasJost@users.noreply.github.com> --- .github/workflows/ci.yml | 4 +- .readthedocs.yml | 1 - doc/api.rst | 2 + pypesto/hierarchical/__init__.py | 27 ++++++++- pypesto/sample/adaptive_metropolis.py | 58 +++++++++++++++++-- pypesto/sample/adaptive_parallel_tempering.py | 19 +++++- pypesto/sample/metropolis.py | 25 +++++++- pypesto/sample/parallel_tempering.py | 21 ++++++- pypesto/visualize/parameters.py | 4 +- pyproject.toml | 4 +- setup.cfg | 6 +- tox.ini | 8 +-- 12 files changed, 157 insertions(+), 22 deletions(-) diff --git a/.github/workflows/ci.yml b/.github/workflows/ci.yml index ec4710d23..f3fb62908 100644 --- a/.github/workflows/ci.yml +++ b/.github/workflows/ci.yml @@ -188,8 +188,8 @@ jobs: runs-on: ubuntu-latest strategy: matrix: - # update to 3.11 when nlopt allows - python-version: ['3.9', '3.10'] + # ipopt does not work on 3.9 (https://github.com/mechmotum/cyipopt/issues/225) + python-version: ['3.11'] steps: - name: Check out repository diff --git a/.readthedocs.yml b/.readthedocs.yml index 4a2a4d3c1..b75a9c82b 100644 --- a/.readthedocs.yml +++ b/.readthedocs.yml @@ -17,7 +17,6 @@ python: path: . extra_requirements: - doc - - amici build: os: "ubuntu-20.04" diff --git a/doc/api.rst b/doc/api.rst index f23a6a0cb..470318cd0 100644 --- a/doc/api.rst +++ b/doc/api.rst @@ -12,6 +12,8 @@ API reference pypesto.history pypesto.logging pypesto.objective + pypesto.objective.aesara + pypesto.objective.jax pypesto.objective.julia pypesto.optimize pypesto.petab diff --git a/pypesto/hierarchical/__init__.py b/pypesto/hierarchical/__init__.py index fc326f25e..74caae7d8 100644 --- a/pypesto/hierarchical/__init__.py +++ b/pypesto/hierarchical/__init__.py @@ -2,7 +2,32 @@ Hierarchical ============ -Hierarchical optimization sub-package. +Contains an implementation of the hierarchical optimization approach, which +decomposes the parameter estimation problem into an outer and an inner problem. +In the outer problem, only dynamic parameters are optimized. +In the inner problem, conditional on the outer solution, static parameters are +optimized. +Static parameters can be parameters affecting directly the model observables, +such as scaling factors, offsets, and noise parameters. + +Hierarchical optimization has the advantage that the outer problem is typically +less complex than the full problem, and thus can be solved more efficiently. +Further, the inner problem can often be solved analytically, which is more +efficient. +Thus, hierarchical optimization can be used to speed up parameter estimation, +finding optimal values more efficiently and reliably. + +The implementation in this package is based on: + +* Loos et al. 2018 (https://doi.org/10.1093/bioinformatics/bty514), + who give an analytic solution for the inner problem for scaling factors and + noise standard deviations, for Gaussian and Laplace noise, using forward + sensitivity analysis (FSA). +* Schmiester et al. 2020 (https://doi.org/10.1093/bioinformatics/btz581), + who give an analytic solution for the inner problem for scaling factors, + offsets and noise standard deviations, for Gaussian and Laplace noise, + using adjoint sensitivity analysis (ASA). ASA allows to calculate gradients + substantially more efficiently in high dimension. """ from .calculator import HierarchicalAmiciCalculator diff --git a/pypesto/sample/adaptive_metropolis.py b/pypesto/sample/adaptive_metropolis.py index cbe9b1825..6b29bd804 100644 --- a/pypesto/sample/adaptive_metropolis.py +++ b/pypesto/sample/adaptive_metropolis.py @@ -8,7 +8,54 @@ class AdaptiveMetropolisSampler(MetropolisSampler): - """Metropolis-Hastings sampler with adaptive proposal covariance.""" + """Metropolis-Hastings sampler with adaptive proposal covariance. + + A core problem of the standard Metropolis-Hastings sampler is + its fixed proposal distribution, which must be manually tuned. + This sampler adapts the proposal distribution during the sampling + process based on previous samples. + It adapts the correlation structure and the scaling factor of the + proposal distribution. + For both parts, there exist a variety of methods, see + + * Ballnus et al. 2017. + Comprehensive benchmarking of Markov chain Monte Carlo methods for + dynamical systems + (https://doi.org/10.1186/s12918-017-0433-1) + * Andrieu et al. 2008. + A tutorial on adaptive MCMC + (https://doi.org/10.1007/s11222-008-9110-y) + + for a review. + + Here, we approximate the covariance matrix via a weighted average of + current and earlier samples, + with a decay factor determining the relative contribution of the + current sample and earlier ones to the weighted average of mean and + covariance. + The scaling factor we aim to converge to a fixed target acceptance rate + of 0.234, as suggested by theoretical results. + The implementation is based on: + + * Lacki et al. 2015. + State-dependent swap strategies and automatic reduction of number of + temperatures in adaptive parallel tempering algorithm + (https://doi.org/10.1007/s11222-015-9579-0) + * Miasojedow et al. 2013. + An adaptive parallel tempering algorithm + (https://doi.org/10.1080/10618600.2013.778779) + + In turn, these are based on adaptive MCMC as discussed in: + + * Haario et al. 2001. + An adaptive Metropolis algorithm + (https://doi.org/10.2307/3318737) + + For reference matlab implementations see: + + * https://github.com/ICB-DCM/PESTO/blob/master/private/performPT.m + * https://github.com/ICB-DCM/PESTO/blob/master/private/updateStatistics.m + """ def __init__(self, options: Dict = None): super().__init__(options) @@ -79,7 +126,7 @@ def _update_proposal( decay_constant=decay_constant, ) - # compute covariance scaling factor + # compute covariance scaling factor based on the target acceptance rate self._cov_scale *= np.exp( (np.exp(log_p_acc) - target_acceptance_rate) / np.power(n_sample_cur + 1, decay_constant) @@ -100,8 +147,11 @@ def update_history_statistics( n_cur_sample: int, decay_constant: float, ) -> Tuple[np.ndarray, np.ndarray]: - """ - Update sampling statistics. + """Update sampling mean and covariance matrix via weighted average. + + Update sampling mean and covariance matrix based on the previous + estimate and the most recent sample via a weighted average, + with gradually decaying update rate. Parameters ---------- diff --git a/pypesto/sample/adaptive_parallel_tempering.py b/pypesto/sample/adaptive_parallel_tempering.py index 788f752b4..2475db3e3 100644 --- a/pypesto/sample/adaptive_parallel_tempering.py +++ b/pypesto/sample/adaptive_parallel_tempering.py @@ -7,7 +7,24 @@ class AdaptiveParallelTemperingSampler(ParallelTemperingSampler): - """Parallel tempering sampler with adaptive temperature adaptation.""" + """Parallel tempering sampler with adaptive temperature adaptation. + + Compared to the base class, this sampler adapts the temperatures + during the sampling process. + This both simplifies the setup as it avoids manual tuning, + and improves the performance as the temperatures are adapted to the + current state of the chains. + + This implementation is based on: + + * Vousden et al. 2016. + Dynamic temperature selection for parallel tempering in Markov chain + Monte Carlo simulations + (https://doi.org/10.1093/mnras/stv2422), + + via a matlab reference implementation + (https://github.com/ICB-DCM/PESTO/blob/master/private/performPT.m). + """ @classmethod def default_options(cls) -> Dict: diff --git a/pypesto/sample/metropolis.py b/pypesto/sample/metropolis.py index 88a8ca41c..09495f534 100644 --- a/pypesto/sample/metropolis.py +++ b/pypesto/sample/metropolis.py @@ -11,7 +11,30 @@ class MetropolisSampler(InternalSampler): - """Simple Metropolis-Hastings sampler with fixed proposal variance.""" + """Simple Metropolis-Hastings sampler with fixed proposal variance. + + The Metropolis-Hastings sampler is a Markov chain Monte Carlo (MCMC) + method generating a sequence of samples from a probability + distribution. + + This class implements a simple Metropolis algorithm with fixed + symmetric Gaussian proposal distribution. + + For the underlying original publication, see: + + * Metropolis et al. 1953. + Equation of State Calculations by Fast Computing Machines + (https://doi.org/10.1063/1.1699114) + * Hastings 1970. + Monte Carlo sampling methods using Markov chains and their + applications + (https://doi.org/10.1093/biomet/57.1.97) + + For reference matlab implementations see: + + * https://en.wikipedia.org/wiki/Metropolis%E2%80%93Hastings_algorithm + * https://github.com/ICB-DCM/PESTO/blob/master/private/performPT.m + """ def __init__(self, options: Dict = None): super().__init__(options) diff --git a/pypesto/sample/parallel_tempering.py b/pypesto/sample/parallel_tempering.py index e82ab9e34..44a9f4453 100644 --- a/pypesto/sample/parallel_tempering.py +++ b/pypesto/sample/parallel_tempering.py @@ -10,7 +10,26 @@ class ParallelTemperingSampler(Sampler): - """Simple parallel tempering sampler.""" + """Simple parallel tempering sampler. + + Parallel tempering is a Markov chain Monte Carlo (MCMC) method that + uses multiple chains with different temperatures to sample from a + probability distribution. + The chains are coupled by swapping samples between them. + This allows to sample from distributions with multiple modes more + efficiently, as high-temperature chains can jump between modes, while + low-temperature chains can sample the modes more precisely. + + This implementation is based on: + + * Vousden et al. 2016. + Dynamic temperature selection for parallel tempering in Markov chain + Monte Carlo simulations + (https://doi.org/10.1093/mnras/stv2422), + + via a matlab-based reference implementation + (https://github.com/ICB-DCM/PESTO/blob/master/private/performPT.m). + """ def __init__( self, diff --git a/pypesto/visualize/parameters.py b/pypesto/visualize/parameters.py index 52509116f..d175c25b6 100644 --- a/pypesto/visualize/parameters.py +++ b/pypesto/visualize/parameters.py @@ -561,7 +561,7 @@ def optimization_scatter( parameter_indices = process_parameter_indices( parameter_indices=parameter_indices, result=result ) - # remove all start indices, that encounter an inf value at the start + # remove all start indices that encounter an inf value at the start # resulting in optimize_result[start]["x"] being None start_indices_finite = start_indices[ [ @@ -570,7 +570,7 @@ def optimization_scatter( ] ] # compare start_indices with start_indices_finite and log a warning - if not np.all(start_indices == start_indices_finite): + if len(start_indices) != len(start_indices_finite): logger.warning( 'Some start indices were removed due to inf values at the start.' ) diff --git a/pyproject.toml b/pyproject.toml index ed59d2789..06bda6a3f 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -4,8 +4,8 @@ [build-system] requires = [ - "setuptools >= 52.0.0", - "wheel >= 0.36.2", + "setuptools >= 68.2.0", + "wheel >= 0.41.2", ] build-backend = "setuptools.build_meta" diff --git a/setup.cfg b/setup.cfg index 6194ed8bd..573406d37 100644 --- a/setup.cfg +++ b/setup.cfg @@ -96,7 +96,7 @@ amici = petab = petab >= 0.2.0 ipopt = - ipopt + cyipopt >= 1.3.0 dlib = dlib >= 19.19.0 nlopt = @@ -143,6 +143,10 @@ doc = ipykernel >= 6.15.1 %(select)s %(fides)s + %(amici)s + %(petab)s + %(aesara)s + %(jax)s example = notebook >= 6.1.4 select = diff --git a/tox.ini b/tox.ini index 9866525c1..9477f90f0 100644 --- a/tox.ini +++ b/tox.ini @@ -95,10 +95,8 @@ description = Test Julia interface [testenv:optimize] -extras = test,dlib,pyswarm,cmaes,nlopt,fides,mpi,pyswarms,petab +extras = test,dlib,ipopt,pyswarm,cmaes,nlopt,fides,mpi,pyswarms,petab commands = - # workaround as ipopt has incomplete build - pip install git+https://github.com/mechmotum/cyipopt.git@master pytest --cov=pypesto --cov-report=xml --cov-append \ test/optimize -s description = @@ -127,8 +125,6 @@ description = allowlist_externals = bash extras = example,amici,petab,pyswarm,pymc3,cmaes,nlopt,fides commands = - # workaround as ipopt has incomplete build - pip install git+https://github.com/mechmotum/cyipopt.git@master bash test/run_notebook.sh 1 description = Run notebooks 1 @@ -171,7 +167,7 @@ description = [testenv:doc] extras = - doc,amici,petab + doc,amici,petab,aesara,jax commands = sphinx-build -W -b html doc/ doc/_build/html description = From 7458d04086eeefec5b6df8266c412eeb279a2921 Mon Sep 17 00:00:00 2001 From: Yannik Schaelte Date: Mon, 2 Oct 2023 15:40:11 +0200 Subject: [PATCH 2/2] update changelog+version --- CHANGELOG.rst | 36 ++++++++++++++++++++++++++++++++---- pypesto/version.py | 2 +- 2 files changed, 33 insertions(+), 5 deletions(-) diff --git a/CHANGELOG.rst b/CHANGELOG.rst index b51746437..516e6465c 100644 --- a/CHANGELOG.rst +++ b/CHANGELOG.rst @@ -6,8 +6,38 @@ Release notes .......... +0.3.2 (2023-10-02) +------------------ + +* Visualize: + * Restrict fval magnitude in waterfall with order_by_id (#1090) + * Hierarchical parameter plot fix (#1106) + * Fix y-limits on waterfall (#1109) +* Sampling: + * Use cloudpickle for pickling dynesty sampler (#1094) +* Optimize + * Small fix on hierarchical initialise (#1095) + * Fix startpoint sampling for hierarchical optimization #1105) + * SacessOptimizer: retry reading, delay deleting (#1110) + * SacessOptimizer: Fix logging with multiprocessing (#1112) + * SacessOptimizer: tmpdir option (#1115) +* Storage: + * fix storage (#1099) +* Examples + * Notebook on differences (#1098) +* Problem + * Add startpoint_method to Problem (#1093) +* General + * Added new entry to bib (#1100) + * PetabJL integration (#1089) + * Other platform tests (#1113) + * Dokumentation fixes (#1120) + * Updated CODEOWNER (#1122) + + 0.3.1 (2023-06-22) -------------------- +------------------ + * Visualize: * Parameter plot w/ hier. pars, noise estimation for splines (#1061) * Sampling: @@ -34,7 +64,7 @@ Release notes 0.3.0 (2023-05-02) -------------------- +------------------ New functionalities compared to 0.2.0: @@ -61,8 +91,6 @@ Not supported functionalities and versions compared to 0.2.0: * **Changed parameter indexing from boolean to int in profiling routines** - - 0.2 series .......... diff --git a/pypesto/version.py b/pypesto/version.py index 260c070a8..f9aa3e110 100644 --- a/pypesto/version.py +++ b/pypesto/version.py @@ -1 +1 @@ -__version__ = "0.3.1" +__version__ = "0.3.2"