imported>Michael Hardy at 17:11, 14 June 2025

2025-06-14T17:11:18Z

New page

{{Bayesian statistics}}
The '''nested sampling algorithm''' is a [[computation]]al approach to the [[Bayesian statistics]] problems of comparing models and generating samples from posterior distributions. It was developed in 2004 by [[physicist]] John Skilling.<ref name="skilling1">{{cite journal |last= Skilling |first= John |title= Nested Sampling |journal= AIP Conference Proceedings |pages= 395–405 |year= 2004 |volume= 735 |doi= 10.1063/1.1835238|bibcode= 2004AIPC..735..395S }}</ref>

==Background==
[[Bayes' theorem]] can be applied to a pair of competing models <math>M_1</math> and <math>M_2</math> for data <math>D</math>, one of which may be true (though which one is unknown) but which both cannot be true simultaneously. The posterior probability for <math>M_1</math> may be calculated as:

: <math>
\begin{align}
P(M_1\mid D) & = \frac{P(D\mid M_1) P(M_1)}{P(D)} \\
& = \frac{P(D\mid M_1) P(M_1)}{P(D\mid M_1) P(M_1) + P(D\mid M_2) P(M_2)} \\
& = \frac{1}{1 + \frac{P(D\mid M_2)}{P(D\mid M_1)} \frac{P(M_2)}{P(M_1)} }
\end{align}
</math>

The prior probabilities <math>M_1</math> and <math>M_2</math> are already known, as they are chosen by the researcher ahead of time. However, the remaining [[Bayes factor]] <math>P(D\mid M_2)/P(D\mid M_1)</math> is not so easy to evaluate, since in general it requires marginalizing nuisance parameters. Generally, <math>M_1</math> has a set of parameters that can be grouped together and called <math>\theta</math>, and <math>M_2</math> has its own vector of parameters that may be of different dimensionality, but is still termed <math>\theta</math>. The marginalization for <math>M_1</math> is

: <math>P(D\mid M_1) = \int d \theta \, P(D\mid \theta,M_1) P(\theta\mid M_1)</math>

and likewise for <math>M_2</math>. 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 <math>P(\theta\mid D,M_1)</math>.<ref name = "skilling2">{{cite journal |last= Skilling |first= John |title= Nested Sampling for General Bayesian Computation |journal= Bayesian Analysis |volume= 1 |issue= 4 |pages= 833–860 |year= 2006 |doi= 10.1214/06-BA127|doi-access= free }}</ref> It is an alternative to methods from the Bayesian literature<ref name="chen">{{cite book |author= Chen, Ming-Hui, Shao, Qi-Man, and Ibrahim, Joseph George |title= Monte Carlo methods in Bayesian computation |publisher= Springer |year= 2000 |isbn= 978-0-387-98935-8 |url= https://books.google.com/books?id=R3GeFfshc7wC}}</ref> 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 <math>Z=P(D\mid M)</math> where <math>M</math> is <math>M_1</math> or <math>M_2</math>:

Start with <math>N</math> points <math>\theta_1,\ldots,\theta_N</math> sampled from prior.
'''for''' <math>i=1</math> to <math>j</math> '''do''' % The number of iterations j is chosen by guesswork.
<math>L_i := \min(</math>current likelihood values of the points<math>)</math>;
<math>X_i := \exp(-i/N);</math>
<math>w_i := X_{i-1} - X_i</math>
<math>Z := Z + L_i\cdot w_i;</math>
Save the point with least likelihood as a sample point with weight <math>w_i</math>.
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 <math>L_i</math>.
'''end'''
'''return''' <math>Z</math>;

At each iteration, <math>X_i</math> is an estimate of the amount of prior mass covered by the hypervolume in parameter space of all points with likelihood greater than <math>\theta_i</math>. The weight factor <math>w_i</math> is an estimate of the amount of prior mass that lies between two nested hypersurfaces <math>\{ \theta \mid P(D\mid\theta,M) = P(D\mid\theta_{i-1},M) \}</math> and <math>\{ \theta \mid P(D\mid\theta,M) = P(D\mid\theta_i,M) \}</math>. The update step <math>Z := Z+L_i w_i</math> computes the sum over <math>i</math> of <math>L_i w_i</math> to numerically approximate the integral

: <math>
\begin{align}
P(D\mid M) &= \int P(D\mid \theta,M) P(\theta\mid M) \,d \theta \\
&= \int P(D\mid \theta,M) \,dP(\theta\mid M)
\end{align}
</math>

In the limit <math>j \to \infty</math>, this estimator has a positive bias of order <math> 1 / N</math><ref>{{cite journal |last= Walter |first= Clement| title=Point-process based Monte Carlo estimation| journal=Statistics and Computing|pages=219–236 |year=2017| volume=27|doi=10.1007/s11222-015-9617-y |arxiv=1412.6368|s2cid= 14639080}}</ref> which can be removed by using <math>(1 - 1/N)</math> instead of the <math>\exp (-1/N)</math> in the above algorithm.

The idea is to subdivide the range of <math>f(\theta) = P(D\mid\theta,M)</math> and estimate, for each interval <math>[f(\theta_{i-1}), f(\theta_i)]</math>, how likely it is a priori that a randomly chosen <math>\theta</math> would map to this interval. This can be thought of as a Bayesian's way to numerically implement [[Lebesgue integration]].<ref name="Jasa">{{cite journal |last1= Jasa|first1= Tomislav |last2= Xiang |first2= Ning|title= Nested sampling applied in Bayesian room-acoustics decay analysis |journal= Journal of the Acoustical Society of America |pages= 3251–3262 |year= 2012 |volume= 132 |issue= 5 |doi=10.1121/1.4754550|pmid= 23145609 |bibcode= 2012ASAJ..132.3251J|s2cid= 20876510 }}</ref>

== Choice of MCMC algorithm ==
The original procedure outlined by Skilling (given above in pseudocode) does not specify what specific Markov chain Monte Carlo algorithm should be used to choose new points with better likelihood.

Skilling's own code examples (such as one in Sivia and Skilling (2006),<ref>{{cite book | last=Sivia | first=Devinderjit | last2=Skilling | first2=John | title=Data Analysis: A Bayesian Tutorial | publisher=Oxford University Press, USA | publication-place=Oxford | date=June 2006 | isbn=978-0-19-856832-2 | page=}}</ref> [https://www.inference.org.uk/bayesys/sivia/lighthouse.c available on Skilling's website]) chooses a random existing point and selects a nearby point chosen by a random distance from the existing point; if the likelihood is better, then the point is accepted, else it is rejected and the process repeated. Mukherjee et al. (2006)<ref name="mukherjee">{{cite journal |last1= Mukherjee |first1= P. |last2= Parkinson |first2= D. |last3= Liddle |first3= A.R. |title= A Nested Sampling Algorithm for Cosmological Model Selection |journal= Astrophysical Journal |volume= 638 |issue= 2 |pages= 51–54 |year= 2006 |bibcode= 2006ApJ...638L..51M |doi= 10.1086/501068|arxiv= astro-ph/0508461|s2cid= 6208051 }}</ref> found higher acceptance rates by selecting points randomly within an ellipsoid drawn around the existing points; this idea was refined into the MultiNest algorithm<ref name="multinest">{{cite journal |last1= Feroz |first1= F. |last2= Hobson |first2= M.P. | last3=Bridges|first3=M.|title= MULTINEST: an efficient and robust Bayesian inference tool for cosmology and particle physics |journal= MNRAS |volume= 398 |issue= 4 |year= 2008 |url= https://arxiv.org/abs/0809.3437 |doi= 10.1111/j.1365-2966.2009.14548.x|arxiv= 0809.3437 }}</ref> which handles multimodal posteriors better by grouping points into likelihood contours and drawing an ellipsoid for each contour.

==Implementations==
Example implementations demonstrating the nested sampling algorithm are publicly available for download, written in several [[programming language]]s.
* Simple examples in [[C (programming language)|C]], [[R (programming language)|R]], or [[Python (programming language)|Python]] are on [http://www.inference.phy.cam.ac.uk/bayesys/ John Skilling's website].
* A [[Haskell (programming language)|Haskell]] port of [http://hackage.haskell.org/package/NestedSampling the above simple codes is on Hackage].
* An example in [[R (programming language)|R]] originally designed for [[Curve fitting|fitting]] [[Spectrum|spectra]] is described on [http://www.mrao.cam.ac.uk/~bn204/galevol/speca/rnested.html Bojan Nikolic's website] and is [https://github.com/bnikolic/ available on GitHub].
* A NestedSampler is part of the [[Python (programming language)|Python]] toolbox BayesicFitting<ref name="kester">{{cite journal |last1= Kester |first1= D. |last2= Mueller |first2= M. |title= BayesicFitting, a PYTHON toolbox for Bayesian fitting and evidence calculation.: Including a Nested Sampling implementation. |journal= Astronomy and Computing |volume= 37 |pages= 100503 |year= 2021 | doi= 10.1016/j.ascom.2021.100503|doi-access= free |arxiv= 2109.11976 |bibcode= 2021A&C....3700503K }}</ref> for generic model fitting and evidence calculation. It is [https://github.com/dokester/BayesicFitting available on GitHub].
* An implementation in [[C++]], named DIAMONDS, [https://github.com/JorisDeRidder/ is on GitHub].
* A highly modular [[Python (programming language)|Python]] parallel example for [[statistical physics]] and [[condensed matter physics]] uses is [https://github.com/js850/nested_sampling on GitHub].
* pymatnest is a package designed for exploring the [[energy landscape]] of different materials, calculating thermodynamic variables at arbitrary temperatures and locating [[phase transitions]] is [https://github.com/libAtoms/pymatnest on GitHub]
* The MultiNest software package is capable of performing nested sampling on multi-modal posterior distributions.<ref name="multinest" /><ref name="feroz">{{cite journal |last1= Feroz |first1= F. |last2= Hobson |first2= M.P. |title= Multimodal nested sampling: an efficient and robust alternative to Markov Chain Monte Carlo methods for astronomical data analyses |journal= MNRAS |volume= 384 |issue= 2 |pages= 449–463 |year= 2008 |url= http://adsabs.harvard.edu/cgi-bin/bib_query?arXiv:0704.3704 |doi= 10.1111/j.1365-2966.2007.12353.x |doi-access= free |bibcode=2008MNRAS.384..449F |arxiv= 0704.3704|s2cid= 14226032 }}</ref> It has interfaces for C++, [[Fortran]] and Python inputs, and is [https://github.com/farhanferoz/MultiNest available on GitHub].
* PolyChord is another nested sampling software package [https://github.com/PolyChord/PolyChordLite available on GitHub]. 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.<ref>{{cite journal |last1=Handley |first1=Will |first2=Mike |last2=Hobson |first3=Anthony |last3=Lasenby |title=polychord: next-generation nested sampling |journal=Monthly Notices of the Royal Astronomical Society |date=2015 |volume=453 |issue=4 |pages=4384–4398 |doi=10.1093/mnras/stv1911 |doi-access=free |bibcode=2015MNRAS.453.4384H |arxiv=1506.00171 |s2cid=118882763 }}</ref> It has interfaces to likelihood functions written in Python, Fortran, C, or C++.
* NestedSamplers.jl, a [[Julia (programming language)|Julia]] package for implementing single- and multi-ellipsoidal nested sampling algorithms is [https://github.com/TuringLang/NestedSamplers.jl on GitHub].
* [https://www.cse-lab.ethz.ch/korali/ Korali] is a high-performance framework for uncertainty quantification, optimization, and deep reinforcement learning, which also implements nested sampling.

==Applications==
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 [[cosmology|cosmological]] [[model selection]] and object detection, as it "uniquely combines accuracy, general applicability and computational feasibility."<ref name="mukherjee" /> A refinement of the algorithm to handle multimodal posteriors has been suggested as a means to detect astronomical objects in extant datasets.<ref name="feroz"/> 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]].<ref>{{cite journal |last1= Mthembu |first1= L. |last2= Marwala |first2= T. |last3= Friswell |first3= M.I. |last4= Adhikari |first4= S. |title= Model selection in finite element model updating using the Bayesian evidence statistic |journal= Mechanical Systems and Signal Processing |volume= 25 |issue= 7 |pages= 2399–2412 |year= 2011 |doi=10.1016/j.ymssp.2011.04.001|bibcode= 2011MSSP...25.2399M }}</ref> This sampling method has also been used in the field of materials modeling. It can be used to learn the [[Partition function (statistical mechanics)|partition function]] from [[statistical mechanics]] and derive [[thermodynamics|thermodynamic]] properties.<ref name="partay">{{cite journal|last1=Partay|first1=Livia B.|year=2010|title=Efficient Sampling of Atomic Configurational Spaces|journal=The Journal of Physical Chemistry B|volume=114|issue=32|pages=10502–10512|doi=10.1021/jp1012973|pmid=20701382|arxiv=0906.3544|s2cid=16834142}}</ref>

==Dynamic nested sampling==

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.<ref>{{cite journal |last1=Higson |first1=Edward |last2=Handley |first2=Will |last3=Hobson |first3=Michael |last4=Lasenby |first4=Anthony |title=Dynamic nested sampling: an improved algorithm for parameter estimation and evidence calculation |journal=Statistics and Computing |date=2019 |volume=29 |issue=5 |pages=891–913 |doi=10.1007/s11222-018-9844-0 |bibcode=2019S&C....29..891H |arxiv=1704.03459 |s2cid=53514669 }}</ref> 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:
* {{Proper name|dynesty}} – a Python implementation of dynamic nested sampling which can be [https://github.com/joshspeagle/dynesty downloaded from GitHub].<ref>{{cite journal |last=Speagle |first=Joshua |title=dynesty: A Dynamic Nested Sampling Package for Estimating Bayesian Posteriors and Evidences |journal=Monthly Notices of the Royal Astronomical Society |year=2020 |volume=493 |issue=3 |pages=3132–3158 |doi=10.1093/mnras/staa278 |doi-access=free |arxiv=1904.02180|s2cid=102354337 }}</ref>
* dyPolyChord: a software package which can be used with Python, C++ and Fortran likelihood and prior distributions.<ref>{{cite journal |last1=Higson |first1=Edward |title=dyPolyChord: dynamic nested sampling with PolyChord |journal=Journal of Open Source Software |date=2018 |volume=3 |issue=29 |page=965 |doi=10.21105/joss.00965 |doi-access=free }}</ref> dyPolyChord is [https://github.com/ejhigson/dyPolyChord available on GitHub].

Dynamic nested sampling has been applied to a variety of scientific problems, including analysis of gravitational waves,<ref>{{cite journal |last1=Ashton |first1=Gregory |title=Bilby: A User-friendly Bayesian Inference Library for Gravitational-wave Astronomy |journal=The Astrophysical Journal Supplement Series |date=2019 |volume=241 |issue=2 |page=13 |doi=10.3847/1538-4365/ab06fc |display-authors=etal|bibcode=2019ApJS..241...27A |arxiv=1811.02042 |s2cid=118677076 |doi-access=free }}</ref> mapping distances in space<ref>{{cite journal |last1=Zucker |first1=Catherine |title=Mapping Distances across the Perseus Molecular Cloud Using {CO} Observations, Stellar Photometry, and Gaia {DR}2 Parallax Measurements |journal=The Astrophysical Journal |date=2018 |volume=869 |issue=1 |page=83 |doi=10.3847/1538-4357/aae97c |display-authors=etal|arxiv=1803.08931 |s2cid=119446622 |doi-access=free }}</ref> and exoplanet detection.<ref>{{cite journal |last1=Günther |first1=Maximilian |title=A super-Earth and two sub-Neptunes transiting the nearby and quiet M dwarf TOI-270 |journal=Nature Astronomy |date=2019 |volume=3 |issue=12 |pages=1099–1108 |doi=10.1038/s41550-019-0845-5 |display-authors=etal|bibcode=2019NatAs...3.1099G |arxiv=1903.06107 |s2cid=119286334 }}</ref>

==See also==
*[[Bayesian model comparison]]
*[[List of algorithms]]

==References==
{{Reflist}}

[[Category:Bayesian statistics]]
[[Category:Model selection]]
[[Category:Randomized algorithms]]

Nested sampling algorithm - Revision history

imported>Michael Hardy at 17:11, 14 June 2025