Nested Sampling#

See Bayesian Model Selection - Nested Sampling for more information on nested sampling.

minimisers.NS.calcEllipsoid(u, VS)#

Calculate a bounding ellipsoid for a set of points u.

Parameters:
  • u (array) – Array of points to contain in the ellipsoid. Nxndims array where N is the number point and ndims is the number of dimensions.

  • VS (float) – The minimum volume that the bounding ellipsoid should have.

Returns:

  • B (array) – Bounding matrix for ellipsoid including scale factor for minimum volume.

  • mu (float) – The centroid of the ellipsoid.

  • VE (float) – The volume of the ellipsoid.

  • flag (bool) – Equals 1 if the ellipsoid cannot be drawn (too few points or bad condition number).

minimisers.NS.drawEllipsoidPoints(B, mu, N)#

Draw points uniformly from an ellipsoid.

Parameters:
  • B (array) – The bounding matrix for the ellipsoid.

  • mu (float) – The centroid of the ellipsoid.

  • N (int) – The number of points to draw.

Returns:

pnts – An array of N points drawn from the N-dimensional ellipsoid defined by bounding matrix B and centroid mu.

Return type:

array

minimisers.NS.drawMCMC(livepoints, cholmat, logLmin, prior, data, likelihood, model, nMCMC, parnames, extraparvals)#

Draw a sample from the prior volume using MCMC chains.

The new point will be found by evolving a random multi-dimensional sample from within the sample array, livepoints, using an MCMC with nMCMC iterations. The MCMC proposals will use a Student’s t-distribution (with 2 degrees of freedom) based on the Cholesky decomposed covariance matrix of the array, cholmat. 10% of the samples will actually be drawn using differential evolution by taking two random points from the current live points.

Parameters:
  • livepoints (array) – The live points in the current iteration of the nested sampler.

  • cholmat (array) – The Cholesky decomposed covariance matrix of livepoints.

  • logLmin (float) – The worst likelihood permitted in this iteration.

  • prior (array) – The prior information for the parameters.

  • data (array) – The problem struct and controls.

  • likelihood (function) – The likelihood function for the problem.

  • model (unknown) – Unused.

  • nMCMC (int) – The number of chains to use.

  • parnames (array) – The names of the parameters.

  • extraparvals (array) – A vector of additional parameters needed by the model.

Returns:

  • sample (array) – The new point sampled by the algorithm.

  • logL (float) – The log-likelihood of the new sample.

minimisers.NS.drawMultiNest(fracvol, Bs, mus, logLmin, prior, data, likelihood, model, parnames, extraparvals)#

Draw a sample from the prior volume using the MultiNest algorithm.

The new point will be found by drawing a random multi-dimensional sample from within the set of optimal ellipsoids constructed using the MultiNest algorithm.

Parameters:
  • fracvol (float) – The cumulative fractional volume of the ellipsoids.

  • Bs (array) – The bounding matrices of the ellipsoids to draw from.

  • mus (array) – The centroids of the ellipsoids to draw from.

  • logLmin (float) – The worst likelihood permitted in this iteration.

  • prior (array) – The prior information for the parameters.

  • data (array) – The problem struct and controls.

  • likelihood (function) – The likelihood function for the problem.

  • model (unknown) – Unused.

  • nMCMC (int) – The number of chains to use.

  • parnames (array) – The names of the parameters.

  • extraparvals (array) – A vector of additional parameters needed by the model.

Returns:

  • sample (array) – The new point sampled by the algorithm.

  • logL (float) – The log-likelihood of the new sample.

minimisers.NS.inEllipsoids(pnt, Bs, mus)#

Calculate how many of the ellipsoids contain the point pnt.

Parameters:
  • pnt (array) – The point to calculate ellipsoid membership for.

  • Bs (array) – The bounding matrices of the ellipsoids.

  • mus (array) – The centroids of the ellipsoids.

Returns:

N – The number of ellipsoids containing pnt.

Return type:

int

minimisers.NS.kmeans(X, K, maxerr)#

Perform K-means clustering for the data in X.

Finds K prototypes representing the samples in data matrix X, where each row of X represents a sample. Iterates until maximum norm difference between prototypes found in successive iterations is < maxerr

This script uses square Euclidean distance, but can be easily modified to use other metrics.

Parameters:
  • X (array) – An array of points, where each row of X is one point.

  • K (int) – The number of clusters to separate the points into.

  • maxerr (float) – The maximum norm difference between cluster prototypes.

Returns:

  • means (array) – matrix with each row a cluster prototype.

  • Nmeans (array) – Number of samples in each cluster.

  • membership (array) – Assigned class for each sample.

minimisers.NS.logPlus(logx, logy)#

Calculate log(x+y) given logx and logy.

It avoids problems of dynamic range when the exponentiated values of x or y are very large or small.

Parameters:
  • logx (float) – The two logarithmic values to add.

  • logy (float) – The two logarithmic values to add.

Returns:

logz – The value of log(x+y).

Return type:

