{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Prior definition\n", "\n", "In this notebook we demonstrate how to include prior knowledge in a parameter inference problem, in particular how to define (log-)priors for parameters. If you want to maximize your posterior distribution, you need to define\n", "\n", "* A (negative log-)likelihood\n", "* A (log-)prior\n", "\n", "The posterior is then built as an `AggregatedObjective`. If you import a problem via `PEtab` and the priors are contained in the parameter table, the definition of priors is done automatically.\n", "\n", "**CAUTION**: The user needs to specify the **negative** _log-likelihood_, while the _log-prior_ is internally multiplied by -1." ] }, { "cell_type": "raw", "metadata": { "raw_mimetype": "text/restructuredtext", "tags": [] }, "source": [ "In this notebook:\n", "\n", "* :class:`pypesto.objective.AggregatedObjective `\n", "* :func:`pypesto.objective.get_parameter_prior_dict `" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "# install if not done yet\n", "# %pip install pypesto --quiet" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "import numpy as np\n", "\n", "import pypesto" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Example: Rosenbrock Banana\n", "\n", "We will use the Rosenbrock Banana\n", "\n", "\\begin{align}\n", "f(x, \\theta) = \\sum_{i=1}^{N} \\underbrace{100 \\cdot(x_{i}-x_{i-1}^2)^2}_{\\text{\"negative log-likelihood\"}} + \\underbrace{(x_{i-1}-1)^2}_{\\text{\"Gaussian log-prior\"}}\n", "\\end{align}\n", "\n", "as an example. Here, we interpret the first term as the _negative log-likelihood_ and the second term as Gaussian _log-prior_ with mean $1$ and standard deviation $1/\\sqrt{2}$.\n", "\n", "Note that the second term is only equivalent to the negative log-distribution of a Gaussian up to a constant." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Define the negative log-likelihood" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "n_x = 5\n", "\n", "\n", "def rosenbrock_part_1(x):\n", " \"\"\"\n", " Calculate obj. fct + gradient of the \"likelihood\" part.\n", " \"\"\"\n", " obj = sum(100.0 * (x[1:] - x[:-1] ** 2.0) ** 2.0)\n", "\n", " grad = np.zeros_like(x)\n", " grad[:-1] += -400 * (x[1:] - x[:-1] ** 2.0) * x[:-1]\n", " grad[1:] += 200 * (x[1:] - x[:-1] ** 2.0)\n", "\n", " return obj, grad\n", "\n", "\n", "neg_log_likelihood = pypesto.Objective(fun=rosenbrock_part_1, grad=True)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Define the log-prior\n", "\n", "A prior on an individual parameter is defined in a `prior_dict`, which contains the following key-value pairs:\n", "\n", "* `index`: Index of the parameter\n", "* `density_fun`: (Log-)posterior. (Scalar function!)\n", "* `density_dx`: d/dx (Log-)posterior (optional)\n", "* `density_ddx`: d^2/dx^2 (Log-)posterior (optional)\n", "\n", "A `prior_dict` can be either obtained by `get_parameter_prior_dict` for several common priors, or defined by the user." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "from pypesto.objective.priors import get_parameter_prior_dict\n", "\n", "# create a list of prior dicts...\n", "prior_list = []\n", "mean = 1\n", "std_dev = 1 / np.sqrt(2)\n", "\n", "for i in range(n_x - 1):\n", " prior_list.append(get_parameter_prior_dict(i, \"normal\", [mean, std_dev]))\n", "\n", "# create the prior\n", "neg_log_prior = pypesto.objective.NegLogParameterPriors(prior_list)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Define the negative log-posterior and the problem\n", "\n", "The negative log-posterior is defined as an `AggregatedObjective`. Since optimization/visualization is not the main focus of this notebook, the reader is referred to other examples for a more in-depth presentation of these. Here, basic visualisation is provided." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "neg_log_posterior = pypesto.objective.AggregatedObjective(\n", " [neg_log_likelihood, neg_log_prior]\n", ")\n", "\n", "lb = -5 * np.ones((n_x, 1))\n", "ub = 5 * np.ones((n_x, 1))\n", "\n", "problem = pypesto.Problem(\n", " objective=neg_log_posterior,\n", " lb=lb,\n", " ub=ub,\n", ")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Optimize" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "import pypesto.optimize as optimize\n", "\n", "result = optimize.minimize(problem=problem, n_starts=10, filename=None)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Some basic visualizations" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "import pypesto.visualize as visualize\n", "\n", "visualize.waterfall(result, size=(15, 6))\n", "\n", "# parallel coordinates plot for best 5 fits\n", "visualize.parameters(result, start_indices=5);" ] } ], "metadata": { "kernelspec": { "display_name": "Python 3 (ipykernel)", "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.3" }, "vscode": { "interpreter": { "hash": "44a9cdcbdccbf05a880e90d2e6fe72470baab4d1b82472d890be0596ed887a6b" } } }, "nbformat": 4, "nbformat_minor": 4 }