Improving the efficiency and robustness of nested sampling using posterior repartitioning

In real-world Bayesian inference applications, prior assumptions regarding the parameters of interest may be unrepresentative of their actual values for a given dataset. In particular, if the likelihood is concentrated far out in the wings of the assumed prior distribution, this can lead to extremely inefficient exploration of the resulting posterior by nested sampling algorithms, with unnecessarily high associated computational costs. Simple solutions such as broadening the prior range in such cases might not be appropriate or possible in real-world applications, for example when one wishes to assume a single standardised prior across the analysis of a large number of datasets for which the true values of the parameters of interest may vary. This work therefore introduces a posterior repartitioning (PR) method for nested sampling algorithms, which addresses the problem by redefining the likelihood and prior while keeping their product fixed, so that the posterior inferences and evidence estimates remain unchanged but the efficiency of the nested sampling process is significantly increased. Numerical results show that the PR method provides a simple yet powerful refinement for nested sampling algorithms to address the issue of unrepresentative priors.


Introduction
Bayesian inference (see e.g. MacKay 2003) provides a comprehensive framework for estimating unknown parameter(s) θ of some model with the assistance both of observed data D and prior knowledge of θ. One is interested in obtaining the posterior distribution of θ, and this can be expressed using Bayes' theorem as: where M represents model (or hypothesis) assumption(s), Pr(θ|D, M) ≡ P(θ) is the posterior probability density, Pr(D|θ, M) ≡ L(θ) is the likelihood, and Pr(θ|M) ≡ π(θ) is the prior of θ. Pr(D|M) ≡ Z is called the evidence (or marginal likelihood). We then have a simplified expression: and where Ψ represents the prior space of θ. The evidence Z is often used for model selection. It is the average of the likelihood over the prior, considering every possible choice of θ, and thus is not a function of the parameters θ. By ignoring the constant Z, the posterior P(θ) is proportional to the product of likelihood L(θ) and prior π(θ). The likelihood L(θ) is fully determined by the observation model (or measurement model / forward model) along with its corresponding noise assumptions. It is common that the structure of the observation model is predefined in real-world applications. By contrast, the prior distribution is often less well defined, and can be chosen in a number of ways, provided it is consistent with any physical requirements on the parameters θ (or quantities derived therefrom). One role of the prior distribution π(θ) is to localise the appropriate region of interest in the parameter space, which assists the inference process. One often chooses a standard distribution (such as Gaussian or uniform) as the prior when limited information is available a priori. In particular, the prior should be representative of the range of values that the parameters might take for the dataset(s) under analysis. An interesting discussion related to prior belief in a broader context can be found in Gelman (2008).
The approach outlined above works well in most scenarios, but it can be problematic if an inappropriate prior is chosen. In particular, if the true values of the parameters θ [or, more meaningfully, the location(s) of the peak(s) of the likelihood] lie very far out in the wings of the prior distribution π(θ), then this can result in very inefficient exploration of the parameter space by NS algorithms. In extreme cases, it can even result in a sampling algorithm failing to converge correctly, usually because of numerical inaccuracies, and incorrect posterior inferences (a toy example will be used to illustrate this problem in later sections).
This paper seeks to address the unrepresentative prior problem. One obvious solution is simply to augment the prior so that it covers a wider range of the parameter space. In some common cases, however, this might not be applicable. This is particularly true when one wishes to assume the same prior across a large number of datasets, for each of which the peak(s) of the likelihood may lie in very different regions of the parameter space. Moreover, in practical implementations, the specialists responsible for defining the prior knowledge, developing the measurement model, building the software, performing the data analysis, and testing the solution are often different people. Thus, there may be a significant overhead in communicating and understanding the full analysis pipeline before a new suitable prior could be agreed upon for a given scenario. This is a common occurrence in the analysis of, for example, production data in the oil and gas industry.
We therefore adopt an approach in this paper that circumvents the above difficulties. In particular, we present a posterior repartitioning (PR) method for addressing the unrepresentative prior problem in the context of NS algorithms (Skilling, 2006) for exploring the parameter space. One important way in which nested sampling differs from other methods is that it makes use of the likelihood L(θ) and prior π(θ) separately in its exploration of the parameter space, in that samples are drawn from the prior π(θ) such that they satisfy some likelihood constraint L(θ) > L * . By contrast, Markov chain Monte Carlo (MCMC) sampling methods or genetic algorithm variants are typically blind to this separation 1 , and deal solely in terms of the product L(θ)π(θ), which is proportional to the posterior P(θ). This difference provides an opportunity in the case of NS to 'repartition' the product L(θ)π(θ) by defining a new effective likelihoodL(θ) and priorπ(θ) (which is typically 'broader' than the original prior), subject to the conditionL(θ)π(θ) = L(θ)π(θ), so that the (unnormalised) posterior remains unchanged. Thus, in principle, the inferences obtained are unaffected by the use of the PR method, but, as we will demonstrate, the approach can yield significant improvements in sampling efficiency and also helps to avoid the convergence problems that can occur in extreme examples of unrepresentative priors. More generally, this approach highlights the intrinsic degeneracy between the 'effective' likelihood and prior in the formulation of Bayesian inference problems, which it may prove advantageous to exploit using NS methods more broadly than in merely addressing the unrepresentative prior problem, although we will defer such considerations to future publications. More discussion about generalised Bayesian prior design is given in Simpson et al (2017). This paper is organized as follows. Section 2 gives a brief summary of NS. Section 3 details the underlying problem, and illustrates it using a simple toy example. Section 4 describes the PR method and its implementation in the widely-used NS algorithm MultiNest. Section 5 shows some numerical results in simple synthetic examples. Section 6 concludes the proposed approach and discusses its advantages and limitations.
2 Nested sampling NS is a sequential sampling method that can efficiently explore the posterior distribution by repeatedly finding a higher likelihood region while keeping the number of samples the same. It consists of the following steps: -A certain number (N live ) of samples of the parameters θ are drawn from the prior distribution π(θ); these are termed 'live points'. -The likelihoods of these samples are computed through the likelihood function L(θ). -The sample with the lowest likelihood is removed and replaced by a sample again drawn from the prior, but constrained to a higher likelihood than that of the discarded sample.
-The above step is repeated until some convergence criteria are met (e.g. the difference in evidence estimates between two iterations falls below a predefined threshold); the final set of samples and the discarded samples are then used to estimate the evidence Z in model selection and obtain posteriorweighted samples for use in parameter estimation.
Pseudo code for the NS algorithm is given below. Note that it is only one of the various possible NS implementations. Other implementations share the same structure but may differ in details, for example in how X i or w i is calculated, or the method used for drawing new samples. See Skilling (2006) for details.