float

minimisers.NS.mchol(G)#

Calculate the Cholesky (LDL) decomposition of a symmetric matrix G.

Includes a potential small perturbation so that we can calculate for non-positive-definite matrices.

Given a symmetric matrix G, find a matrix E of “small” norm, lower-triangular matrix L, and diagonal matrix D such that G+E is positive definite, and

\[G+E = LDL'\]

Also, calculate a direction pneg, such that if G is not PD, then

\[\text{pneg}' G \text{pneg} < 0\]

Note that if G is positive-definite, then the routine will return pneg=[].

Parameters:

G (array) – The matrix to calculate the Cholesky decomposition for.

Returns:

  • L (array) – The lower-triangular matrix such that G+E = LDL’.

  • D (array) – The diagonal matrix such that G+E = LDL’.

  • E (array) – The perturbation matrix such that G+E = LDL’.

  • pneg (array) – The descent direction when G is not positive-definite.

minimisers.NS.nest2pos(nest_samples, nLive)#

Convert nested samples to samples from the posterior distribution.

Parameters:
  • nest_samples (array) – The samples from the nested sampler.

  • nLive (int) – The number of live points.

Returns:

The samples from the posterior distribution for each of nest_samples.

Return type:

post_samples

minimisers.NS.nestedSampler(data, nLive, nMCMC, tolerance, flike, model, prior, parnames)#

Perform nested sampling on a given likelihood function and set of data with priors.

Parameters:
  • data (array) – The problem struct and controls.

  • nLive (int) – The number of live points to use.

  • nMCMC (int) – If 0, use MultiNest to draw points. if >0, use MCMC and differential evolution with nMCMC chains.

  • tolerance (float) – The tolerance of maximum log-likelihood changes between iterations.

  • flike (function) – The log-likelihood function to use.

  • model (function) – The function handle of a model function to pass to the likelihood function.

  • prior (array) – A cell array of parameter names, prior types, prior parameters, and boundary handling.

  • parnames (array) – The names of the parameters to optimise.

Returns:

  • logZ (float) – The final log-evidence calculated.

  • nest_samples (array) – The full set of points sampled during the run.

  • post_samples (array) – The posterior values of each point in nest_samples.

  • H (float) – The Shannon entropy of the evidence.

minimisers.NS.nsIntraFun(data, p)#

Calculate the log-likelihood of a data point during nested sampling.

Parameters:
  • data (array) – The problem struct and controls.

  • p (array) – The point in parameter space to calculate likelihood for.

minimisers.NS.optimalEllipsoids(u, VS)#

Optimally partition the samples u into a set of subclusters enclosed by ellipsoids.

Parameters:
  • u (array) – The points to enclose in subclusters.

  • VS (float) – The minimum volume of the set of ellipsoids.

Returns:

  • Bs (array) – An array of bounding matrices for the ellipsoids enclosing the subclusters, scaled to have at least the minimum volume required by the subclusters. ( (K x ndims) x ndims )

  • mus (array) – An array of centroids for the bounding ellipsoids. (K x ndims)

  • VEs (array) – An array of volumes for the bounding ellipsoids. (K x 1)

  • ns (array) – An array containing the number of points for each subcluster. (K x 1)

Notes

This is Algorithm 1 in the paper:

Feroz, F.; Hobson, M.P.; Bridges, M. (2008), “MULTINEST: an efficient and robust Bayesian inference tool for cosmology and particle physics”. DOI: 10.1111/j.1365-2966.2009.14548.x, URL: https://arxiv.org/abs/0809.3437

minimisers.NS.rescaleParameters(prior, params)#

Convert a uniform value from the unit hypercube to a value from the parameter priors.

Parameters:
  • prior (array) – The prior information for the parameters.

  • params (array) – The values of the sample from the unit hypercube.

Returns:

scaled – The values of params scaled to the priors.

Return type:

array

minimisers.NS.runNestedSampler(problemStruct, controls)#

Run the nested sampling algorithm for a given problem and controls.

Parameters:
  • problemStruct (struct) – the Project struct.

  • controls (struct) – the Controls struct.

Returns:

  • problemStruct (struct) – the output project struct.

  • result (struct) – the calculation and optimisation results object.

  • bayesResults (struct) – Additional Bayesian results from the algorithm.

minimisers.NS.splitEllipsoid(u, VS)#

Optimally bound a set of points by two ellipsoids.

It uses the k-means algorthim to split the points into two sub-clusters and uses an optimisation scheme to re-assign points, if necessary, between the sub-clusters. This is based on the description in Algorithm 1 of the MultiNest paper by Feroz, Hobson, and Bridges.

Parameters:
  • u (array) – The points to bound by two ellipsoids.

  • VS (float) – The minimum volume of the ellipsoids.

Returns:

  • u1, u2 (array) – The points in each ellipsoid.

  • VE1, VE2 (float) – The volume of each ellipsoid.

  • nosplit (bool) – Equals 1 if splitting into two ellipsoids is not possible.