diff --git a/contributors.txt b/contributors.txt index f550fe1b..f392ed49 100644 --- a/contributors.txt +++ b/contributors.txt @@ -23,3 +23,5 @@ Lea Hirsch -- advising RT on development of RV/astrometry joint fitting code Eric Nielsen -- came up with idea for OFTI algorithm; advised SB on development of OFTI code that sampler.py inherits from Roberto Tejada -- developing RV/astrometry joint fitting code + +Aaditya Thokal -- wrote MultiNest tutorial and installation documentation \ No newline at end of file diff --git a/docs/installation.rst b/docs/installation.rst index e7ab5c92..0d4010c8 100644 --- a/docs/installation.rst +++ b/docs/installation.rst @@ -26,6 +26,42 @@ Next, install ``orbitize``: $ pip install orbitize +MultiNest (Optional) +++++++++++++++++++++ + +The MultiNest sampler requires both the ``PyMultiNest`` Python package and the +native MultiNest library. + +First, install ``PyMultiNest``: + +.. code-block:: bash + + $ pip install pymultinest + +Next, build and install the native MultiNest library by following the official +installation instructions: + +https://johannesbuchner.github.io/PyMultiNest/install.html + +After installation, ensure that the directory containing the MultiNest shared +library is included in your library search path. + +On Linux: + +.. code-block:: bash + + $ export LD_LIBRARY_PATH=/path/to/MultiNest/lib:$LD_LIBRARY_PATH + +On macOS: + +.. code-block:: bash + + $ export DYLD_LIBRARY_PATH=/path/to/MultiNest/lib:$DYLD_LIBRARY_PATH + +For quick testing, approximately **200 live points** are generally sufficient. +For publication-quality analyses, we recommend using **1000 or more live +points**, depending on the complexity of the orbital model. + We recommend installing and running ``orbitize`` in a ``conda`` virtual environment. Install ``anaconda`` or ``miniconda`` `here `_, then see instructions diff --git a/docs/tutorials.rst b/docs/tutorials.rst index 46798550..31116237 100644 --- a/docs/tutorials.rst +++ b/docs/tutorials.rst @@ -57,6 +57,7 @@ us if you are still confused). tutorials/Hipparcos_IAD.ipynb tutorials/HGCA_tutorial.ipynb tutorials/abs_astrometry.ipynb + tutorials/multinest_tutorial.ipynb diff --git a/docs/tutorials/multinest_tutorial.ipynb b/docs/tutorials/multinest_tutorial.ipynb new file mode 100644 index 00000000..f652d4fe --- /dev/null +++ b/docs/tutorials/multinest_tutorial.ipynb @@ -0,0 +1,328 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "id": "2dc6f93f", + "metadata": {}, + "source": [ + "# MultiNest Sampler Tutorial\n", + "\n", + "Based on the MultiNest implementation by Tomas Stolker (2024).\n", + "\n", + "MultiNest is a Bayesian nested sampling algorithm that efficiently explores complex posterior distributions while simultaneously estimating the Bayesian evidence. In contrast to traditional Markov Chain Monte Carlo (MCMC) methods, MultiNest is particularly well suited for multimodal or highly degenerate parameter spaces and enables Bayesian model comparison through evidence estimation.\n", + "\n", + "The `orbitize!` implementation interfaces with the MultiNest library through `PyMultiNest`. Because the native MultiNest library must be compiled separately, users should first complete the MultiNest installation described in the `orbitize!` installation guide before running this tutorial.\n", + "\n", + "This tutorial uses the bundled GJ504.csv example dataset distributed with orbitize!, which contains relative astrometric measurements of GJ 504 b.\n", + "\n", + "- loading astrometric observations\n", + "- constructing an `orbitize.System`\n", + "- configuring the MultiNest sampler\n", + "- selecting an appropriate number of live points\n", + "- executing the sampler\n", + "- saving the posterior samples" + ] + }, + { + "cell_type": "markdown", + "id": "0a5e3672", + "metadata": {}, + "source": [ + "## Prerequisites\n", + "\n", + "Before running this tutorial, ensure that:\n", + "\n", + "1. `orbitize!` is installed.\n", + "2. MultiNest and `PyMultiNest` are installed (see the installation guide).\n", + "3. `astrometry_hd135344ab.csv` is located in the same directory as this notebook." + ] + }, + { + "cell_type": "markdown", + "id": "234df74a", + "metadata": {}, + "source": [ + "## Load the Astrometric Data\n", + "\n", + "This tutorial uses relative astrometric measurements of **GJ 504 b**.\n", + "\n", + "The observations are stored in `astrometry_hd135344ab.csv` and consist of measured right ascension and declination offsets together with their uncertainties and instrument identifiers.\n", + "\n", + "These measurements will be used to construct an `orbitize.System`, which defines the orbital model and likelihood function explored by MultiNest.\n", + "\n", + "Ensure that `astrometry_hd135344ab.csv` is located in the same directory as this notebook before continuing.\n", + "\n", + "The data file is included in the tutorials/ folder of the orbitize! repository." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "13abbe86", + "metadata": {}, + "outputs": [], + "source": [ + "from orbitize import read_input\n", + "from orbitize import sampler\n", + "from orbitize import system" + ] + }, + { + "cell_type": "code", + "execution_count": 20, + "id": "1b26ffb0", + "metadata": {}, + "outputs": [], + "source": [ + "import os\n", + "import orbitize\n", + "\n", + "data_file = os.path.join(orbitize.DATADIR, \"GJ504.csv\")\n", + "data_table = read_input.read_file(data_file)" + ] + }, + { + "cell_type": "markdown", + "id": "60d4332a", + "metadata": {}, + "source": [ + "The loaded table contains the observation epochs, measured relative astrometry, associated uncertainties, and instrument identifiers.\n", + "\n", + "These measurements define the likelihood function that will be explored by the MultiNest sampler." + ] + }, + { + "cell_type": "markdown", + "id": "8f403e9a", + "metadata": {}, + "source": [ + "## Construct the Orbit System\n", + "\n", + "Next we create an `orbitize.System` object.\n", + "\n", + "This object stores the observational data together with the stellar properties and orbital model that will be sampled by MultiNest." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "082fe819", + "metadata": {}, + "outputs": [], + "source": [ + "orb_system = system.System(\n", + " num_secondary_bodies=1,\n", + " data_table=data_table,\n", + " stellar_or_system_mass=2.2,\n", + " plx=7.41,\n", + " mass_err=0.1,\n", + " plx_err=0.04,\n", + " restrict_angle_ranges=True,\n", + " tau_ref_epoch=58849,\n", + " fit_secondary_mass=False,\n", + " hipparcos_IAD=None,\n", + " gaia=None,\n", + " fitting_basis=\"Standard\",\n", + " use_rebound=False,\n", + ")" + ] + }, + { + "cell_type": "markdown", + "id": "9e82fb21", + "metadata": {}, + "source": [ + "## Configure the MultiNest Sampler\n", + "\n", + "The principal MultiNest hyperparameter is the number of **live points**.\n", + "\n", + "Increasing the number of live points generally improves posterior accuracy and Bayesian evidence estimation, but also increases computational cost.\n", + "\n", + "200 live points are generally sufficient for testing and tutorial examples.\n", + "\n", + "1000 or more live points are recommended for publication-quality analyses, following the recommendation of Tomas Stolker." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "61fbfc69", + "metadata": {}, + "outputs": [], + "source": [ + "n_live_points = 200\n", + "\n", + "orb_sampler = sampler.MultiNest(orb_system)" + ] + }, + { + "cell_type": "markdown", + "id": "3cd4f936", + "metadata": {}, + "source": [ + "> **Note**\n", + ">\n", + "> MultiNest requires both the native MultiNest library and the\n", + "> `PyMultiNest` Python package.\n", + ">\n", + "> On Linux, ensure that the directory containing `libmultinest.so`\n", + "> is included in `LD_LIBRARY_PATH`.\n", + ">\n", + "> On macOS, use `DYLD_LIBRARY_PATH`.\n", + ">\n", + "> If these libraries are unavailable, the sampling step below will fail.\n", + "> Please complete the installation instructions before continuing." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "131662a8", + "metadata": {}, + "outputs": [], + "source": [ + "orb_sampler.run_sampler(\n", + " n_live_points=n_live_points,\n", + " output_basename=\"./multinest/\",\n", + " hdf5_file=\"orbitize.hdf5\",\n", + " multinest_kwargs={\n", + " \"resume\": False\n", + " },\n", + ")" + ] + }, + { + "cell_type": "markdown", + "id": "2bed2d25", + "metadata": {}, + "source": [ + "## Output\n", + "\n", + "After completion, MultiNest writes the posterior samples to the specified HDF5 file together with the standard MultiNest output directory.\n", + "\n", + "These files can later be analysed using the standard `orbitize!` results utilities without rerunning the sampler." + ] + }, + { + "cell_type": "markdown", + "id": "41776f5c", + "metadata": {}, + "source": [ + "## Tips for Choosing the Number of Live Points\n", + "\n", + "The number of live points controls the balance between computational cost and sampling quality.\n", + "\n", + "| Live Points | Recommended Usage |\n", + "|-------------|-------------------|\n", + "| 200 | Quick tutorial / testing |\n", + "| 500 | General scientific analyses |\n", + "| 1000+ | Publication-quality posterior estimation |\n", + "\n", + "Increasing the number of live points generally improves posterior estimation and Bayesian evidence calculation, although runtime increases approximately linearly." + ] + }, + { + "cell_type": "markdown", + "id": "cabd6581", + "metadata": {}, + "source": [ + "## Tips for Faster Convergence\n", + "\n", + "Efficient nested sampling depends strongly on the choice of prior distributions and sampler configuration.\n", + "\n", + "The following recommendations have consistently improved convergence when fitting orbital models with `orbitize!`." + ] + }, + { + "cell_type": "markdown", + "id": "09e3ad22", + "metadata": {}, + "source": [ + "## Saving Results\n", + "\n", + "The posterior samples are automatically written to the HDF5 file specified with the `hdf5_file` argument in `run_sampler()`.\n", + "\n", + "This is the recommended approach because the sampler safely handles writing the results, including MPI-enabled executions where writing the file manually after sampling can cause issues.\n", + "\n", + "The saved HDF5 file can later be reloaded with `orbitize!` without rerunning the sampler." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "9968c35b", + "metadata": {}, + "outputs": [], + "source": [ + "# Uncomment after MultiNest has been installed.\n", + "\n", + "# results = orb_sampler.run_sampler(\n", + "# n_live_points=n_live_points,\n", + "# output_basename=\"./multinest/\",\n", + "# hdf5_file=\"orbitize.hdf5\",\n", + "# multinest_kwargs={\"resume\": False},\n", + "# )" + ] + }, + { + "cell_type": "markdown", + "id": "6eb2308f", + "metadata": {}, + "source": [ + "The following command runs the MultiNest sampler. It is commented out because documentation builds typically do not include the MultiNest shared library." + ] + }, + { + "cell_type": "markdown", + "id": "c3917a68", + "metadata": {}, + "source": [ + "## When should I use MultiNest?\n", + "\n", + "MultiNest is particularly useful when Bayesian evidence is required for model comparison or when the posterior distribution is expected to contain multiple modes.\n", + "\n", + "For simpler orbit-fitting problems where evidence estimation is unnecessary, OFTI or MCMC may provide a more computationally efficient alternative." + ] + }, + { + "cell_type": "markdown", + "id": "b511a3df", + "metadata": {}, + "source": [ + "## Summary\n", + "\n", + "In this tutorial we:\n", + "\n", + "- loaded astrometric observations of GJ 504 b\n", + "- constructed an `orbitize.System`\n", + "- configured the MultiNest sampler\n", + "- discussed the role of the number of live points\n", + "- executed the MultiNest sampler\n", + "- saved the resulting posterior samples\n", + "\n", + "MultiNest provides an efficient Bayesian nested sampling framework capable of simultaneously estimating posterior distributions and Bayesian evidence, making it particularly valuable for orbit-fitting problems with complex or multimodal parameter spaces." + ] + } + ], + "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.12.1" + } + }, + "nbformat": 4, + "nbformat_minor": 5 +}