Source-linked AI summary
MontePython 3: boosted MCMC sampler and other features
Thejs Brinckmann, Julien Lesgourgues
TL;DR
MontePython addresses the difficulty of efficient cosmological parameter inference when proposal covariances are poorly known and posteriors may be non-Gaussian. It introduces jumping-factor adaptation and inverse-Fisher proposal calculations, and evaluates them in difficult runs. The authors report reduced convergence times and simpler convergence procedures, while noting limitations for Fisher-matrix use with real data.
Problem
Efficient Metropolis-Hastings sampling is difficult when the initial covariance matrix is poorly matched to constraining data or when posteriors are non-Gaussian.
Method
MontePython 3 adds automatic jumping-factor adaptation, inverse-Fisher proposal densities, and supporting minimum-finding and iterative Fisher-matrix procedures.
Results
Using superupdate and an inverse Fisher matrix as input covariance matrix reduces convergence time in most cases and simplifies obtaining convergence.
Takeaways & Limitations
The authors recommend Fisher-matrix computation for forecasts lacking an accurate covariance matrix and generally recommend superupdate for Metropolis-Hastings runs.
Takeaways & Limitations
For real-data parameter extraction, the approximate Fisher matrix is only intermittently positive definite and invertible because the true best-fit point is not known accurately.
Abstract
from arXiv · showhide
MontePython is a parameter inference package for cosmology. We present the latest development of the code over the past couple of years. We explain, in particular, two new ingredients both contributing to improve the performance of Metropolis-Hastings sampling: an adaptation algorithm for the jumping factor, and a calculation of the inverse Fisher matrix, which can be used as a proposal density. We present several examples to show that these features speed up convergence and can save many hundreds of CPU-hours in the case of difficult runs, with a poor prior knowledge of the covariance matrix. We also summarise all the functionalities of MontePython in the current release, including new likelihoods and plotting options.
1. Introduction
MontePython 3 is a modular Python package for cosmological parameter inference, sampling, chain analysis, and plotting. This release adds performance features for Metropolis-Hastings and expands interfaces, likelihoods, and plotting functionality.
- Package scope: MontePython is a modular Python MCMC package for cosmological parameter inference.Its modular design facilitates code modification and implementation of modules for other Python sampling packages.
- Package scope: The `run` mode performs sampling, while the `info` mode analyzes MCMC chains and plots results.
- Interfaces: MontePython is interfaced with CLASS and its extensions, with potential extension to CAMB and PyCosmo.The text lists HiCLASS, SONG, ExoCLASS, CLASS SZ, and CLASSgal among available CLASS-related branches.
- Version 3 developments: The release introduces jumping-factor adaptation and inverse-Fisher proposal densities to improve Metropolis-Hastings performance.Examples report faster convergence and savings of many hundreds of CPU-hours for difficult runs with poor covariance-matrix knowledge.
- Version 3 developments: Version 3.0 includes new plotting options, additional likelihoods, and expanded cosmological parameter definitions, sampling, and analysis options.
2. Metropolis-Hastings sampling strategy
MontePython’s Metropolis-Hastings strategy uses Gaussian proposals, fast sampling, and adaptive updates to improve exploration efficiency. Version 3 adds superupdate, which adapts the jumping factor alongside covariance updates, especially for poorly specified or non-Gaussian problems.
- Available methods: MontePython supports Metropolis-Hastings, Nested Sampling, Cosmo Hammer, Fisher sampling, and post-processing methods.
- Proposal density: Metropolis-Hastings proposes Gaussian parameter jumps using a covariance matrix C and jumping parameter c.The covariance matrix encodes parameter correlations and relative standard deviations, while the proposal covariance is cC.
- Proposal density: For Gaussian posteriors, an approximately optimal proposal uses a covariance matrix close to the posterior covariance and jumping factor 2.4.The associated target acceptance rate is approximately 0.25.
- Fast sampling: Cholesky-based fast sampling separates parameter blocks so nuisance-parameter jumps can leave slower cosmological parameters unchanged.Parameters are ordered from slower to faster computational blocks, and oversampling factors govern repeated jumps within blocks.
- Adaptive updates: The `--update` option periodically updates C, while `--superupdate` additionally adapts the jumping factor j.The default covariance-update interval is 50 cycles, and the recommended superupdate interval is 20 cycles.
- Adaptive updates: Superupdate is intended for non-Gaussian posteriors and poor initial proposal densities, where manually selecting a suitable jumping factor is difficult.It is complementary to covariance updating and is normally used together with `--update`.
- Adaptive updates: The jumping factor is bounded below at 10% of its initial value to reduce the risk of chains becoming stuck in local minima.A lower desired jumping parameter should instead be supplied explicitly through `--f`.
3. Fisher matrix
MontePython v3 introduces a direct likelihood-based Fisher-matrix calculation whose inverse can initialize MCMC proposals and accelerate sampling. The method uses finite differences around a best-fit point, but its reliability for real-data extraction depends on locating that point accurately.
- Motivation and approach: The Fisher matrix is computed from numerical derivatives of the likelihood at a χ2 minimum, where its inverse gives the Gaussian covariance approximation.At the maximum-likelihood point, the Fisher matrix is positive definite and invertible.
- Applications: The inverse Fisher matrix can serve as the input covariance matrix for a Metropolis-Hastings run, providing a useful approximate starting proposal that the update mechanism can refine.High accuracy is unnecessary for this initialization because the MCMC update mechanism can improve an approximate covariance estimate.
- Applications: For mock-data forecasts, the fiducial model is the known maximum-likelihood point, so the Fisher method works particularly well.The forecast procedure uses fiducial spectra as observed spectra without generating a random realization.
- Limitations: For real-data extraction, Fisher calculations remain unreliable when the best-fit point is only approximate, and the implemented minimization algorithms are not robust with many nuisance parameters.The paper reports that the method may still work when the approximate and true best fits are sufficiently close, but this does not happen often.
- Motivation and approach: MontePython uses a direct likelihood-based Fisher evaluation rather than derivatives of observables such as Cℓ or P(k).The paper describes the two approaches as mathematically equivalent, while noting practical advantages of working directly with the likelihood.
- Computational cost: For a typical Planck run, the Fisher calculation requires 2049 likelihood evaluations, while reordered evaluation reduces Boltzmann-solver calls to Ncalls = 85.Without the special reordering, the corresponding count is Ncalls = 697.
4. Illustration of performance
The performance tests compare covariance updates, jumping-factor adaptation, and inverse-Fisher initialization across mock and current-data cosmological fits. Superupdate generally improves convergence, while Fisher initialization is especially valuable for difficult high-dimensional forecasts but can fail or degrade performance with non-Gaussian nuisance parameters.
- Mock forecasts: The tests compare update, superupdate, and superupdate + Fisher for mock datasets fitted with ΛCDM and extensions reaching 12 free parameters.The mock forecasts use fake Planck and DESI likelihoods and assess convergence with the Gelman-Rubin criterion.
- Mock forecasts: Poor initial covariance knowledge makes the proposal density too wide and requires learning new jump scales and correlations for added parameters.This difficulty is strongest when DESI is more constraining than the Planck-derived covariance and when extra parameters lack covariance information.
- Mock forecasts: For 6- and 8-parameter models, all methods reach R −1 = 0.03 or better, while superupdate and superupdate + Fisher achieve values up to 2 times smaller.For 10-, 11-, and 12-parameter models at 12 hours, only runs starting from a Fisher matrix obtain any convergence; at 48 hours, superupdate methods are 2.6 to 3.7 times better than update alone.
- Mock forecasts: In the representative 10-parameter run, Fisher starts near the right covariance and jumping factor, superupdate adapts through an initial low-jump phase, and update needs about 3500 additional steps.The update delay corresponds to about 12 hours on 48 cores, while superupdate saves about 600 core-hours for that run.
- Current data: The inverse Fisher matrix can be counterproductive for current data with many non-Gaussian nuisance parameters, degrading performance by a factor of three in the small-dataset run.In the extended-data case, the full Planck TTTEEE + BAO + LSS Fisher matrix was not invertible, so the analysis switched to Planck-lite TTTEEE.
- Conclusions: The authors conclude that superupdate is generally recommended, with benefits ranging from negligible to large depending on run difficulty.Update and superupdate are nearly equivalent when a new run already has excellent covariance and jumping-factor knowledge.
5. Summary and conclusions
MontePython v3 adds convergence-oriented tools and recommends Fisher preprocessing and superupdate for Metropolis-Hastings runs. The release also supports extensive parameter reparameterisations, while Fisher calculation remains non-robust in difficult non-Gaussian cases.
- Using superupdate and an inverse Fisher covariance matrix reduces convergence time in most cases and simplifies achieving convergence.The authors report that fewer trial-and-error runs are needed.
- The Fisher matrix computation is quick and is recommended before forecast runs when no accurate covariance matrix is available.
- Superupdate generally performs at least as well as update alone because it optimizes the acceptance rate.
- Fisher calculation is not entirely robust in v3.0.0 for difficult cases with many non-Gaussian parameters, where numerical errors or poor best-fit estimates can make the matrix non-invertible.
- MontePython supports reparameterising input variables before they are passed to CLASS, including transformations for cosmological and dark-energy-related parameters.The implementation can be extended through the data.py module.
Appendix A.1. Neutrino hierarchy
MontePython can sample total neutrino mass while converting each sampled value into individual masses under normal or inverted hierarchy assumptions. CLASS receives the individual masses, with standard neutrino settings fixed unless extra relics are studied.
- MontePython samples total neutrino mass while arranging individual masses under Normal or Inverted Hierarchy assumptions.Normal Hierarchy has two less massive and one more massive neutrino; Inverted Hierarchy reverses those counts.
- For each sampled total mass, the individual neutrino masses are calculated by solving a system of equations before being passed to CLASS.
- The mass-splitting inputs use current central values obtained from neutrino oscillation experiments.
- Standard active-neutrino settings fix N_ur=0.00641, N_ncdm=3, and T_ncdm='0.71611,0.71611,0.71611' unless extra relativistic relics are being studied.
Appendix A.2. Degenerate massive ν’s and varying Neff
MontePython supports degenerate massive-neutrino approximations and can vary Neff through the neutrino-species degeneracy. It also provides dark-energy reparameterisations, including sampling w0 + wa while restricting it to negative values.
- Degenerate massive neutrinos: Three degenerate massive neutrinos can be used as a computationally faster, though not fully realistic, approximation.
- Degenerate massive neutrinos: With M_tot, MontePython divides the total neutrino mass by three and passes the resulting particle mass to CLASS using one species with degeneracy three.
- Varying Neff: Neff can be varied by using deg_ncdm as a varying cosmological parameter instead of fixing it.
- Degenerate massive neutrinos: Using only one or two degenerate massive neutrinos with the rest massless is slightly inaccurate for the precision of current experiments.
- Dark-energy reparameterisation: For CPL dark energy, MontePython can sample w0 + wa and w0 to restrict w0 + wa to negative values.
Appendix A.4. Sterile ν parametrization
MontePython parameterises sterile neutrinos using an effective sterile mass and their contribution to Neff, then converts the effective mass into the physical mass passed to CLASS.
- The sterile-neutrino parameterisation samples effective mass m_s,eff and sterile-neutrino contribution deg_ncdm_2 to Neff.The effective mass avoids the region where sterile mass is arbitrarily large while its Neff contribution is arbitrarily small.
- The implementation uses two ncdm species at neutrino temperature and requires a fixed phantom varying parameter for normal-neutrino degeneracy.
- The effective sterile mass is converted to physical sterile mass by dividing by ΔN_s within data.py before being passed to CLASS.The conversion assumes the sterile species is ncdm species number 2 and Dodelson-Widrow-like.
Appendix B. Sampling options
MontePython v3 offers multiple sampling and post-processing methods, with Metropolis-Hastings options for adapting proposal distributions and jumping factors, plus Fisher-matrix computation.
- Sampling methods: MontePython supports Metropolis-Hastings, Nested Sampling, Cosmo Hammer, Fisher sampling, Importance Sampling, and derived-parameter reprocessing methods.These methods are selected through the code’s sampling-method options.
- Additional controls: MontePython also exposes controls for temperature sampling, minimization, parameter-space boundaries, and alternative samplers through MultiNest and CosmoHammer.The temperature option samples from P^1/T, while minimization can re-evaluate the starting point.
- Metropolis-Hastings options: The default Metropolis-Hastings method supports proposal-distribution updates and the superupdate approach, which also adapts the jumping factor.Superupdate is deactivated by default, with a recommended adaptation delay of 20 cycles.
- Adaptation controls: The adaptive option updates the covariance matrix and jumping factor during sampling, but is suitable only for single-chain runs.The jumping-factor adaptation start is controlled separately by --adaptive-ts.
- Fisher sampling: Fisher sampling computes a Fisher matrix that can be used for proposal construction, with options for asymmetric steps, step iterations, and likelihood-change tolerances.The Fisher-related controls include --fisher-asymmetric, --fisher-step-it, --fisher-delta, and --fisher-tol.
Appendix C. Analyze and plotting options
MontePython v3 significantly improves plotting options and presentation while remaining compatible with GetDist, whose advanced plotting tools users may prefer.
- Plotting overview: MontePython v3 significantly improves the range and general presentation of its plotting options.Users may instead use GetDist, and MontePython writes .paramnames files for compatibility.
- Plotting overview: MontePython’s analysis and plotting tools automatically evaluate burn-in, remove non-Markovian chain parts, and read parameter metadata from log.param.These conveniences complement compatibility with GetDist.
Appendix C.1. Chain analysis
MontePython’s chain-analysis and plotting options automate burn-in handling, expose configurable chain statistics, and support extensive customization of one- and two-dimensional plots, including Fisher contours.
- Chain analysis: MontePython defines burn-in as chain points preceding the first effective χ2 value below a configurable threshold.The threshold is based on χ2_min + 6 by default, with the value adjustable through LOG LKL CUTOFF.
- Chain analysis: Chain-analysis options can retain non-Markovian samples, keep a chosen post-burn-in fraction, compute covariance matrices, set histogram bins, temper posteriors, and suppress output.The defaults include retaining the full post-burn-in fraction and using 20 histogram bins.
- Chain analysis: Adaptive proposal updates make earlier chain steps non-Markovian, so MontePython normally analyzes only steps after adaptation stops.The --keep-non-markovian option can override this behavior, while --keep-fraction can discard an initial fraction of retained samples.
- Plot customization: Plotting controls cover one- and two-dimensional output, smoothing, interpolation, layout, typography, legends, tick formatting, line styles, and file formats.Options include posterior and Gaussian smoothing, contour selection, subplot counts, font and tick sizes, and PDF, EPS, or PNG output.
- Plot customization: Additional output controls can disable plots, restrict results to one-dimensional distributions, export separate files, and remove means from one-dimensional plots.These options include --no-plot, --no-plot-2d, --all, and --no-mean.
- Plot customization: Contour rendering can use line-only contours, while one-dimensional posterior smoothing supports no smoothing, cubic interpolation, or polynomial fits to ln(P).Two-dimensional interpolation can also refine histograms before Gaussian smoothing and contour extraction.
- Fisher plotting: When an inverse Fisher matrix is available, MontePython can overlay Fisher contours on MCMC contours and center them on best-fit parameters.Fisher overlays are enabled with --plot-fisher, while --center-fisher changes their center from input starting values.
Appendix C.3. More advanced plot customisation
MontePython supports advanced plot customization through .plot files, including parameter transformations, relabeling, rescaling, selection, executable plotting code, and legend controls.
- Advanced customization: Advanced customization is supplied through .plot files invoked with the --extra option.These files can override parameters defined by the plotting code and accept additional Python lines.
- Parameter transformations: Plot files can redefine chain parameters algebraically while the chains are read, enabling derived quantities such as Ωm to replace an existing column.The example replaces ωcdm with (0.01*ωb + ωcdm)/(H0/100.)^2.
- Parameter labels: Users can rename parameters for redefinitions or more readable LaTeX plot labels, while MontePython also performs basic automatic formatting.The info.to_change functionality supports further label customization.
- Plot selection: Plot files can change overall parameter rescaling factors and specify exactly which parameters appear in plots, including removal of nuisance parameters.These controls use info.new_scales and info.to_plot.
- New v3 functionality: New v3 functionality adds controls for one- and two-dimensional plot legends and legend names.The relevant settings are info.plot_legend_1d, info.plot_legend_2d, and info.legendnames.
- Plot-specific code: Users can add Python code for plot-specific lines, bands, arrows, labels, and other custom annotations before finalizing selected plots.The code is placed in small .py files and can target specific one- or two-dimensional plots with conditional statements.
Appendix D. Likelihoods
Appendix D catalogs the likelihoods included in MontePython v3.0.0, documents their citation references and update status, and distinguishes dependency requirements across current-data and forecast likelihoods.
- Likelihood catalog: MontePython v3.0.0 provides a comprehensive list of likelihoods, with references to the papers that should be cited when they are used.The references identify where each likelihood was published or first used with MontePython.
- Likelihood catalog: The package's likelihoods can be modified or newly created, with implementation guidance available in the MontePython online documentation.
- Dependency codes: Dependency codes distinguish self-contained likelihoods from those requiring external data, downloaded wrappers, or automatically generated mock data.SC denotes self-contained; D, W, and M identify the corresponding external-data, wrapper, and mock-data requirements.
- Table organization: Tables D.3 and D.4 list current-data likelihoods divided alphabetically into entries a–c and c–z.
- Table organization: Table D.5 lists forecast likelihoods, with Euclid and SKA entries subject to updates linked to publication status.