jupyter_notebooks/Nested_Sampling_classes.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Nested Samplers available in Gleipnir\n",
    "\n",
    "Gleipnir currently has five classes that can be used to create Nested Sampling objects: one for using a built-in Nesting Sampling implementation and four interfaces to external Nested Sampling packages. In this notebook, we will explore each of these Nested Samplers and how to use them. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Gleipnir currently has six Nested Sampling classes:\n",
    "\n",
    "  * NestedSampling - Gleipnir's built-in implementation of Nested Sampling.\n",
    "  * MultiNestNestedSampling - Wrapper on top of PyMultiNest for running MultiNest\n",
    "  * PolyChordNestedSampling - Wrapper on top of pypolychord for running PolyChord\n",
    "  * dyPolyChordNestedSampling - Wrapper on top of dyPolyChord for running dyPolyChord\n",
    "  * DNest4NestedSampling - Wrapper on top of dnest4 (python bindings of DNest4) for running DNest4\n",
    "  * NestleNestedSampling - Wrapper on top of Nestle for running Nestle\n",
    "  \n",
    "We'll cover each one in more detail in the following sections.  "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Model: Egg Carton likelihood\n",
    "\n",
    "For the purposes of this tutorial, we will again use the Egg Carton likelihood landscape model (as in the ([Intro to Nested Sampling Notebook](https://github.com/LoLab-VU/Gleipnir/blob/master/jupyter_notebooks/Intro_to_Nested_Sampling_with_Gleipnir.ipynb)).  \n",
    "\n",
    "The model is typically two-dimensional (two parameters) and the landscape generated by the likelihood function is a multi-modal egg carton-like shape; see slide 15 of [this pdf](http://www.nbi.dk/~koskinen/Teaching/AdvancedMethodsInAppliedStatistics2016/Lecture14_MultiNest.pdf) for a visualization of the likelihood landscape. The parameters are each defined on \\[0:10pi\\] with uniform priors. \n",
    "\n",
    "Here is the Egg Carton loglikelihood function which returns the natural logarithm of the likelihood for a given parameter vector:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Import NumPy\n",
    "import numpy as np\n",
    "# Define the loglikelihood function.\n",
    "def loglikelihood(parameter_vector):\n",
    "    chi = (np.cos(parameter_vector)).prod()\n",
    "    return (2. + chi)**5"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Sampled Parameters\n",
    "\n",
    "Now that we have our loglikelihood function, let's look at how to define parameters for sampling during the Nested Sampling run. The parameters that are sampled are defined by a list of SampledParameter class instances. The SampledParameter class object stores data on the name of the parameter and the parameter's prior probability distribution. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Import the SampledParameters class.\n",
    "from gleipnir.sampled_parameter import SampledParameter"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "A new SampledParameter needs two arguments: a name and a object defining the prior.\n",
    "\n",
    "For priors we can use frozen RV objects from scipy.stats; in special cases you could also write your own prior distribution class objects, but for most purposes scipy.stats distributions will be sufficient."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Let's import the uniform distribution.\n",
    "from scipy.stats import uniform"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Now we'll create our list sampled parameters.\n",
    "# There are two parameters 'x' and 'y',each with a uniform prior on [0:10pi]. \n",
    "sampled_parameters = list()\n",
    "sampled_parameters.append(SampledParameter(name='x', prior=uniform(loc=0.0,scale=10.0*np.pi)))\n",
    "sampled_parameters.append(SampledParameter(name='y', prior=uniform(loc=0.0,scale=10.0*np.pi)))  "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now we've defined our list of paramters that are to be sampled and their prior probability distributions.  We'll use these sampled parameters for each of the Nested Samplers. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 1. NestedSampling\n",
    "The first Nested Sampler is Gleipnir's built-in implementation of the Nested Sampling algorithm, NestedSampling. It is imported from the gleipnir.nestedsampling module. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "from gleipnir.nestedsampling import NestedSampling"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This Nested Sampler uses a plug-in style approach to sampling (i.e., replacing dead points) and stopping criterion (i.e., when to exit the Nested Sampling run). These are passed in to the NestedSampling class initialization as optional keyword arguments. Although there are defaults set for these if we don't pass anything in, for the puposes of this tutorial we will go over each and explicitly set them."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Sampler\n",
    "The sampler is used in the classic Nested Sampling algorithm to replace dead points during each nested iteration. The samplers are imported from the nestedsampling.samplers sub-module:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Import the sampler we want to use during the Neseted Sampling run.\n",
    "from gleipnir.nestedsampling.samplers import MetropolisComponentWiseHardNSRejection"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Currently, Gleipnir just has the one sampler: MetropolisComponentWiseHardNSRejection. This sampler uses a [Metropololis Monte Carlo](http://xbeams.chem.yale.edu/~batista/vaa/node42.html) scheme adapted for Nested Sampling with a hard Nested likelihood level rejection criterion. During the Nested Sampling iteration the most recent dead point is replaced with a survivor that is  then modified using the sampler. More samplers may be added in the future.\n",
    "\n",
    "Now we can initialize the sampler:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Initialize the sampler.\n",
    "sampler = MetropolisComponentWiseHardNSRejection(iterations=10, tuning_cycles=1)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Here the iterations=10 are the number of total component wise update cycles, so the value of 10 will yield 20 component-wise Monte Carlo trial moves (i.e., iterations*(number of sampled parameters)). tuning_cycles=1 sets the number of trial move size tuning cycles; each tuning_cycle is 20 iterations. Note that MetropolisComponentWiseHardNSRejection(10, tuning_cyles=1) is the default value that will be used if we don't explicitly set the sampler."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Stopping Criterion\n",
    "\n",
    "The stopping criterion sets the method used to determine when terminate the Nested Sampling iterations. They are imported from the nestedsampling.stopping_criterion modulue. There are currenlty two criterion which can be used:\n",
    "  * NumberOfIterations - Stop after a fixed number of nested iterations.\n",
    "  * RemainingPriorMass - Stop after the remaining fraction of prior mass is less than or equal to a preset threshold.\n",
    "Here, we'll use the fixed number of iterations:  "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Import the stopping criterion object. In this case, we'll use a fixed number of iterations.\n",
    "from gleipnir.nestedsampling.stopping_criterion import NumberOfIterations\n",
    "# Initialize the stopping criterion -- We'll stop after 1000 Nested Sampling iterations.\n",
    "stopping_criterion = NumberOfIterations(1000)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Note that NumberOfIterations(1000) is the default that will be used by the Nested Sampler if an explicit value is not set."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now that we've got the sampled parameters, sampler, and stopping criterion, all we need now is the Nested Sampling population size. Let's go ahead and set it:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Set the NS population size.\n",
    "population_size=500"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now we create the instance of the NestedSampling class:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [],
   "source": [
    "NS = NestedSampling(sampled_parameters,\n",
    "                    loglikelihood,\n",
    "                    population_size,\n",
    "                    sampler=sampler,\n",
    "                    stopping_criterion=stopping_criterion)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Then the Nested Sampling is launched with the run function:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [],
   "source": [
    "log_evidence, log_evidence_error = NS.run(verbose=False)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Then we can check the output:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      " Ln(Evidence): 235.14699033558188+-0.06321396463872749\n"
     ]
    }
   ],
   "source": [
    "print(\" Ln(Evidence): {}+-{}\".format(log_evidence, log_evidence_error))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In addition to the log_evidence and log_evidence_error estimates returned by the run function, the NestedSampling Nested Sampler has the following accesible properties:\n",
    "  * evidence\n",
    "  * evidence_error\n",
    "  * information"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [],
   "source": [
    "evidence = NS.evidence\n",
    "evidence_error = NS.evidence_error\n",
    "information = NS.information"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The NestedSampling Nested Sampler also has the following functions:\n",
    "  * posteriors - Estimates of the posterior marginal probability distributions of each parameter. \n",
    "  * posterior_moments - Estimates of the first four moments of each posterior marginal probability distribution of the parameters.  \n",
    "  * akaike_ic - Estimate of the Akaike Information Criterion.\n",
    "  * bayesian_ic - Estimate of the Bayesian Information Criterion.\n",
    "  * deviance_ic - Estimate of the Deviance Information Criterion.\n",
    "  * best_fit_likelihood - Get the maximum likelihood parameter vector.\n",
    "  * best_fit_posterior - Get highest posterior probability parameter vector."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 2. MultiNestNestedSampling\n",
    "\n",
    "Gleipnir also provides a Nested Sampler class, MultiNestNestedSampling, that is an interface to [MultiNest](https://academic.oup.com/mnras/article/398/4/1601/981502); this class is a wrapper on top of [PyMultiNest](https://github.com/JohannesBuchner/PyMultiNest). Use of MultiNest via Gleipnir requires separate building and installation of both MultiNest and the Python wrapper PyMultiNest; see Gleinir's [README](https://github.com/LoLab-VU/Gleipnir#multinest) for the links to those instructions. \n",
    "\n",
    "It is imported from the multinest module:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [],
   "source": [
    "from gleipnir.multinest import MultiNestNestedSampling"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Then the instance is created with a similar input pattern as the NestedSampling, but only requires the sampled parameters, loglikelihood, and population size:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [],
   "source": [
    "MNNS = MultiNestNestedSampling(sampled_parameters, loglikelihood, population_size)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "However, the MulitNestNestedSampling object does have some extra keyword arguments that can be passed in to alter the behavior of the MultiNest run (internally passed along to PyMultiNest):\n",
    "  * importance_nested_sampling (bool): Should MultiNest use\n",
    "    Importance Nested Sampling (INS). Default: True\n",
    "  * constant_efficiency_mode (bool): Should MultiNest run in\n",
    "    constant sampling efficiency mode. Default: False\n",
    "  * sampling_efficiency (float): Set the MultiNest sampling\n",
    "    efficiency. 0.3 is recommended for evidence evaluation,\n",
    "    while 0.8 is recommended for parameter estimation.\n",
    "    Default: 0.8\n",
    "  * resume (bool): Resume from a previous MultiNest run (using\n",
    "    the last saved checkpoint in the MultiNest output files).\n",
    "    Default: True\n",
    "  * write_output (bool): Specify whether MultiNest should write\n",
    "    to output files. True is required for additional\n",
    "    analysis. Default: True\n",
    "  * multimodal (bool): Set whether MultiNest performs mode\n",
    "    separation. Default: True\n",
    "  * max_mode (int): Set the maximum number of modes allowed in\n",
    "    mode separation (if multimodal=True). Default: 100\n",
    "  * mode_tolerance (float): A lower bound for which MultiNest will\n",
    "    use to separate mode samples and statistics with\n",
    "    log-evidence value greater the given value.\n",
    "    Default: -1e90\n",
    "  * n_clustering_params (int): If multimodal=True, set the number\n",
    "    of parameters to use in clustering during mode separation.\n",
    "    If None, then MultiNest will use all the paramters for\n",
    "    clustering during mode separation. If\n",
    "    n<(number of sampled parameters), then MultiNest will only\n",
    "    use a subset composed of the first n parameters for\n",
    "    clustering during mode separation. Default: None\n",
    "  * null_log_evidence (float): If multimodal=True, a lower bound\n",
    "    for which MultiNest can use to separte mode samples and\n",
    "    statistics with a local log-evidence value greater than the\n",
    "    given bound. Default: -1.e90\n",
    "  * log_zero (float): Set a threshold value for which points with\n",
    "    a loglikelihood less than the given value will be ignored\n",
    "    by MultiNest. Default: -1e100\n",
    "  * max_iter (int): Set the maximum number of nested sampling\n",
    "    iterations performed by MultiNest. If 0, then it is\n",
    "    unlimited and MultiNest will stop using a different\n",
    "    criterion. Default: 0"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "By default, MultiNest will output a set of files with the root 'multinest_run_'. If you want to change the file root you can do so via the multinest_file_root property:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Default: multinest_run_\n",
      "Changed to: run_eggcarton_multinest_\n"
     ]
    }
   ],
   "source": [
    "print(\"Default: {}\".format(MNNS.multinest_file_root))\n",
    "MNNS.multinest_file_root = 'run_eggcarton_multinest_'\n",
    "print(\"Changed to: {}\".format(MNNS.multinest_file_root))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now the MultiNest output files with begin with 'run_eggcarton_multinest'"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Once the instance has been defined, you start the Nested Sampling using the run function:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "  analysing data from run_eggcarton_multinest_.txt\n"
     ]
    }
   ],
   "source": [
    "log_evidence, log_evidence_error = MNNS.run(verbose=False)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      " Ln(Evidence): 235.96857829537618+-0.10992362784591564\n"
     ]
    }
   ],
   "source": [
    "print(\" Ln(Evidence): {}+-{}\".format(log_evidence, log_evidence_error))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In addition to the log_evidence and log_evidence_error estimates returned by the run function, the MultiNestNestedSampling Nested Sampler has the following accesible properties:\n",
    "  * evidence\n",
    "  * evidence_error\n",
    "  \n",
    "The MultiNestNestedSampling Nested Sampler also has the following functions:\n",
    "  * posteriors - Estimates of the posterior marginal probability distributions of each parameter. \n",
    "  * posterior_moments - Estimates of the first four moments of each posterior marginal probability distribution of the parameters.  \n",
    "  * akaike_ic - Estimate of the Akaike Information Criterion.\n",
    "  * bayesian_ic - Estimate of the Bayesian Information Criterion.\n",
    "  * deviance_ic - Estimate of the Deviance Information Criterion.  \n",
    "  * best_fit_likelihood - Get the maximum likelihood parameter vector.\n",
    "  * best_fit_posterior - Get highest posterior probability parameter vector."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 3. PolyChordNestedSampling\n",
    "\n",
    "Gleipnir also provides a Nested Sampler class, PolyChordNestedSampling, that is an interface to the public version of [PolyChord](https://github.com/PolyChord/PolyChordLite); this class is a wrapper on top of pypolychord, which is itself the python wrapper of PolyChord.  Use of PolyChord via Gleipnir requires separate building and installation of pypolychord; see Gleinir's [README](https://github.com/LoLab-VU/Gleipnir#polychord) for the links to those instructions. \n",
    "\n",
    "It is imported from the polychord module:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [],
   "source": [
    "from gleipnir.polychord import PolyChordNestedSampling"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Then the instance is created with a similar input pattern as the NestedSampling, but like the MultiNestNestedSampling object it only requires the sampled parameters, loglikelihood, and population size:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {},
   "outputs": [],
   "source": [
    "PCNS = PolyChordNestedSampling(sampled_parameters, loglikelihood, population_size)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The PolyChordNestedSampling object does not have any extra keyword arguments to worry about."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "By default, PolyChord will output a set of files with the root 'polychord_run_'. If you want to change the file root you can do so via the polychord_file_root property:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Default: polychord_run\n",
      "Changed to: run_eggcarton_polychord_\n"
     ]
    }
   ],
   "source": [
    "print(\"Default: {}\".format(PCNS.polychord_file_root))\n",
    "PCNS.polychord_file_root = 'run_eggcarton_polychord_'\n",
    "print(\"Changed to: {}\".format(PCNS.polychord_file_root))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now the PolyChord output files will begin with 'run_eggcarton_polychord'"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Once the instance has been defined, you start the Nested Sampling using the run function:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {},
   "outputs": [],
   "source": [
    "%%capture\n",
    "# Run. There is no verbose setting for the PCNS object's run function. \n",
    "log_evidence, log_evidence_error = PCNS.run()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      " Ln(Evidence): 236.210847680415+-0.109259064530656\n"
     ]
    }
   ],
   "source": [
    "print(\" Ln(Evidence): {}+-{}\".format(log_evidence, log_evidence_error))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In addition to the log_evidence and log_evidence_error estimates returned by the run function, the PolyChordNestedSampling Nested Sampler has the following accesible properties:\n",
    "  * evidence\n",
    "  * evidence_error\n",
    "  \n",
    "The PolyChordNestedSampling Nested Sampler also has the following functions:\n",
    "  * posteriors - Estimates of the posterior marginal probability distributions of each parameter. \n",
    "  * posterior_moments - Estimates of the first four moments of each posterior marginal probability distribution of the parameters.  \n",
    "  * akaike_ic - Estimate of the Akaike Information Criterion.\n",
    "  * bayesian_ic - Estimate of the Bayesian Information Criterion.\n",
    "  * deviance_ic - Estimate of the Deviance Information Criterion. \n",
    "  * best_fit_likelihood - Get the maximum likelihood parameter vector.\n",
    "  * best_fit_posterior - Get highest posterior probability parameter vector.  "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 4. dyPolyChordNestedSampling\n",
    "\n",
    "Gleipnir also provides a Nested Sampler class, dyPolyChordNestedSampling, that is an interface to [dyPolyChord](https://github.com/ejhigson/dyPolyChord); this class is a wrapper on top of dyPolyChord, which itself uses the python bindings of PolyChord (i.e., pypolychord). Use of dyPolyChord via Gleipnir requires separate building and installation of pypolychord; see Gleinir's [README](https://github.com/LoLab-VU/Gleipnir#polychord) for the links to those instructions. \n",
    "\n",
    "It is imported from the dypolychord module:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 24,
   "metadata": {},
   "outputs": [],
   "source": [
    "from gleipnir.dypolychord import dyPolyChordNestedSampling"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Then the instance is created with the common input pattern:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 25,
   "metadata": {},
   "outputs": [],
   "source": [
    "dyPCNS = dyPolyChordNestedSampling(sampled_parameters, loglikelihood, population_size)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Note that the population_size is used as the dyPolyChord parameter `nlive_const`, which is the total computational budget for allocating population sizes during the dynamic Nested Sampling run. The dyPolyChordNestedSampling object does have some extra keyword arguments that can be passed in to alter the behavior of the dyPolyChord run:\n",
    "   * initial_population_size (int): Set the initial population size for first dyPolyChord exploratory run. Should be < population_size. Default: population_size/2\n",
    "   * dynamic_goal (float or int): Set the dynamic Nested Sampling goal, which is a number in (0, 1) that determines how to allocate computational effort between parameter estimation and evidence calculation; i.e., with a value of 0 dyPolyChord will focus on parameter estimation only, and opposingly with a value of  1 dyPolyChord will focus on evidence calculation only. Default: 0.5\n",
    "\n",
    "\n",
    "  "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "By default, dyPolyChord will output a set of files with the root 'dypolychord_run_' in a 'dypolychord_chains' directory. If you want to change the file root you can do so via the dypolychord_file_root property:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 26,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Default: dypolychord_run\n",
      "Changed to: run_eggcarton_dypolychord_\n"
     ]
    }
   ],
   "source": [
    "print(\"Default: {}\".format(dyPCNS.dypolychord_file_root))\n",
    "dyPCNS.dypolychord_file_root = 'run_eggcarton_dypolychord_'\n",
    "print(\"Changed to: {}\".format(dyPCNS.dypolychord_file_root))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now the dyPolyChord output files will begin with 'run_eggcarton_dypolychord'"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Once the instance has been defined, you start the Nested Sampling using the run function:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 27,
   "metadata": {},
   "outputs": [],
   "source": [
    "%%capture\n",
    "# Run. There is no verbose setting for the dyPCNS object's run function. \n",
    "log_evidence, log_evidence_error = dyPCNS.run()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 28,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      " Ln(Evidence): 235.86150386710474+-None\n"
     ]
    }
   ],
   "source": [
    "print(\" Ln(Evidence): {}+-{}\".format(log_evidence, log_evidence_error))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Note that dyPolyChord doesn't return a rough estimate of the error in the evidence, it instead suggests using a bootstrap analysis of multiple dyPolyChord runs via nestcheck to estimate the error. At the moment, Gleipnir does not implement an interface for this analysis, so users must do that on their own if they want the bootstrapped error estimates.\n",
    "\n",
    "In addition to the log_evidence estimate returned by the run function, the dyPolyChordNestedSampling Nested Sampler has the following accesible properties:\n",
    "  * evidence\n",
    "  * evidence_error --> Just returns None\n",
    "  \n",
    "The dyPolyChordNestedSampling Nested Sampler also has the following functions:\n",
    "  * posteriors - Estimates of the posterior marginal probability distributions of each parameter. \n",
    "  * posterior_moments - Estimates of the first four moments of each posterior marginal probability distribution of the parameters.\n",
    "  * akaike_ic - Estimate of the Akaike Information Criterion.\n",
    "  * bayesian_ic - Estimate of the Bayesian Information Criterion.\n",
    "  * deviance_ic - Estimate of the Deviance Information Criterion. \n",
    "  * best_fit_likelihood - Get the maximum likelihood parameter vector.\n",
    "  * best_fit_posterior - Get highest posterior probability parameter vector.  "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 5. DNest4NestedSampling\n",
    "\n",
    "Gleipnir also provides a Nested Sampler class, DNest4NestedSampling, that is an interface to [DNest4](https://github.com/eggplantbren/DNest4); this class is a wrapper on top of dnest4, which is itself the python wrapper of DNest4. Use of DNest4 via Gleipnir requires separate building and installation of DNest4 and its Python wrapper; see Gleinir's [README](https://github.com/LoLab-VU/Gleipnir#dnest4) for the links to those instructions. \n",
    "\n",
    "It is imported from the dnest4 module:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 29,
   "metadata": {},
   "outputs": [],
   "source": [
    "from gleipnir.dnest4 import DNest4NestedSampling"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Then the instance is created with a similar input pattern as the NestedSampling, but only requires the sampled parameters, loglikelihood, and population size:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 30,
   "metadata": {},
   "outputs": [],
   "source": [
    "DNS = DNest4NestedSampling(sampled_parameters, loglikelihood, population_size)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "However, the DNest4NestedSampling object does have some extra keyword arguments that can be passed in to alter the behavior of the DNest4 run:\n",
    "  * n_diffusive_levels (int): Set the maximum number of\n",
    "      diffusive likelihood levels for DNest4 to use.\n",
    "      Default: 20\n",
    "  * dnest4_backend (str : \"memory\" or \"csv\"): Set which\n",
    "      DNest4 backend to use. \"memory\" means outputs are \n",
    "      just kept in memory. \"csv\" means outputs are \n",
    "      written to disk in files with a csv format. Default: \"memory\"\n",
    "  * num_steps (int): The number of nested iterations to run. If None, \n",
    "      will run forever. \n",
    "      Default: 1000\n",
    "  * new_level_interval (int): The number of nested iterations to run before\n",
    "     creating a new diffusive likelihood level. Default: 10000\n",
    "  * lam (float): Set the backtracking scale length. Default: 5.0\n",
    "  * beta (float): Set the strength of effect to force the histogram\n",
    "      to equal bin counts. Default: 100"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's re-initialize the sampler and set the num_steps variable (the default value is 1000)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 31,
   "metadata": {},
   "outputs": [],
   "source": [
    "DNS = DNest4NestedSampling(sampled_parameters, loglikelihood, population_size, num_steps=500)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "If the dnest4_backend keyword argument is set to \"csv\" (the default is \"memory\", which has no output files) then the DNest4 will output a set of files with the root 'dnest4_run_'. If you want to change the file root you can do so via the dnest4_file_root property:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 32,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Default: ./dnest4_run_\n",
      "Changed to: run_eggcarton_dnest4_\n"
     ]
    }
   ],
   "source": [
    "print(\"Default: {}\".format(DNS.dnest4_file_root))\n",
    "DNS.dnest4_file_root = 'run_eggcarton_dnest4_'\n",
    "print(\"Changed to: {}\".format(DNS.dnest4_file_root))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now any DNest4 output files would begin with 'run_eggcarton_dnest4'"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Once the instance has been defined, you start the Nested Sampling using the run function:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 33,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "500\n",
      "500\n"
     ]
    }
   ],
   "source": [
    "log_evidence, log_evidence_error = DNS.run(verbose=False)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 34,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      " Ln(Evidence): 235.8798566270103+-0.1106471957938031\n"
     ]
    }
   ],
   "source": [
    "print(\" Ln(Evidence): {}+-{}\".format(log_evidence, log_evidence_error))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In addition to the log_evidence and log_evidence_error estimates returned by the run function, the DNest4NestedSampling Nested Sampler has the following accesible properties:\n",
    "  * evidence\n",
    "  * evidence_error\n",
    "  * information\n",
    "  \n",
    "The DNest4NestedSampling Nested Sampler also has the following functions:\n",
    "  * posteriors - Estimates of the posterior marginal probability distributions of each parameter. \n",
    "  * posterior_moments - Estimates of the first four moments of each posterior marginal probability distribution of the parameters.\n",
    "  * akaike_ic - Estimate of the Akaike Information Criterion.\n",
    "  * bayesian_ic - Estimate of the Bayesian Information Criterion.\n",
    "  * deviance_ic - Estimate of the Deviance Information Criterion.  \n",
    "  * best_fit_likelihood - Get the maximum likelihood parameter vector.\n",
    "  * best_fit_posterior - Get highest posterior probability parameter vector.  "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 6. NestleNestedSampling\n",
    "\n",
    "Gleipnir also provides a Nested Sampler class, NestleNestedSampling, that is an interface to [Nestle](https://github.com/kbarbary/nestle); this class is a wrapper on top of nestle, which is the Nestle python package.  \n",
    "\n",
    "It is imported from the gleipnir.nestle module:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 35,
   "metadata": {},
   "outputs": [],
   "source": [
    "from gleipnir.nestle import NestleNestedSampling"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Then the instance is created with a similar input pattern as the NestedSampling with the sampled parameters, loglikelihood, and population size:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 36,
   "metadata": {},
   "outputs": [],
   "source": [
    "NNS = NestleNestedSampling(sampled_parameters, loglikelihood, population_size)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "However, the NestleNestedSampling object does have some extra keyword arguments that can be passed in to alter the behavior of the Nestle run:\n",
    "  * method (str): Method to use to generate new sample points.\n",
    "      Options are 'classic' (MCMC), 'single' \n",
    "      (single-ellipsoidal), and 'multi' (multi-ellipsoidal).   \n",
    "      Default: 'single'\n",
    "  * update_interval (int): Only update the new point generator\n",
    "      every this many likelihood calls. \n",
    "      Default: round(0.6*npoints)\n",
    "  * maxiter (int): Set a maximum number of nested iterations that\n",
    "      are allowed to be run before terminating.\n",
    "      Default: None -> no limit\n",
    "  * maxcall (int): Set a maximum number of likelihood evaluations that\n",
    "      are allowed to be run before terminating the nested sampling\n",
    "      run. Default: None -> no limit\n",
    "  * dlogz (float): If specified, set a termination threshold for\n",
    "      the contribution of the remaining prior volume to the\n",
    "      log-evidence. Terminates when the following condition is\n",
    "      true: log(Z + Z_remain) - log(Z) < dlogz\n",
    "      to output files. True is required for additional\n",
    "      analysis. Default: None -> If neither this nor\n",
    "      decline_factor are set, then dlogz will be set to\n",
    "      dlogz=0.5.\n",
    "  * decline_factor (float): If specified, set a termination\n",
    "      threshold based on the weight (prior mass times likelihood)\n",
    "      of the most recent dead points. Terminates if the weight\n",
    "      has been declining for decline_factor * nsamples iterations.\n",
    "      Nestle documentation suggests a value of 1.0 works well.\n",
    "      Mutually exclusive from dlogz. Default: None -> If neither\n",
    "      this nor dlogz are set, then dlogz will be set to dlogz=0.5.    "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 38,
   "metadata": {},
   "outputs": [],
   "source": [
    "NNS = NestleNestedSampling(sampled_parameters, loglikelihood, population_size, num_steps=500)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Once the instance has been defined, you start the Nested Sampling using the run function:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 39,
   "metadata": {},
   "outputs": [],
   "source": [
    "log_evidence, log_evidence_error = NNS.run(verbose=False)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 40,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      " Ln(Evidence): 235.81611128234087+-0.1111767110157801\n"
     ]
    }
   ],
   "source": [
    "print(\" Ln(Evidence): {}+-{}\".format(log_evidence, log_evidence_error))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In addition to the log_evidence and log_evidence_error estimates returned by the run function, the NestleNestedSampling Nested Sampler has the following accesible properties:\n",
    "  * evidence\n",
    "  * evidence_error\n",
    "  * information\n",
    "  \n",
    "The NestleNestedSampling Nested Sampler also has the following functions:\n",
    "  * posteriors - Estimates of the posterior marginal probability distributions of each parameter. \n",
    "  * posterior_moments - Estimates of the first four moments of each posterior marginal probability distributions of the parameters.\n",
    "  * akaike_ic - Estimate of the Akaike Information Criterion.\n",
    "  * bayesian_ic - Estimate of the Bayesian Information Criterion.\n",
    "  * deviance_ic - Estimate of the Deviance Information Criterion.  \n",
    "  * best_fit_likelihood - Get the maximum likelihood parameter vector.\n",
    "  * best_fit_posterior - Get highest posterior probability parameter vector.  "
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.6.8"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}