Skip to content

Commit a272255

Browse files
dweindldweindl
committed
Refactor priors
Introduce `Distribution` classes for handling prior distributions and sampling from them. Later on this can be extended to noise models for measurements.
1 parent 4a551a7 commit a272255

File tree

8 files changed

+725
-108
lines changed

8 files changed

+725
-108
lines changed

‎doc/example.rst‎

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -10,6 +10,7 @@ The following examples should help to get a better idea of how to use the PEtab
1010

1111
example/example_petablint.ipynb
1212
example/example_visualization.ipynb
13+
example/distributions.ipynb
1314

1415
Examples of systems biology parameter estimation problems specified in PEtab
1516
can be found in the `systems biology benchmark model collection <https://github.com/Benchmarking-Initiative/Benchmark-Models-PEtab>`_.

‎doc/example/distributions.ipynb‎

Lines changed: 197 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,197 @@
1+
{
2+
"cells": [
3+
{
4+
"metadata": {},
5+
"cell_type": "markdown",
6+
"source": [
7+
"# Prior distributions in PEtab\n",
8+
"\n",
9+
"This notebook gives a brief overview of the prior distributions in PEtab and how they are represented in the PEtab library.\n",
10+
"\n",
11+
"Prior distributions are used to specify the prior knowledge about the parameters.\n",
12+
"Parameter priors are specified in the parameter table. A prior is defined by its type and its parameters.\n",
13+
"Each prior type has a specific set of parameters. For example, the normal distribution has two parameters: the mean and the standard deviation.\n",
14+
"\n",
15+
"There are two types of priors in PEtab - objective priors and initialization priors:\n",
16+
"\n",
17+
"* *Objective priors* are used to specify the prior knowledge about the parameters that are to be estimated. They will enter the objective function of the optimization problem. They are specified in the `objectivePriorType` and `objectivePriorParameters` columns of the parameter table.\n",
18+
"* *Initialization priors* can be used as a hint for the optimization algorithm. They will not enter the objective function. They are specified in the `initializationPriorType` and `initializationPriorParameters` columns of the parameter table.\n",
19+
"\n",
20+
"\n"
21+
],
22+
"id": "372289411a2aa7b3"
23+
},
24+
{
25+
"metadata": {
26+
"collapsed": true
27+
},
28+
"cell_type": "code",
29+
"source": [
30+
"import matplotlib.pyplot as plt\n",
31+
"import numpy as np\n",
32+
"import seaborn as sns\n",
33+
"\n",
34+
"from petab.v1.C import *\n",
35+
"from petab.v1.distributions import *\n",
36+
"\n",
37+
"sns.set_style(None)\n",
38+
"\n",
39+
"\n",
40+
"def plot(distr: Distribution, ax=None):\n",
41+
" \"\"\"Visualize a distribution.\"\"\"\n",
42+
" if ax is None:\n",
43+
" fig, ax = plt.subplots()\n",
44+
"\n",
45+
" sample = distr.sample(10000)\n",
46+
"\n",
47+
" # pdf\n",
48+
" xmin = min(sample.min(), distr.lb_scaled if distr.bounds is not None else sample.min())\n",
49+
" xmax = max(sample.max(), distr.ub_scaled if distr.bounds is not None else sample.max())\n",
50+
" x = np.linspace(xmin, xmax, 500)\n",
51+
" y = distr.pdf(x)\n",
52+
" ax.plot(x, y, color='red', label='pdf')\n",
53+
"\n",
54+
" sns.histplot(sample, stat='density', ax=ax, label=\"sample\")\n",
55+
"\n",
56+
" # bounds\n",
57+
" if distr.bounds is not None:\n",
58+
" for bound in (distr.lb_scaled, distr.ub_scaled):\n",
59+
" if bound is not None and np.isfinite(bound):\n",
60+
" ax.axvline(bound, color='black', linestyle='--', label='bound')\n",
61+
"\n",
62+
" ax.set_title(str(distr))\n",
63+
" ax.set_xlabel('Parameter value on the parameter scale')\n",
64+
" ax.grid(False)\n",
65+
" handles, labels = ax.get_legend_handles_labels()\n",
66+
" unique_labels = dict(zip(labels, handles))\n",
67+
" ax.legend(unique_labels.values(), unique_labels.keys())\n",
68+
" plt.show()"
69+
],
70+
"id": "initial_id",
71+
"outputs": [],
72+
"execution_count": null
73+
},
74+
{
75+
"metadata": {},
76+
"cell_type": "markdown",
77+
"source": "The basic distributions are the uniform, normal, Laplace, log-normal, and log-laplace distributions:\n",
78+
"id": "db36a4a93622ccb8"
79+
},
80+
{
81+
"metadata": {},
82+
"cell_type": "code",
83+
"source": [
84+
"plot(Uniform(0, 1))\n",
85+
"plot(Normal(0, 1))\n",
86+
"plot(Laplace(0, 1))\n",
87+
"plot(LogNormal(0, 1))\n",
88+
"plot(LogLaplace(1, 0.5))"
89+
],
90+
"id": "4f09e50a3db06d9f",
91+
"outputs": [],
92+
"execution_count": null
93+
},
94+
{
95+
"metadata": {},
96+
"cell_type": "markdown",
97+
"source": "If a parameter scale is specified, the sample is transformed accordingly (but not the distribution parameters):\n",
98+
"id": "dab4b2d1e0f312d8"
99+
},
100+
{
101+
"metadata": {},
102+
"cell_type": "code",
103+
"source": [
104+
"plot(Normal(10, 2, transformation=LIN))\n",
105+
"plot(Normal(10, 2, transformation=LOG))\n",
106+
"# Note that the log-normal is different from the log-transformed normal distribution:\n",
107+
"plot(LogNormal(10, 2, transformation=LIN))"
108+
],
109+
"id": "f6192c226f179ef9",
110+
"outputs": [],
111+
"execution_count": null
112+
},
113+
{
114+
"metadata": {},
115+
"cell_type": "markdown",
116+
"source": "Prior distributions can also be defined on the parameter scale by using the types `parameterScaleUniform`, `parameterScaleNormal` or `parameterScaleLaplace`. In these cases, a sample from the given distribution is used directly, without applying any transformation according to `parameterScale` (this implies, that for `parameterScale=lin`, there is no difference between `parameterScaleUniform` and `uniform`):",
117+
"id": "263c9fd31156a4d5"
118+
},
119+
{
120+
"metadata": {},
121+
"cell_type": "code",
122+
"source": [
123+
"plot(Uniform(0, 1, transformation=LOG10))\n",
124+
"plot(ParameterScaleUniform(0, 1, transformation=LOG10))\n",
125+
"\n",
126+
"plot(Uniform(0, 1, transformation=LIN))\n",
127+
"plot(ParameterScaleUniform(0, 1, transformation=LIN))\n"
128+
],
129+
"id": "5ca940bc24312fc6",
130+
"outputs": [],
131+
"execution_count": null
132+
},
133+
{
134+
"metadata": {},
135+
"cell_type": "markdown",
136+
"source": "To prevent the sampled parameters from exceeding the bounds, the sampled parameters are clipped to the bounds. The bounds are defined in the parameter table. Note that the current implementation does not support sampling from a truncated distribution. Instead, the samples are clipped to the bounds. This may introduce unwanted bias, and thus, should only be used with caution (i.e., the bounds should be chosen wide enough):",
137+
"id": "b1a8b17d765db826"
138+
},
139+
{
140+
"metadata": {},
141+
"cell_type": "code",
142+
"source": [
143+
"plot(Normal(0, 1, bounds=(-4, 4))) # negligible clipping-bias at 4 sigma\n",
144+
"plot(Uniform(0, 1, bounds=(0.1, 0.9))) # significant clipping-bias"
145+
],
146+
"id": "4ac42b1eed759bdd",
147+
"outputs": [],
148+
"execution_count": null
149+
},
150+
{
151+
"metadata": {},
152+
"cell_type": "markdown",
153+
"source": "Further distribution examples:",
154+
"id": "45ffce1341483f24"
155+
},
156+
{
157+
"metadata": {},
158+
"cell_type": "code",
159+
"source": [
160+
"plot(Normal(10, 1, bounds=(6, 14), transformation=\"log10\"))\n",
161+
"plot(ParameterScaleNormal(10, 1, bounds=(10**6, 10**14), transformation=\"log10\"))\n"
162+
],
163+
"id": "581e1ac431860419",
164+
"outputs": [],
165+
"execution_count": null
166+
},
167+
{
168+
"metadata": {},
169+
"cell_type": "code",
170+
"source": "",
171+
"id": "802a64be56a6c94f",
172+
"outputs": [],
173+
"execution_count": null
174+
}
175+
],
176+
"metadata": {
177+
"kernelspec": {
178+
"display_name": "Python 3",
179+
"language": "python",
180+
"name": "python3"
181+
},
182+
"language_info": {
183+
"codemirror_mode": {
184+
"name": "ipython",
185+
"version": 2
186+
},
187+
"file_extension": ".py",
188+
"mimetype": "text/x-python",
189+
"name": "python",
190+
"nbconvert_exporter": "python",
191+
"pygments_lexer": "ipython2",
192+
"version": "2.7.6"
193+
}
194+
},
195+
"nbformat": 4,
196+
"nbformat_minor": 5
197+
}

‎doc/modules.rst‎

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -14,6 +14,7 @@ API Reference
1414
petab.v1.composite_problem
1515
petab.v1.conditions
1616
petab.v1.core
17+
petab.v1.distributions
1718
petab.v1.lint
1819
petab.v1.measurements
1920
petab.v1.models

0 commit comments

Comments
 (0)