4
• Find the lowest likelihood in live sample and save it in L i .
• Replace the individual sample with likelihood L i by a newly drawn sample from restricted prior space Ψ i such that θ ∈ Ψ i satisfies L(θ) > L i .

8
• If max{L(θ n )}X i < exp(tol)Z, then stop. 9 end for 10 Increment Z by N live n=1 L(θ n )X I /N live . 11 Assign the sample replaced at iteration i the importance weight p i = L i w i /Z.
In Algorithm 1, X 0 represents the whole prior volume of prior space Ψ , and {X i } I i=1 are the constrained prior volumes at each iteration. The number of iterations I depends on a pre-defined convergence criterion tol on the accuracy of the final log-evidence value and on the complexity of the problem.
Among the various implementations of the NS algorithm, two widely used packages are MultiNest (Feroz et al, 2009(Feroz et al, , 2013 and PolyChord (Handley et al, 2015). MultiNest draws the new sample at each iteration using rejection sampling from within a multi-ellipsoid bound approximation to the iso-likelihood surface defined by the discarded point; the bound is constructed from the samples present at that iteration. PolyChord draws the new sample at each iteration using a number of successive slice-sampling steps taken in random directions. Please see Feroz et al (2009) andHandley et al (2015) for more details.

Unrepresentative prior problem
We describe a prior π(θ) as unrepresentative in the analysis of a particular dataset, if the true values of the parameters [or, more precisely, the peak(s) of the likelihood L(θ)] for that dataset lie very far into the wings of π(θ). In real-world applications, this can occur for a number of reasons, for example: (i) limited prior knowledge may be available, resulting in a simple tractable distribution being chosen as the prior, which could be unrepresentative; (ii) one may wish to adopt the same prior across a large number of datasets that might correspond to different true values of the parameters of interest, and for some of these datasets the prior may be unrepresentative. In any case, as we illustrate below in a simple example, an unrepresentative prior may result in very inefficient exploration of the parameter space, or failure of the sampling algorithm to converge correctly in extreme cases. This can be particularly damaging in applications where one wishes to perform analyses on many thousands (or even millions) of different datasets, since those (typically few) datasets for which the prior is unrepresentative can absorb a large fraction of the computational resources. Indeed, the authors have observed this phenomenon in practice in an industrial geophysical application consisting of only ∼ 1000 different datasets.
It is also worth mentioning that one could, of course, encounter the even more extreme case where the true parameter values, or likelihood peak(s), for some dataset(s) lie outside an assumed prior having compact support. This case, which one might describe as an unsuitable prior, is not addressed by our PR method, and is not considered here.

A univariate toy example
One may demonstrate the unrepresentative prior problem using a simple one-dimensional toy example. Suppose one makes N independent measurements (or observations) X = [x 1 , · · · , x n , · · · , x N ] of some quantity θ, such that where ξ denotes the simulated measurement noise, which is Gaussian distributed ξ ∼ N (µ ξ , σ 2 ξ ) with mean µ ξ and variance σ 2 ξ . For simplicity, we will assume the measurement process is unbiased, so that µ ξ = 0, and that the variance σ 2 ξ of the noise is known a priori (although it is a simple matter to relax these two assumptions).
The likelihood L(θ) is therefore simply the product of N Gaussian densities: For the purposes of illustration, we will assume the prior π(θ) also to be a Gaussian, with mean µ π = 0 and standard deviation σ π = 4, such that a priori one expects θ to lie in the range [−10, 10] with probability of approximately 0.99. Since the likelihood and prior are both Gaussian in θ, then so too is the posterior P(θ).
To illustrate the problem of an unrepresentative prior, we consider three cases in which the true value θ * of the unknown parameter is given, respectively, by: (1) θ * = 5, (2) θ * = 30 and (3) θ * = 40. Thus, case (1) corresponds to a straightforward situation in which the true value θ * lies comfortably within the prior, whereas cases (2) and (3) represent the more unusual eventuality in which the true value lies well into the wings of the prior distribution. In our simple synthetic example, one expects cases (2) and (3) to occur only extremely rarely. In real-world applications, however, the prior distribution is typically constructed on a case-bycase basis by analysts, and may not necessarily support a standard frequentist's interpretation of the probability of 'extreme' events. In fact, such situations are regularly encountered in real-world applications, when a large number of datasets are analysed. In each of the three cases considered, we set the variance of the simulated measurement noise to be σ ξ = 1 and the number of measurements is N = 20. Note that the width of the likelihood in (5) is proportional to 1/ √ N , so the unrepresentative prior problem becomes more acute as N increases.
Figures 1 (a), (c) and (e) show the prior, likelihood and posterior distributions for the cases (1), (2) and (3), respectively. One sees that as the true value θ * increases and lies further into the wings of the prior, the posterior lies progressively further to the left of the likelihood, as expected. As a result, in cases (2) and (3), the peak of the posterior (red dashed curve) is displaced to the left of the true value (black dashed line). This can be clearly observed in the zoomed-in plots within sub-figures (c) and (e). Figures 1 (b), (d) and (f) show histograms (blue bins) of the posterior samples obtained using MultiNest for cases (1), (2) and (3), respectively, together with the corresponding true analytical posterior distributions (red solid curves). In each case, the MultiNest sampling parameters were set to N live = 2000, efr = 0.8 and tol = 0.5 (see Feroz et al 2009 for details), and the algorithm was run to convergence. A natural estimatorθ and uncertainty ∆θ, respectively, for the value of the unknown parameter are provided by the mean and standard deviation of the posterior samples in each case, and are given in Table 1.
In case (1), one sees that the samples obtained are indeed consistent with being drawn from the true posterior, as expected. The meanθ and standard deviation ∆θ of the samples listed in Table 1 agree well with the mean µ P and standard deviation σ P of the true posterior distribution. In this case, MultiNest converged relatively quickly, requiring a total of 13529 likelihood evaluations. On repeating the entire analysis a total of 10 times, one obtains statistically consistent results in each case.
In case (2), one sees that the samples obtained are again consistent with being drawn from the true posterior. Indeed, from Table 1, one may verify that the mean and standard deviation of the samples agree well with those of the true posterior distribution. In this case, however, the convergence of MultiNest is much slower, requiring about 6 times the number of likelihood evaluations needed in case (1). This is a result of the true value lying far out in the wings of the prior distribution. Recall that NS begins by drawing N live samples from the prior and at each subsequent iteration replaces the sample having the lowest likelihood with a sample again drawn from the prior but constrained to have a higher likelihood. Thus, as the iterations progress, the collection of N live 'live points' gradually migrates from the prior to the peak of the likelihood. When the likelihood is concentrated very far out in the wings of the prior, this process can become very slow, even if one is able to draw each new sample from the constrained prior using standard methods (sometimes termed perfect nested sampling). In practice, this is usually not possible, so algorithms such as MultiNest and Poly-Chord use other methods that may require several likelihood evaluations before a new sample is accepted. Depending on the method used, an unrepresentative prior can also result in a significant drop in sampling efficiency, thereby increasing the required number of likelihood evaluations still further. On repeating the entire analysis a total of 10 times, once again obtains statistically consistent results in each case.
In case (3), one sees that the samples obtained are clearly inconsistent with being drawn from the true posterior. Indeed, the samples are concentrated at just a single value of θ. This behaviour may be understood by again considering the operation of NS. The algorithm begins by drawing N live = 2000 samples from the prior, which is a Gaussian with mean µ π = 0 and standard deviation σ π = 4. Thus, one would expect approximately only one such sample to lie outside the range [ −14, 14]. Moreover, since the likelihood is a Gaussian centred near the true value θ * = 40 with standard deviation ∼ 0.25, the live points will typically all lie in a region over which the likelihood is very small and flat (although, in this particular example, the values of the log-likelihood for the live points -which is the quantity used in the numerical calculations -are still distinguishable to machine precision). When the point with the lowest likelihood value is discarded, it must be replaced at the next NS iteration by another drawn from the prior, but with a larger likelihood. How this replacement sample is obtained de-pends on the particular NS implementation being used. As discussed in Section 2, MultiNest draws candidate replacement samples at each iteration using rejection sampling from within a multi-ellipsoid bound approximation to the iso-likelihood surface defined by the discarded point, which in just one dimension reduces simply to a range in θ. Since this bound is constructed from the samples present at that iteration, it will typically not extend far beyond the locations of the live points having the extreme values of the parameter θ. Thus, there is very limited opportunity to sample candidate replacement points from much larger values of θ, where the likelihood is significantly higher. Hence, as the NS iterations proceed, the migration of points from the prior towards the likelihood is extremely slow. Indeed, in this case, the migration is sufficiently slow that the algorithm terminates (in this case after 96512 likelihood evaluations) before reaching the main body of the likelihood and produces a set of posterior-weighted samples from the discarded points (see Feroz et al 2009 for details). Since this weighting is proportional to the likelihood, in this extreme case the recovered posterior is merely a 'spike' corresponding to the sample with the largest likelihood, as observed in Figure 1 (f). In short, the algorithm has catastrophically failed. On repeating the entire analysis a total of 10 times, one finds similar pathological behaviour in each case.
One may, of course, seek to improve the performance of NS in such cases in a number of ways. Firstly, one may adjust the convergence criterion (tol in Multi-Nest) so that many more NS iterations are performed, although there is no guarantee in any given problem that this will be sufficient to prevent premature convergence. Perhaps more useful is to ensure that there is a greater opportunity at each NS iteration of drawing candidate replacement points from larger values of θ, where the likelihood is larger. This may be achieved in a variety of ways. In MultiNest, for example, one may reduce the efr parameter so that the volume of the ellipsoidal bound (or the θ-range in this one-dimensional problem) becomes larger. Alternatively, as in other NS implementations, one may draw candidate replacement points using either MCMC sampling (Feroz and Hobson, 2008) or slice-sampling (Handley et al, 2015) and increase the number of steps taken before a candidate point is chosen.
All the of above approaches may mitigate the problem to some degree in particular cases (as we have verified in further numerical tests), but only at the cost of a simultaneous dramatic drop in sampling efficiency caused precisely by the changes made in obtaining candidate replacement points. Moreover, in more extreme cases these measures fail completely. In particular, if Table 1 MultiNest performance in the toy example illustrated in Figure 1.
Case (1) Case (2) Case (3)  the prior and the likelihood are extremely widely separated, the differences in the values of the log-likelihood of the live samples may fall below the machine accuracy used to perform the calculations. Thus, the original set of prior-distributed samples are likely to have loglikelihood values that are indistinguishable to machine precision. Thus, the 'lowest likelihood' sample to be discarded will be chosen effectively at random. Moreover, in seeking a replacement sample that is drawn from the prior but having a larger likelihood, the algorithm is very unlikely to obtain a sample for which the likelihood value is genuinely larger to machine precision.
Even if such a sample is obtained, then the above problems will re-occur in the next iteration when seeking to replace the next discarded sample, and so on. In this scenario, the sampling efficiency again drops dramatically, but more importantly the algorithm essentially becomes stuck and will catastrophically fail because of accumulated numerical inaccuracies.

Simple 'solutions'
A number of potential simple 'solutions' to the unrepresentative prior problem are immediately apparent. For example, one might consider the following: modify the prior distribution across one's analysis, either by increasing its standard deviation σ π , or even by adopting a different functional form, so that it should comfortably encompass the likelihood for all datasets; perform the analysis using the original prior for all the datasets, identify the datasets for which it is unrepresentative by monitoring the sampling efficiency and examining the final set of posterior samples for pathologies, and then modify the prior as above for these datasets.
Unfortunately, neither of these approaches is appropriate or realistic. The former approach is inapplicable since the prior may be representative for the vast majority of the datasets under analysis, and one should use this information in deriving inferences. Also, the former solution sacrifices the overall speed and computational efficiency, as the augmented prior is applied to all cases but not only the problematic ones. Choosing a proper trade-off between the efficiency and the coverage of prior is difficult when a large number of experiments need to be examined.
The latter solution requires one to identify various outlier cases (as the outlier cases could be very different from one to another), and also perform re-runs of those identified. It becomes a non-trivial computational problem when a single algorithm run requires a considerable amount of run time, or when the results of the outlier cases are needed for the next step computation, i.e. the whole process waits for the outlier cases to proceed. This could be trivial for some applications and could be very difficult for others in which many different outlier cases exist.

Posterior repartitioning method
The posterior repartitioning (PR) method addresses the unrepresentative prior problem in the context of NS algorithms (Skilling, 2006) for exploring the parameter space, without sacrificing computational speed or changing the inference obtained.

General expressions
In general, the 'repartition' of the product L(θ)π(θ) can be expressed as: whereL(θ) andπ(θ) are the new effectivelikelihood and prior, respectively. As a result, the (unnormalised) posterior remains unchanged. The modified priorπ(θ) can be any tractable distribution, which we assume to be appropriately normalised to unit volume. The possibility of repartitioning the posterior in NS was first mentioned in Feroz et al (2009), but equation (6) can also be viewed as the vanilla case (when the importance weight function equals to 1) of nested importance sampling proposed in Chopin and Robert (2010). One general advantage of NS is that the evidence (or marginal likelihood), which is intractable in most cases, can be accurately approximated. This is achieved by first defining V (l ) as the prior volume within the iso-likelihood contour L(θ) = l , namely where l is a real number that gradually rises from zero to the maximum of L(θ) as the NS iterations progress, so that V (l ) monotonically decreases from unity to zero. After PR, π(θ) is replaced byπ(θ), and the evidence can be calculated as It is worth noting, however, that in the case whereπ(θ) is not properly normalised, the 'modified evidence' Z obtained after PR is simply related to the original evidence by Provided one can evaluate the volume of the modified priorπ(θ), one may therefore straightforwardly recover the original evidence, if required. For many simple choices ofπ(θ), this is possible analytically, but may require numerical integration in general. It should be noted, however, that the normalistion of the modified prior is irrelevant for obtaining posterior samples. We now discuss some particular special choices forπ(θ).

Power posterior repartitioning
Rather than introducing a completely new prior distribution into the problem, a sensible choice is often simply to takeπ(θ) to be the original prior π(θ) raised to some power, and then renormalised to unit volume, such that where β ∈ [0, 1] and Z π (β) ≡ π(θ) β dθ. By altering the value of β, the modified prior can be chosen from a range between the original prior (β = 1) and the uniform distribution (β = 0). As long as the equality in equation (6) holds, the PR method can be applied separately for multiple unknown parameters with different forms of prior distributions. Figure 2 illustrates how the prior changes for different values of β in a one-dimensional problem. As the parameter β decreases from 1 to 0, the prior distribution evolves from a Gaussian centred on zero with standard deviation σ π = 4 to a uniform distribution, where the normalisation depends on the assumed support [−50, 50] of the unknown parameter θ. Indeed, the uniform modified priorπ(θ) ∼ U(a, b) is a special case, but often a useful choice. One advantage of this choice is that the range [a, b] can be easily set such that it accommodates the range of θ values required to overcome the unrepresentative prior problem, and the modified prior is trivially normalised. It can cause the sampling to be inefficient, however, since it essentially maximally broadens the search space (within the desired range).
The above approach is easily extended to multivariate problems with parameter vector θ = (θ 1 , θ 2 , · · · , θ N ) T . It is worth noting in particular the case where the original prior is a multivariate Gaussian, such that π(θ) = N (µ, Σ), where µ is the vector of means for each variable and Σ is the covariance matrix. The power modified priorπ(θ) is then given simply by N (µ, β −1 Σ) over the assumed supported region R of the parameter space, and There is unfortunately no robust universal guideline for choosing an appropriate value for β, since this depends on the dimensionality and complexity of the posterior and on the initial prior distribution assumed. Nonetheless, as demonstrated in the numerical examples presented in Section 5, there is a straightforward approach for employing the PR method in more realistic problems, in which the true posterior is not known. Namely, starting from β = 1 (which corresponds to the original prior), one can obtain inferences for progressively smaller values of β, according to some pre-defined or dynamic 'annealing schedule', until the results converge to a statistically consistent solution. The precise nature of the annealing schedule is unimportant, although either linearly or exponentially decreasing values of β seem the most natural approaches.

More general posterior repartitioning
Raising the original prior to some power β merely provides a convenient way of defining the modified prior, since it essentially just broadens the original prior by some specified amount. In general, however,π(θ) can be any tractable distribution. For example, there is no requirement for the modified prior to be centred at the same parameter value as the original prior. One could, therefore, choose a modified prior that broadens and/or shifts the original one, or a modified prior that has a different form from the original. Note that, in this generalised setting, the modified prior should at least be non-zero everywhere that the original prior is non-zero.

Diagnostics of the unrepresentative prior problem
This paper focuses primarily on how to mitigate the unrepresentative prior problem using PR. Another critical issue, however, is how one may determine when the prior is unrepresentative in the course the analysis of some (large number of) dataset(s). We comment briefly on this issue here.
Diagnosing the unrepresentative prior problem beforehand is generally difficult. Thus, designing a practical engineering-oriented solution is helpful in addressing most such problems. The goal of this diagnostic is to identify abnormal cases amongst a number of datasets during the analysis procedure. We assume that at least a few 'reliable' (sometimes called 'gold standard') datasets, which do not suffer from the unrepresentative prior problem, have been analysed before the diagnostics. The reliability threshold of a dataset varies depending on different scenarios, but (ideally) a gold standard dataset should: (1) be recognised as such by field experts; (2) have all of its noise sources clearly identified and characterised; (3) yield parameter estimates that are consistent with true values either known a priori or determined by other means. These provide us with some rough but reliable information and prior knowledge, such as runtime, convergence rate, and the shape of posterior distribution. We denote this information as the available knowledge for the problem of interest.
One may then employ a diagnostic scheme of the type illustrated in Figure 3, which is composed of two parts: on-the-fly diagnostics and after-run diagnostics. On-the-fly diagnostics involve monitoring the runtime and convergence status during the analysis of each dataset. Specifically, runtime monitoring involves simply checking whether the runtime of an individual analysis is greatly different from those of the available knowledge. Similarly, convergence rate checks compare the speed Fig. 3 A flow chart of a designed diagnostic process. The two main steps of the diagnostic process are highlighted in dark blue. The process starts by running a sampling algorithm for Bayesian parameter estimation (the top small block), and proceeds with two hierarchical diagnostics steps to evaluate the trail of interest. 'Available knowledge' is defined as reliable experimental information and prior knowledge that one could obtain in advance.
of convergence between the current run and the available knowledge. If both results are consistent with those in the available knowledge, the diagnostic process proceeds to after-run diagnostics. Note that the quantitative consistency check can be defined in various ways. A simple method is to set a threshold for the difference between available knowledge and individual runs. For instance, the result from an individual run can be considered as a reliable one if the error between the individual run result and the mean of the available knowledge is within a certain threshold. Such criteria should be carefully discussed by field experts on a case-by-case basis.
After-run diagnostics compare the computed posterior with the available knowledge. One plausible afterrun diagnostic is to evaluate some 'distance' measure between the assumed prior and the posterior distribu-tion resulting from the analysis. An obvious choice is to employ the Kullback-Leibler (KL) divergence (see, e.g., Bishop 2006). The KL divergence quantifies the difference between two probability distributions by calculating their relative entropy. A larger KL divergence indicates a greater difference between the two distributions. The KL divergence is, however, an asymmetric measure and its value is not bounded. To overcome these drawbacks, one could also consider the Jensen-Shannon divergence (Endres and Schindelin, 2003), which is a symmetric variant of the KL divergence. The posterior may also be compared with the available knowledge in the outlier check step.
Finally, we note that a diagnostic analysis is valid when it is performed using the same algorithm specifications. For instance, N live , efr, and tol settings should be the same in MultiNest when performing diagnostic analysis. In any case, once a reasonable diagnostic metric is constructed, the abnormal trials can be identified according to some predetermined criteria and examined, and the proposed PR scheme can be applied on a case-by-case basis. A simple illustration of this process is presented in the bivariate example case in the next section.

Numerical examples
We begin by illustrating the PR method in two numerical examples, one univariate and the other a bivariate Gaussian posterior. Our investigation is then extended to higher dimensional (from 3 to 15 dimensions) Gaussian posteriors, to explore its stability to the 'curse of dimensionality'. Finally, we consider a bivariate non-Gaussian example. In particular, we compare the performance of the MultiNest sampler before and after applying PR.
We use the open-source MultiNest package (Feroz et al, 2009) and set efficiency parameter efr = 0.8, convergence tolerance parameter tol = 0.5, multi-modal parameter mmode = False, random seed control parameter seed = −1, and the constant efficiency mode ceff = False for all the following examples. The number of live samples N live varies in different cases. We keep the other MultiNest tuning options in their default values. See (Feroz et al, 2009) and its corresponding MultiNest Fortran package for details of these default settings.
In some of the multi-dimensional cases, we also compare the MultiNest performance with MCMC. Specifically, a standard Metropolis-Hastings sampler is implemented and applied to the same numerical examples. Other MCMC samplers such as No-U-Turn Sampler (NUTS), and slice samplers give similar performance in the numerical examples. One popular Python implementation of these samplers can be found in PyMC3 (Salvatier et al, 2016) package. In some cases, we also compare the performance of importance sampling (Neal, 2001;Tokdar and Kass, 2010;Martino et al, 2018), using a standard IS implementation from Python package 'pypmc' (Jahn et al, 2018).

Toy univariate example revisited
Here we re-use case (3) of the toy example discussed in Section 3.1, for which MultiNest was shown to fail without applying PR. In this case, the true value of the unknown parameter is θ * = 40 and the number of observations is set to N = 20 (see Figure 4(a)).
We use power prior redefinition and consider the β values 0, 0.2, 0.4, 0.6, 0.8 and 1; note that β = 1 is equivalent to the original method implemented in the toy example, and β = 0 corresponds to using a uniform distribution as the modified prior. The range of the uniform prior for β = 0 is set as θ ∈ [0, 50] in this example. Figure 4 shows the performance of MultiNest assisted by the PR method. Panels (b) to (f) show the MultiNest performance with decreasing β. One sees that as β decreases, the posterior samples obtained approximate the true posterior with increasing accuracy, although in this extreme example one requires β = 0.4 or lower to obtain consistent results.
To evaluate the performance of the PR method further, MultiNest was run on 10 realisations for each value of β. The resulting histograms of MultiNest's posterior samples were then fitted with a standard Gaussian distribution. For each value of β, the average of the means of the fitted Gaussian distributions and the root mean squared error (RMSE) between these estimates and the true value are presented in Table 2, along with the average number of likelihood calls for MultiNest to converge; since the time spent for each likelihood calculation is similar, this quantity is proportional to the runtime. The RMSE clearly decreases as β decreases from unity to zero, which demonstrates that a wider prior allows MultiNest to obtain more accurate results, even in this extreme example of an unrepresentative prior. Also, one sees that the averaged number of likelihood evaluations also decreases significantly with β, so that the computational efficiency is also increased as the effective prior widens. Table 2 A numerical comparison of the results in the univariate toy example of the PR method for different values of β (where β = 1 corresponds to the standard method). The quantityμ denotes the averaged mean value of the fitted Gaussian distribution to the posterior histogram over 10 realisations. RMSE denotes the root mean squared error between the ground truth value andμ. N like is the averaged number of likelihood evaluations, and Z denotes the averaged estimated log-evidence and its uncertainty given by MultiNest. These results illustrate the general procedure mentioned at the end of Section 4.2, in which one obtains inferences for progressively smaller values of β, according to some pre-defined or dynamic 'annealing schedule', until the results converge to a statistically consistent solution. This is explored further in the example considered in the next section. Before moving on, however, it is also of interest to investigate the evidence values obtained with and without the PR method. For completeness, we reconsider all three cases of the toy example discussed in Section 3.1, namely: (1) θ * = 5; (2) θ * = 30; and θ * = 40. In each case, we calculate the mean and standard deviation of the log-evidence reported by MultiNest over 20 realisations of the data for β = 0, 0.2, 0.4, 0.6, 0.8, 1, respectively. The results are shown in Figure 5, in which the blue solid line and the light blue shaded area indicate, respectively, the average and standard deviation of the log-evidence values produced by MultiNest without PR (β = 1), and the red marker black cap errorbar shows the corresponding quantities produced using PR with other β-values.
For case (1), the red dashed curve fluctuates around the benchmark blue line at −22.0182 (β = 1 case), and the evidence estimates have similar size uncertainties, as one would expect. For case (2), however, one sees that the mean evidence values do change slightly as β is decreased from unity, converging on a final value for β < 0.8 that is ∼ 0.1 log-units larger than its mean value for β = 1. This indicates that case (2) also suffers (to a small extent) from the unrepresentative prior issue, despite this not being evident from the posterior samples plotted in Figure 1(d). For case (3), as expected, one sees that the mean log-evidence values change vastly as β is decreased from unity, converging for β < 0.6 on a value that is ∼ 500 log-units higher than for β = 1. This large difference means that the error-bars are not visible in this case, so the mean and standard deviation of the log-evidence for each β value are also reported in the last column of the Table 2.
These results demonstrate that the PR method also works effectively for evidence approximation in nested sampling, as well as producing posterior samples. Indeed, it also suggests that the evidence might be a useful statistic to monitor for convergence as one gradually lowers the value of β in the PR method.

Bivariate example
As our second example we consider a bivariate generalisation of our previous example, since it is straightforward to visualise. The bivariate case can easily be extended to higher dimensionality.
In each case, the MultiNest sampling parameters were set to N live = 100, efr = 0.8 and tol = 0.5 (see Feroz et al 2009 for details), and the algorithm was run to convergence. The results obtained without apply the PR method (which is equivalent to setting β = 1) are shown in Figure 6. One sees that the MultiNest samples are consistent with the true posterior distribution for case (1), but the sampler fails in cases (2) and (3) in which the ground truth lies far into the wings of the prior.
The MultiNest posterior samples obtained using the PR method, with β = 0.6, 0.3, 0.1, respectively, are shown in Figure 7 for case (2) and case (3). In each case, one sees that as β decreases the samples become consistent with the true posterior. In practice, it is thus necessary to reduce the value of β until the inferences converge to a sufficient accuracy. Table 3 summarises the inference accuracy and the computational efficiency for all three cases for Multi-Nest without PR (which corresponds to β = 1) and with PR for β = 0.4, 0.2, 0.1, 0.05, 10 −5 . One sees that for case (1) θ * = (0.5, 0.5), applying PR to Multi-Nest has only a weak effect on the RMSE performance and the number of likelihood evaluations, with both changing by about a factor of about two (in opposite directions) across the range of β values considered. For case (2) θ * = (1.5, 1.5) and case (3) θ * = (2, 2), however, MultiNest without PR suffers from the unrepresentative prior problem and the corresponding RMSE and number of likelihood evaluations are considerably higher than in case (1). Nonetheless, by combining MultiNest with the PR method, the RMSE and number of likelihood evaluations can be made consistent across the three cases considered. One sees that the RMSE decreases as β decreases and the maximum accuracy is obtained when β = 10 −5 (for which the modified prior is very close to uniform). This should be contrasted with the total number of likelihood evaluations, which increases as β decreases. Indeed, it is clear that the minimum number of likelihood evaluations are required for intermediate values of β. These results show that a reasonable compromise between accuracy and computational efficiency is obtained for β = 0.05 in this problem, which also provides the best consistency for both RMSE and the number of likelihood evaluations across all three cases.

Other sampling algorithms
Since our focus here is to introduce the PR method to improve NS performance in problems with unrepresentative priors, a full comparison between MultiNest and other competing sampling algorithms is beyond the scope of this paper. Nonetheless, we report here on some results of a brief such comparison on the same bivariate example. In particular, we perform a comparison of MultiNest with PR, in terms both of the RMSE and the number of likelihood evaluations, with our own implementation of standard MCMC sampling using the Metropolis-Hastings algorithm with a Gaussian proposal distribution and also with importance sampling (IS), using a standard IS implementation from Python package 'pypmc' (Jahn et al, 2018). The results obtained using MCMC and IS are shown in the final two columns of Table 3. For β = 0.05 Multi-Nest achieves relatively better and consistent RMSE and N like performance across the three cases. The number of posterior samples from the competing algorithms (MCMC and IS) are fixed to values around 1100 in order to obtain a similar number of likelihood evaluations as required by MultiNest for β = 0.05. For case (1), the performance of IS is comparable to that of MCMC. However, in cases (2) and (3), IS is comparable to Multi-Nest with β = 1, so it is clear that IS also suffers from the unrepresentative prior problem.
The detailed comparison of different sampling algorithms is a broad topic that has been widely discussed and explored in the literature. For example, importance sampling was formulated as a special case of bridge sampling, and was compared in (Gronau et al, 2017). An importance nested sampling was proposed to incorporate importance sampling into NS evidence calculation step in (Feroz et al, 2013). A comparison between NS and MCMC was discussed in (Allison and Dunkley, 2013). A review of importance sampling is presented in (Tokdar and Kass, 2010).

Diagnostics for bivariate example
We take the opportunity here to illustrate the diagnostics process discussed in Section 4.4 using the bivariate example. Since case (1) does not suffer from the unrepresentative prior problem, it can be treated as a reliable example and we assume that the 'available knowledge' is gained by analysing this case. As shown in Table  3, the number of likelihood evaluations (which is proportional to the runtime) for MultiNest without PR (β = 1) increases significantly from case (1) to case (3). Thus, the unrepresentative prior problem can be identified on-the-fly by monitoring the runtime. An onthe-fly convergence rate check may also be straightforwardly applied using existing rate of convergence methods (Süli and Mayers, 2003) to the problem. In either case, one may identify that case (2) and case (3) differ significantly from the available knowledge, and hence the PR method should be applied.
However, for some sampling methods, on-the-fly diagnostic of monitoring the runtime would fail in the case (adopted here) in which the number of likelihood evaluations is fixed. In this case, one must therefore rely on an after-run diagnostic, such as the KL divergence, which quantifies the differences between the assumed prior and the corresponding posterior obtained in the analysis. Figure 8 shows MultiNest samples from the prior and the posterior for case (1) and case (2), respectively, of the bivariate example. By computing the standard KL divergence, we find a value (termed KL score) of 169 for case (1) (the available knowledge) and 1833 for case (2). It is clear that the KL score for the unrepresentative prior problem is much larger than normal case, and so case (2) could be flagged as an outlier according to some predefined criterion on KL score. Similarly considerations apply to case (3).

Higher-dimensional examples
In order to investigate the performance of PR in higher dimensionality, we reconsider case (2) in the bivariate example, but extend the dimensionality over the range 3 to 15 dimensions. In particular, we consider the performance with β = 0.05, 0.1, and 0.4. Each of the experiments is repeated 20 times, and the test results are evaluated by calculating the mean and standard deviation of the RMSE over these 20 realisations. As shown in Figure 9 (a), with an increase of dimensionality, the RMSE error-bar undergoes an obvious increase for both β = 0.05 and 0.1 cases. For the case β = 0.4, the RMSE increases at lower dimensionality, but then remains at a stable level for higher dimensionality. Overall, the RMSE performance in higher dimensions is consistent with that in the bivariate example in terms of its order of magnitude, which demonstrates that the PR method is stable and effective for problems with higher dimensionality. Figure 9 (b) shows a set of equivalent plots for the number of likelihood evaluations. This clearly shows that for a smaller β value MultiNest makes a larger number of likelihood evaluations. This is not surprising as a smaller β corresponds to a broader modified prior space. We note that the number of likelihood evaluations required for β = 0.05 is almost twice that for β = 0.1. Figure 10 shows the RMSE comparison between Multi-Nest with PR (β = 0.05) and MCMC methods for the same higher dimensional examples. Note that the RMSE is computed using a comparable number of likelihood evaluations for the two methods for each dimensionality. As can be observed from the figure, MCMC remains stable and accurate (albeit with a slight increase in RMSE with dimension), but has a higher RMSE than MultiNest with PR across the dimensionalities considered. By contrast, for MultiNest with PR, the RMSE increases more noticably with the number of dimensions, as might be expected from a NS algorithm that is based on a form rejection sampling.

Non-Gaussian bivariate example
As our final numerical example, we consider a non-Gaussian bivariate likelihood function. In particular, we adapt the Gaussian bivariate likelihood considered in Table 4 The performance of MultiNest with and without PR method (for various values of β, and both 100 live samples) for all three cases of the non-Gaussian bivariate example. The top half of the table is a comparison of RMSE, and the second half is the number of likelihood evaluations (N like ) per individual algorithm run. For β = 0.05 (highlighted in bold) MultiNest achieves relatively better and consistent RMSE and N like performance across the three cases.

RMSE
MN (β = 1) β = 0.4 β = 0.2 β = 0.1 β = 0.05 β = 10 −5 Case (1)  Section 5.2 by replacing the product of Gaussian distributions in each dimension by a product of Laplace distributions Laplace(µ, b), so that in each dimension the Gaussian form (5) is re-written as: where x n is the nth measurement (although N = 1 in this example), which acts as the location parameter similar to in a Gaussian distribution, and b is the scale parameter in the Laplace distribution analogous to σ ξ in (5). We choose the Laplace distribution as our non-Gaussian test example since: (1) it is valid for both positive and negative values of the parameter θ, unlike Beta/Gamma distributions; and (2) a Laplace distribution with a small b-value has a similar tail to that of a Gaussian (i.e. it is not heavy-tailed), which facilitates easier comparison. The prior distribution is identical to that used Section 5.2, i.e. the same Gaussian distribution. Indeed, all of the other experimental settings are kept the same as those in Section 5.2, and we again consider MultiNest with and without PR method in all three cases.
The results of the analysis are given in Table 4 for runs with N live = 100. Comparing the N live = 100 results with the corresponding ones given in Table 3 for the Gaussian bivariate example, ones sees that the trends for both RMSE and N like are similar to those in the Gaussian bivariate cases, but are in general higher for the Laplace distribution. This is because the peak of the Laplace distribution is sharper than that of a Gaussian. Again reasonable results are obtained for β = 0.05. Figure 11 shows the RMSE resulting from different N live values for β = 0.05 and 10 −5 , respectively. Comparing these RMSE values with those given in Table 3, which were obtained for the Gaussian bivariate example with N live = 100, one sees that higher N live values are required for the Laplace distribution to achieve similar levels of accuracy.

Conclusions
This paper addresses the unrepresentative prior problem in Bayesian inference problems using NS, by introducing the posterior repartitioning method.
The key advantages of the method are that: (i) it is general in nature and can be applied to any such inference problem; (ii) it is simple to implement; and (iii) the posterior distribution is unaltered and hence so too are the inferences. The method is demonstrated in univariate and bivariate numerical examples on Gaussian posteriors, and its performance is further validated and compared with MCMC sampling methods in examples up to 15 dimensions. The method is also tested on a non-Gaussian bivariate example. In all cases, we demonstrate that NS algorithms, assisted by the PR method, can achieve accurate posterior estimation and evidence approximation in problems with an unrepresentative prior.
The proposed scheme does, however, have some limitations: (i) if the prior and likelihood are extremely widely separated, the sampling can still be inefficient and slow, because of the large augmented search space for very small β; (ii) the approach cannot be readily applied to problems with discrete parameters; and (iii) the normalisation of the modified prior will in general not be possible analytically, but require numerical integration.