Planck (via clik)

This plugin is an interface between the Planck likelihood code clik and CosmoSlik. You need clik already installed on your machine, which you can get from here.

You also need to download the "clik files" for whichever likelihoods you would like to use. You can find these here under "Likelihoods" / "Notes".

Quickstart

CosmoSlik provides several plugins which wrap clik and have all the necessary nuisance parameters set up for particular data files. You can use them in your script by adding something like the following to your __init__,

# set up cosmological params and solver
self.cosmo = models.cosmology("lcdm")
self.cmb = models.classy()

# load Planck clik file and set up nuisance parameters
self.clik = likelihoods.planck.planck_2015_highl_TT(
    clik_file="plik_dx11dr2_HM_v18_TT.clik/",
)

then compute the likelihood in __call__ by calling clik with a parameter cmb of the kind returned by CAMB or CLASS,

# compute likelihood
self.clik(self.cmb(**self.cosmo))

The generic clik wrapper

Using the SlikPlugin named clik, we can load up any generic clik file. Supposing we've downloaded the file plik_lite_v18_TT.clik, we can load it in via,


In [1]:
%pylab inline
sys.path = sys.path[1:]
from cosmoslik import *


Populating the interactive namespace from numpy and matplotlib

In [2]:
clik = likelihoods.planck.clik(
    clik_file="plik_lite_v18_TT.clik/",
    A_Planck=1
)
clik


Out[2]:
{'A_Planck': 1,
 'auto_reject_errors': False,
 'clik': <clik.lkl.clik at 0x7f03a8eda510>}

Note that we gave it a parameter A_Planck. Most clik files have extra nuisance parameters, which you can list (for a given file) with,


In [3]:
clik.clik.get_extra_parameter_names()


Out[3]:
('A_Planck',)

You should attach parametes with these names to the clik object as we have done above (usually in a script these will be sampled parameters).

With the clik object created, we can call it to compute the likelihood. The function expects a parameter cmb of the kind returned by CAMB or CLASS.


In [4]:
cmb = models.classy(lmax=3000)()
cmb


Out[4]:
{'BB': array([  0.00000000e+00,   0.00000000e+00,   1.77707779e-06, ...,
          1.41341956e-02,   1.41160893e-02,   1.40979982e-02]),
 'EE': array([ 0.        ,  0.        ,  0.05413355, ...,  0.92829409,
         0.92452331,  0.92077857]),
 'PP': array([  0.00000000e+00,   0.00000000e+00,   6.55525164e+04, ...,
          4.37564298e-04,   4.36858020e-04,   4.36153030e-04]),
 'TE': array([ 0.        ,  0.        ,  3.57593976, ..., -1.57751022,
        -1.57159761, -1.56562643]),
 'TP': array([  0.00000000e+00,   0.00000000e+00,   3.60597308e+03, ...,
          4.68131564e-05,   4.67227124e-05,   4.66374675e-05]),
 'TT': array([    0.        ,     0.        ,  1110.41116627, ...,    26.09287144,
           26.04596539,    25.99892254])}

Here's the negative log likelihood:


In [5]:
clik(cmb)


Out[5]:
201.12250756838722

Putting it all together, a simple script which runs this likelihood would look like:


In [6]:
class planck(SlikPlugin):

    def __init__(self, **kwargs):
        super().__init__()
        
        # load Planck clik file and set up nuisance parameters
        self.clik = likelihoods.planck.clik(
            clik_file="plik_lite_v18_TT.clik/",
            
            # sample over nuisance parameter
            A_Planck=param(start=1, scale=0.0025, gaussian_prior=(1,0.0025))
        )
        
        # set up cosmological params and solver
        self.cosmo = models.cosmology("lcdm")
        self.cmb = models.classy(lmax=3000)
        
        self.sampler = samplers.metropolis_hastings(self)

    def __call__(self):
        # compute likelihood
        return self.clik(self.cmb(**self.cosmo))

In [7]:
s = Slik(planck())
lnl, e = s.evaluate(**s.get_start())
lnl


Out[7]:
956.87241321269414

