CosmoCoffee

[Antony Lewis]

This paper discusses the use of Nested Sampling for cosmological model selection and parameter analysis. This is certainly a neat method, and nice to see it being put to use.

I hope to discuss some of these things at the Sussex conference next week. A couple of things to consider:

A. Intermediate steps in this method require sampling from the prior P(x) conditional on the likelihood L(x) being greater than a particular value L_j. To implement this sampling the authors seem to

1. Construct an ellipse than encloses the current set of samples
2. Expand the ellipse by some constant to allow for the fact that true likelihood contours are not elliptical
3. Importance sample within this expanded ellipse (i.e. take samples until one sample lies in the L(x) > L_j region).

This is in contrast to the original Skilling method where he suggests using MCMC for this sampling step. My worry about the method used here as a general method is the following. In n-dimensions, the volume of a sphere of radius r is ~ rⁿ. If the sphere is expanded by a factor f, the fraction of the new volume contained in the target volume is (1/f)ⁿ, i.e. the acceptance rate is approximately 1/fⁿ and the number of likelihood calculations to get the new sample is ~ fⁿ. For n=5, f=1.5 this gives an acceptance rate ~ 0.13, similar to the 20% they report. For n=7, f=1.8 it gives a rate of 0.016. Since fⁿ grows exponentially with dimension, this sampling method would seem to get exponentially bad in high dimensions. i.e. perhaps for n>~8 using MCMC would be a much better idea (though maybe fine for n<=8 considered in this paper). Also expanding the elliptical region is of course not guaranteed to enclose all points in the L(x) > L_j region, so the method is not formally correct, and can only be tested in each case by expanding my a much larger factor (very slow!) and checking for stability of the result.

B. Is it useful to use different 'priors'? e.g.



so one could perhaps apply this method sampling from P'(θ) rather than the original P(θ) [and using the appropriately modified likelihood]. i.e. perhaps one could use some approximation to L^β that can be quickly normalized.

[Hans Kristian Eriksen]

[Antony Lewis]

I should think you could use almost anything for simple cases, e.g. just the standard CosmoMC sampler (see the CosmoMC paper, notes and readme for references; basically Metropolis-Hastings with Gaussian jumps in random 2D subspaces).

On the other hand I'm not at all sure how well nested sampling would work on a general distribution. I'd have thought sampling from the sharply-truncated prior could give a very odd disjoint shape to sample from even if the posterior distribution is only weakly multimodal itself - in general very difficult. An advantage of using theromodynamic integration is that it least it is pretty robust to odd distribution shapes.

[Hans Kristian Eriksen]

[Antony Lewis]

I think that is right - MCMC will certainly involve lots of likelihood evaluations. One could however hope that it would scale with dimension less badly than exponentially.

On the other hand David Mackay pointed out that for a n-D ellipse, almost all the probability volume is around the surface, so it is indeed also very hard for MCMC to find a new sample. He suggested starting chains at one of the just-rejected points - which would probably be valid for some class of unimodal distributions, though certainly not generally.

In general it seems to me that nested sampling is difficult in high dimensions, and may indeed not be better than thermodynamic integration - especially for odd-shaped distributions.

[Pia Mukherjee]

We tried MCMC-ing (from any one of the remaining points, or from the last rejected point) and found that in general there remained a bias in the answer. In comparison, the sampling from the spheroid scheme was efficient and gave a very accurate and unbiased answer. True, this scheme for finding new points will become less efficient with increasing dimentionality.

As regards comparison with thermodynamic integration (TI), we seemed to require orders of magnitude fewer likelihood evaluations to reach a given final accuracy (based on TI requiring ~10⁷ likelihood evaluations to get to δ(ln E) of 0.1 for a 5 or 6 parameter model. This is based on estimates from Beltran et al. and private communications. Whether TI can be made more efficient, through more tuned thermodynamic scheduling, I dont know.

[Andrew Liddle]

Hi there,

I think it is pretty generic that MCMC won't work to generate new live points, which is what we found in testing various methods.

The challenge is to sample **uniformly** from the remaining prior, which becomes a small part of the parameter volume as the code progresses. With an MCMC method, the samples will only find (or remain in) the high likelihood region if the sampler is given knowledge of what the likelihood is, but as soon as you do that the sampling stops being uniform and biases the answer.

It would be interesting to explore thermodynamic integration more extensively for large-dimensional problems, where I agree that nested sampling may scale badly. In Beltran et al astro-ph/0501477 we found that Therm Int was very inefficient for cosmological problems, and we subsequently found the same even for simple gaussian likelihoods. But we focussed on one particular annealing schedule and significant optimization may be possible.

best,

Andrew

[Antony Lewis]

[Andrew Liddle]

Dear Antony,

Yes, you're right. That should work but the number of MCMC steps needed to converge might indeed be quite large. Worth testing.
Your proposal is somewhat different from the failed MCMC attempts we tried.

best,

Andrew