Nested sampling algorithm

The nested sampling algorithm is a computational approach to the Bayesian statistics problems of comparing models and generating samples from posterior distributions. It was developed in 2004 by physicist John Skilling.[1]

Background

edit

Bayes' theorem can be applied to a pair of competing models   and   for data  , one of which may be true (though which one is unknown) but which both cannot be true simultaneously. The posterior probability for   may be calculated as:

 

The prior probabilities   and   are already known, as they are chosen by the researcher ahead of time. However, the remaining Bayes factor   is not so easy to evaluate, since in general it requires marginalizing nuisance parameters. Generally,   has a set of parameters that can be grouped together and called  , and   has its own vector of parameters that may be of different dimensionality, but is still termed  . The marginalization for   is

 

and likewise for  . This integral is often analytically intractable, and in these cases it is necessary to employ a numerical algorithm to find an approximation. The nested sampling algorithm was developed by John Skilling specifically to approximate these marginalization integrals, and it has the added benefit of generating samples from the posterior distribution  .[2] It is an alternative to methods from the Bayesian literature[3] such as bridge sampling and defensive importance sampling.

Here is a simple version of the nested sampling algorithm, followed by a description of how it computes the marginal probability density   where   is   or  :

Start with   points   sampled from prior.
for   to   do        % The number of iterations j is chosen by guesswork.
     current likelihood values of the points ;
     
     
     
    Save the point with least likelihood as a sample point with weight  .
    Update the point with least likelihood with some Markov chain Monte Carlo steps according to the prior, accepting only steps that
    keep the likelihood above  .
end
return  ;

At each iteration,   is an estimate of the amount of prior mass covered by the hypervolume in parameter space of all points with likelihood greater than  . The weight factor   is an estimate of the amount of prior mass that lies between two nested hypersurfaces   and  . The update step   computes the sum over   of   to numerically approximate the integral

 

In the limit  , this estimator has a positive bias of order  [4] which can be removed by using   instead of the   in the above algorithm.

The idea is to subdivide the range of   and estimate, for each interval  , how likely it is a priori that a randomly chosen   would map to this interval. This can be thought of as a Bayesian's way to numerically implement Lebesgue integration.[5]

Implementations

edit

Example implementations demonstrating the nested sampling algorithm are publicly available for download, written in several programming languages.

  • Simple examples in C, R, or Python are on John Skilling's website.[6]
  • A Haskell port of the above simple codes is on Hackage.[7]
  • An example in R originally designed for fitting spectra is described at [8] and is on GitHub.[9]
  • A NestedSampler is part of the Python toolbox BayesicFitting [10] for generic model fitting and evidence calculation. It is on Github [11]
  • An example in C++, named Diamonds, is on GitHub.[12]
  • A highly modular Python parallel example for statistical physics and condensed matter physics uses is on GitHub.[13]
  • pymatnest is a Python package designed for exploring the energy landscape of different materials, calculating thermodynamic variables at arbitrary temperatures and locating phase transitions is on GitHub.[14]
  • The MultiNest software package is capable of performing nested sampling on multi-modal posterior distributions.[15] It has interfaces for C++, Fortran and Python inputs, and is available on GitHub.[16]
  • PolyChord is another nested sampling software package available on GitHub.[17] PolyChord's computational efficiency scales better with an increase in the number of parameters than MultiNest, meaning PolyChord can be more efficient for high dimensional problems.[18]
  • NestedSamplers.jl a Julia package for implementing single- and multi-ellipsoidal nested sampling algorithms is on GitHub.[19]
  • Korali is a high-performance framework for uncertainty quantification, optimization, and deep reinforcement learning, which also implements nested sampling.

Applications

edit

Since nested sampling was proposed in 2004, it has been used in many aspects of the field of astronomy. One paper suggested using nested sampling for cosmological model selection and object detection, as it "uniquely combines accuracy, general applicability and computational feasibility."[20] A refinement of the algorithm to handle multimodal posteriors has been suggested as a means to detect astronomical objects in extant datasets.[15] Other applications of nested sampling are in the field of finite element updating where the algorithm is used to choose an optimal finite element model, and this was applied to structural dynamics.[21] This sampling method has also been used in the field of materials modeling. It can be used to learn the partition function from statistical mechanics and derive thermodynamic properties. [22]

Dynamic nested sampling

edit

Dynamic nested sampling is a generalisation of the nested sampling algorithm in which the number of samples taken in different regions of the parameter space is dynamically adjusted to maximise calculation accuracy.[23] This can lead to large improvements in accuracy and computational efficiency when compared to the original nested sampling algorithm, in which the allocation of samples cannot be changed and often many samples are taken in regions which have little effect on calculation accuracy.

Publicly available dynamic nested sampling software packages include:

  • dynesty - a Python implementation of dynamic nested sampling which can be downloaded from GitHub.[24][25]
  • dyPolyChord: a software package which can be used with Python, C++ and Fortran likelihood and prior distributions.[26] dyPolyChord is available on GitHub.[27]

Dynamic nested sampling has been applied to a variety of scientific problems, including analysis of gravitational waves,[28] mapping distances in space[29] and exoplanet detection.[30]

See also

edit

References