Ready-to-go wrappers for specific clik files

The previous example was easy because there was one single nuisance parameter, A_Planck. Other clik files have many more nuisance parameters, which must all be sampled over and in some cases have the right priors applied (which you can read about here), otherwise you will not get the right answer.

This is, of course, a huge pain.

For this reason, CosmoSlik comes with several SlikPlugins already containing the correct sampled nuisance parameters for many of these clik files, making writing a script extremely easy. For example, here is the source code for one such plugin, planck_2015_highl_TT:

param = param_shortcut('start','scale')

class planck_2015_highl_TT(clik):

    def __init__(
        self,
        clik_file,
        A_cib_217        = param(60,  10,     range=(0,200)),
        A_planck         = param(1,   0.0025, range=(0.9,1.1), gaussian_prior=(1,0.0025)),
        A_sz             = param(5,   3,      range=(0,10)),
        calib_100T       = param(1,   0.001,  range=(0,3),     gaussian_prior=(0.999,0.001)),
        calib_217T       = param(1,   0.002,  range=(0,3),     gaussian_prior=(0.995,0.002)),
        cib_index        = -1.3,   
        gal545_A_100     = param(7,   2,      range=(0,50),    gaussian_prior=(7,2)),
        gal545_A_143     = param(9,   2,      range=(0,50),    gaussian_prior=(9,2)),
        gal545_A_143_217 = param(21,  8.5,    range=(0,100),   gaussian_prior=(21,8.5)),
        gal545_A_217     = param(80,  20,     range=(0,400),   gaussian_prior=(80,20)),
        ksz_norm         = param(2,   3,      range=(0,10)),
        ps_A_100_100     = param(250, 30,     range=(0,4000)),
        ps_A_143_143     = param(45,  10,     range=(0,4000)),
        ps_A_143_217     = param(40,  10,     range=(0,4000)),
        ps_A_217_217     = param(90,  15,     range=(0,4000)),
        xi_sz_cib        = param(0.5, 0.3,    range=(0,1)),
    ):
        super().__init__(**arguments())

As you can see, all the sampled parameters as automatically set, including ranges and priors. The script to use this likelihood is then extremely simple:


In [8]:
class planck(SlikPlugin):

    def __init__(self):
        super().__init__()
        
        # load Planck clik file and set up nuisance parameters
        self.clik = likelihoods.planck.planck_2015_highl_TT(
            clik_file="plik_dx11dr2_HM_v18_TT.clik/",
        )
        
        # set up cosmological params and solver
        self.cosmo = models.cosmology("lcdm")
        self.cmb = models.classy(lmax=3000)
        
        self.sampler = samplers.metropolis_hastings(self)

    def __call__(self):
        # compute likelihood
        return self.clik(self.cmb(**self.cosmo))

In [9]:
s = Slik(planck())
lnl, e = s.evaluate(**s.get_start())
lnl


Out[9]:
2672.4139027829988

Common calibration parameters

Despite that the Planck likelihood is broken up into different pieces, they sometimes share the same calibration parameters. To apply this correctly in your script, just define one single sampled calibration parameter, then in your __call__, set it across all the different likelihoods.


In [10]:
class planck(SlikPlugin):

    def __init__(self):
        super().__init__()
        
        # set up low and high L likelihood
        self.highl = likelihoods.planck.planck_2015_highl_TT(
            clik_file="plik_dx11dr2_HM_v18_TT.clik/",
        )
        self.lowl = likelihoods.planck.planck_2015_lowl_TT(
            clik_file="commander_rc2_v1.1_l2_29_B.clik/",
            A_planck=None, #turn off this cal parameter, use the one from self.highl
        )
        
        # set up cosmological params and solver
        self.cosmo = models.cosmology("lcdm")
        self.cmb = models.classy(lmax=3000)
        
        self.sampler = samplers.metropolis_hastings(self)

    def __call__(self):
        # set the calibration parameters the same
        self.lowl.A_planck = self.highl.A_planck 
        
        # compute likelihood
        cmb = self.cmb(**self.cosmo)
        return self.lowl(cmb) + self.highl(cmb)