Example python likelihoods and priors¶
Some example PolyChord python likelihoods and priors. You can also add your own!
python_likelihoods¶
Loglikelihood functions for use with PolyChord’s python interface.
PolyChord >= v1.14 requires likelihoods to be callables with parameter and return signatures:
Parameters¶
theta: float or 1d numpy array
Returns¶
- logl: float
Loglikelihood.
- phi: list of length nderived
Any derived parameters.
We use classes with the loglikelihood defined in the __call__ property, as this provides convenient way of storing other information such as hyperparameter values. These objects can be used in the same way as functions due to python’s “duck typing” (alternatively you can define likelihoods using functions).
-
class
dyPolyChord.python_likelihoods.
Gaussian
(sigma=1.0, nderived=0)[source]¶ Symmetric Gaussian loglikelihood centered on the origin.
-
class
dyPolyChord.python_likelihoods.
GaussianMix
(sep=4, weights=(0.4, 0.3, 0.2, 0.1), sigma=1, nderived=0)[source]¶ Gaussian mixture likelihood in \(\ge 2\) dimensions with up to 4 compoents.
Each component has the same standard deviation \(\sigma\), and they their centres respectively have \((\theta_1, \theta_2)\) coordinates:
(0, sep), (0, -sep), (sep, 0), (-sep, 0).
-
__call__
(self, theta)[source]¶ Calculate loglikelihood(theta), as well as any derived parameters.
N.B. this loglikelihood is not normalised.
- Parameters
- theta: float or 1d numpy array
- Returns
- logl: float
Loglikelihood
- phi: list of length nderived
Any derived parameters.
-
__init__
(self, sep=4, weights=(0.4, 0.3, 0.2, 0.1), sigma=1, nderived=0)[source]¶ Define likelihood hyperparameter values.
- Parameters
- sep: float
Distance from each Gaussian to the origin.
- weights: iterable of floats
weights of each Gaussian component.
- sigma: float
Stanard deviation of Gaussian components.
- nderived: int, optional
Number of derived parameters.
-
-
class
dyPolyChord.python_likelihoods.
GaussianShell
(sigma=0.2, rshell=2, nderived=0)[source]¶ Gaussian Shell loglikelihood centred on the origin.
-
__call__
(self, theta)[source]¶ Calculate loglikelihood(theta), as well as any derived parameters.
N.B. this loglikelihood is not normalised.
- Parameters
- theta: float or 1d numpy array
- Returns
- logl: float
Loglikelihood
- phi: list of length nderived
Any derived parameters.
-
__init__
(self, sigma=0.2, rshell=2, nderived=0)[source]¶ Set up likelihood object’s hyperparameter values.
- Parameters
- sigma: float, optional
Standard deviation of Gaussian (the same for each parameter).
- rshell: float, optional
Distance of shell peak from origin.
- nderived: int, optional
Number of derived parameters.
-
-
class
dyPolyChord.python_likelihoods.
LogGammaMix
(nderived=0)[source]¶ The loggamma mix used in Beaujean and Caldwell (2013) and Feroz et al. (2013).
-
class
dyPolyChord.python_likelihoods.
Rastrigin
(a=10, nderived=0)[source]¶ Rastrigin loglikelihood as described in “PolyChord: next-generation nested sampling” (Handley et al., 2015).
-
class
dyPolyChord.python_likelihoods.
Rosenbrock
(a=1, b=100, nderived=0)[source]¶ Rosenbrock loglikelihood as described in “PolyChord: next-generation nested sampling” (Handley et al., 2015).
-
dyPolyChord.python_likelihoods.
log_gaussian_pdf
(theta, sigma=1, mu=0, ndim=None)[source]¶ Log of uncorrelated Gaussian pdf.
- Parameters
- theta: float or 1d numpy array
- sigma: float, optional
- mu: float, optional
- ndim: int, optional
- Returns
- logl: float
Loglikelihood.
-
dyPolyChord.python_likelihoods.
log_loggamma_pdf
(theta, alpha=1, beta=1)[source]¶ Loglikelihood for multidimensional loggamma distribution, with each component of theta having an independent loggamma distribution.
Always returns a float.
- Parameters
- theta: float or numpy array
- alpha: float, optional
- beta: float, optional
- Returns
- logl: float
Loglikelihood.
-
dyPolyChord.python_likelihoods.
log_loggamma_pdf_1d
(theta, alpha=1, beta=1)[source]¶ 1d gamma distribution, with each component of theta independently having PDF:
\[f(x) = \frac{ \mathrm{e}^{\beta x} \mathrm{e}^{-\mathrm{e}^x / \alpha} }{\alpha^{\beta} \Gamma(\beta)}\]This function returns \(\log f(x)\).
- Parameters
- theta: float or array
Where to evaluate \(\log f(x)\). Values \(\in (-\inf, \inf)\).
- alpha: float, > 0
Scale parameter
- beta: float, > 0
Shape parameter
- Returns
- logl: same type as theta
python_priors¶
Python priors for use with PolyChord.
PolyChord v1.14 requires priors to be callables with parameter and return types:
Parameters¶
- hypercube: float or 1d numpy array
Parameter positions in the prior hypercube.
Returns¶
- theta: float or 1d numpy array
Corresponding physical parameter coordinates.
Input hypercube values numpy array is mapped to physical space using the inverse CDF (cumulative distribution function) of each parameter. See the PolyChord papers for more details.
This module use classes with the prior defined in the cube_to_physical function and called with __call__ property. This provides convenient way of storing other information such as hyperparameter values. The objects be used in the same way as functions due to python’s “duck typing” (or alternatively you can just define prior functions).
The BlockPrior class allows convenient use of different priors on different parameters.
Inheritance of the BasePrior class allows priors to:
1. have parameters’ values sorted to give an enforced order. Useful when the parameter space is symmetric under exchange of variables as this allows the space to be explored to be contracted by a factor of N! (where N is the number of such parameters);
adaptively select the number of parameters to use.
You can ignore these if you don’t need them.
-
class
dyPolyChord.python_priors.
BasePrior
(adaptive=False, sort=False, nfunc_min=1)[source]¶ Base class for Priors.
-
__call__
(self, cube)[source]¶ Evaluate prior on hypercube coordinates.
- Parameters
- cube: 1d numpy array
Point coordinate on unit hypercube (in probabily space). Note this variable cannot be edited else PolyChord throws an error.
- Returns
- theta: 1d numpy array
Physical parameter values for prior.
-
__init__
(self, adaptive=False, sort=False, nfunc_min=1)[source]¶ Set up prior object’s hyperparameter values.
- Parameters
- adaptive: bool, optional
- sort: bool, optional
- nfunc_min: int, optional
-
cube_to_physical
(self, cube)[source]¶ Map hypercube values to physical parameter values.
- Parameters
- cube: 1d numpy array
Point coordinate on unit hypercube (in probabily space). See the PolyChord papers for more details.
- Returns
- theta: 1d numpy array
Physical parameter values corresponding to hypercube.
-
-
class
dyPolyChord.python_priors.
BlockPrior
(prior_blocks, block_sizes)[source]¶ Prior object which applies a list of priors to different blocks within the parameters.
-
class
dyPolyChord.python_priors.
Exponential
(lambd=1.0, **kwargs)[source]¶ Exponential prior.
-
class
dyPolyChord.python_priors.
Gaussian
(sigma=10.0, half=False, mu=0.0, **kwargs)[source]¶ Symmetric Gaussian prior centred on the origin.
-
__init__
(self, sigma=10.0, half=False, mu=0.0, **kwargs)[source]¶ Set up prior object’s hyperparameter values.
- Parameters
- sigma: float
Standard deviation of Gaussian prior in each parameter.
- half: bool, optional
Half-Gaussian prior - nonzero only in the region greater than mu. Note that in this case mu is no longer the mean and sigma is no longer the standard deviation of the prior.
- mu: float, optional
Mean of Gaussian prior.
- kwargs: dict, optional
See BasePrior.__init__ for more infomation.
-
cube_to_physical
(self, cube)[source]¶ Map hypercube values to physical parameter values.
- Parameters
- cube: 1d numpy array
Point coordinate on unit hypercube (in probabily space). See the PolyChord papers for more details.
- Returns
- theta: 1d numpy array
Physical parameter values corresponding to hypercube.
-
-
class
dyPolyChord.python_priors.
PowerUniform
(minimum=0.1, maximum=2.0, power=-2, **kwargs)[source]¶ Uniform in theta ** power
-
class
dyPolyChord.python_priors.
Uniform
(minimum=0.0, maximum=1.0, **kwargs)[source]¶ Uniform prior.
-
dyPolyChord.python_priors.
adaptive_transform
(cube, sort=True, nfunc_min=1)[source]¶ Tranform first parameter (nfunc) to uniform in (nfunc_min, nfunc_max) and, if required, perform forced identifiability transform on the next nfunc parameters only.
- Parameters
- cube: 1d numpy array
Point coordinate on unit hypercube (in probabily space).
- Returns
- ad_cube: 1d numpy array
First element is physical coordinate of nfunc parameter, other elements are cube coordinates with any forced identifiability transform already applied.
-
dyPolyChord.python_priors.
forced_identifiability
(cube)[source]¶ Transform hypercube coordinates to enforce identifiability.
For more details see: “PolyChord: next-generation nested sampling” (Handley et al. 2015).
- Parameters
- cube: 1d numpy array
Point coordinate on unit hypercube (in probabily space).
- Returns
- ordered_cube: 1d numpy array