edit
  1. ^ Skilling, John (2004). "Nested Sampling". AIP Conference Proceedings. 735: 395–405. Bibcode:2004AIPC..735..395S. doi:10.1063/1.1835238.
  2. ^ Skilling, John (2006). "Nested Sampling for General Bayesian Computation". Bayesian Analysis. 1 (4): 833–860. doi:10.1214/06-BA127.
  3. ^ Chen, Ming-Hui, Shao, Qi-Man, and Ibrahim, Joseph George (2000). Monte Carlo methods in Bayesian computation. Springer. ISBN 978-0-387-98935-8.{{cite book}}: CS1 maint: multiple names: authors list (link)
  4. ^ Walter, Clement (2017). "Point-process based Monte Carlo estimation". Statistics and Computing. 27: 219–236. arXiv:1412.6368. doi:10.1007/s11222-015-9617-y. S2CID 14639080.
  5. ^ Jasa, Tomislav; Xiang, Ning (2012). "Nested sampling applied in Bayesian room-acoustics decay analysis". Journal of the Acoustical Society of America. 132 (5): 3251–3262. Bibcode:2012ASAJ..132.3251J. doi:10.1121/1.4754550. PMID 23145609. S2CID 20876510.
  6. ^ John Skilling website
  7. ^ Nested sampling algorithm in Haskell at Hackage
  8. ^ Nested sampling algorithm in R on Bojan Nikolic website
  9. ^ Nested sampling algorithm in R on GitHub
  10. ^ Kester, D.; Mueller, M. (2021). "BayesicFitting, a PYTHON toolbox for Bayesian fitting and evidence calculation.: Including a Nested Sampling implementation". Astronomy and Computing. 37: 100503. arXiv:2109.11976. doi:10.1016/j.ascom.2021.100503.
  11. ^ Python toolbox containing a Nested sampling algorithm on GitHub
  12. ^ Nested sampling algorithm in C++ on GitHub
  13. ^ Nested sampling algorithm in Python on GitHub
  14. ^ Nested sampling algorithm for materials simulation on GitHub
  15. ^ a b Feroz, F.; Hobson, M.P. (2008). "Multimodal nested sampling: an efficient and robust alternative to Markov Chain Monte Carlo methods for astronomical data analyses". MNRAS. 384 (2): 449–463. arXiv:0704.3704. Bibcode:2008MNRAS.384..449F. doi:10.1111/j.1365-2966.2007.12353.x. S2CID 14226032.
  16. ^ The MultiNest nested sampling software package on GitHub
  17. ^ The PolyChord nested sampling software package on GitHub
  18. ^ Handley, Will; Hobson, Mike; Lasenby, Anthony (2015). "polychord: next-generation nested sampling". Monthly Notices of the Royal Astronomical Society. 453 (4): 4384–4398. arXiv:1506.00171. Bibcode:2015MNRAS.453.4384H. doi:10.1093/mnras/stv1911. S2CID 118882763.
  19. ^ Implementations of single and multi-ellipsoid nested sampling in Julia on GitHub
  20. ^ Mukherjee, P.; Parkinson, D.; Liddle, A.R. (2006). "A Nested Sampling Algorithm for Cosmological Model Selection". Astrophysical Journal. 638 (2): 51–54. arXiv:astro-ph/0508461. Bibcode:2006ApJ...638L..51M. doi:10.1086/501068. S2CID 6208051.
  21. ^ Mthembu, L.; Marwala, T.; Friswell, M.I.; Adhikari, S. (2011). "Model selection in finite element model updating using the Bayesian evidence statistic". Mechanical Systems and Signal Processing. 25 (7): 2399–2412. Bibcode:2011MSSP...25.2399M. doi:10.1016/j.ymssp.2011.04.001.
  22. ^ Partay, Livia B. (2010). "Efficient Sampling of Atomic Configurational Spaces". The Journal of Physical Chemistry B. 114 (32): 10502–10512. arXiv:0906.3544. doi:10.1021/jp1012973. PMID 20701382. S2CID 16834142.
  23. ^ Higson, Edward; Handley, Will; Hobson, Michael; Lasenby, Anthony (2019). "Dynamic nested sampling: an improved algorithm for parameter estimation and evidence calculation". Statistics and Computing. 29 (5): 891–913. arXiv:1704.03459. Bibcode:2019S&C....29..891H. doi:10.1007/s11222-018-9844-0. S2CID 53514669.
  24. ^ The dynesty nested sampling software package on GitHub
  25. ^ Speagle, Joshua (2020). "dynesty: A Dynamic Nested Sampling Package for Estimating Bayesian Posteriors and Evidences". Monthly Notices of the Royal Astronomical Society. 493 (3): 3132–3158. arXiv:1904.02180. doi:10.1093/mnras/staa278. S2CID 102354337.
  26. ^ Higson, Edward (2018). "dyPolyChord: dynamic nested sampling with PolyChord". Journal of Open Source Software. 3 (29): 965. doi:10.21105/joss.00965.
  27. ^ The dyPolyChord dynamic nested sampling software package on GitHub
  28. ^ Ashton, Gregory; et al. (2019). "Bilby: A User-friendly Bayesian Inference Library for Gravitational-wave Astronomy". The Astrophysical Journal Supplement Series. 241 (2): 13. arXiv:1811.02042. Bibcode:2019ApJS..241...27A. doi:10.3847/1538-4365/ab06fc. S2CID 118677076.
  29. ^ Zucker, Catherine; et al. (2018). "Mapping Distances across the Perseus Molecular Cloud Using {CO} Observations, Stellar Photometry, and Gaia {DR}2 Parallax Measurements". The Astrophysical Journal. 869 (1): 83. arXiv:1803.08931. doi:10.3847/1538-4357/aae97c. S2CID 119446622.
  30. ^ Günther, Maximilian; et al. (2019). "A super-Earth and two sub-Neptunes transiting the nearby and quiet M dwarf TOI-270". Nature Astronomy. 3 (12): 1099–1108. arXiv:1903.06107. Bibcode:2019NatAs...3.1099G. doi:10.1038/s41550-019-0845-5. S2CID 119286334.