# Playing it safe: information constrains collective betting strategies

Philipp Fleig  Department of Cellular Biophysics, Max Planck Institute for
Medical Research, Jahnstraße 29, 69120 Heidelberg, Germany  Department of
Physics & Astronomy, University of Pennsylvania, Philadelphia, PA 19104, USA
Vijay Balasubramanian  Department of Physics & Astronomy, University of
Pennsylvania, Philadelphia, PA 19104, USA  Santa Fe Institute, Santa Fe, NM
87501, USA

###### Abstract

Every interaction of a living organism with its environment involves the
placement of a bet. Armed with partial knowledge about a stochastic world, the
organism must decide its next step or near-term strategy, an act that
implicitly or explicitly involves the assumption of a model of the world.
Better information about environmental statistics can improve the bet quality,
but in practice resources for information gathering are always limited. We
argue that theories of optimal inference dictate that “complex” models are
harder to infer with bounded information and lead to larger prediction errors.
Thus, we propose a principle of playing it safe where, given finite
information gathering capacity, biological systems should be biased towards
simpler models of the world, and thereby to less risky betting strategies. In
the framework of Bayesian inference, we show that there is an optimally safe
adaptation strategy determined by the Bayesian prior. We then demonstrate
that, in the context of stochastic phenotypic switching by bacteria,
implementation of our principle of “playing it safe” increases fitness
(population growth rate) of the bacterial collective. We suggest that the
principle applies broadly to problems of adaptation, learning and evolution,
and illuminates the types of environments in which organisms are able to
thrive.

††preprint: APS/123-QED

## Introduction

Risk is an inherent part of life. Whether in pathogen detection [1], phenotype
selection [2, 3], biochemical and evolutionary mechanisms [4, 5], foraging and
exploration-exploitation strategies [6, 7, 8, 9], or reproduction [10, 11],
biological functions are partly shaped by the need to reduce risk. Broadly,
risk arises from uncertain interactions of an organism with the world around
it [12]. First, both the environment and the typical sensory apparatus are
intrinsically noisy. Thus an organism only has probabilistic information about
the state of the world. Second, any finite system can only observe some
aspects of the world, while others, perhaps most, remain hidden from
observation. Thus, an organism inherently functions with partial information,
bounding its capacity for rational decision-making. In other words, living
systems must necessarily play betting games, producing responses, e.g.
expressing a phenotype, detecting an odor, showing an immune response, that
lead to outcomes that will probably be beneficial given the limited and
uncertain information available for decision-making.

Every betting game has an associated risk. To improve the betting strategy, it
is vital to know the odds of the game. An organism can mitigate the risk of
its strategy, and minimize the likelihood of wrong bets, by gathering
information to assess the odds and adapting to the learned statistics of the
environment. Such learning and adaptation is particularly important if rare,
but potentially catastrophic, events can occur.

Living systems can adapt in this way to the environment over many timescales,
over the course of evolution, in response to environmental cues, and through
ongoing life-long learning. Given the stochastic nature of the environment
from the perspective of an agent inhabiting it, essentially all adaptation
processes can be viewed as schemes for inferring a probability distribution
[13, 14]. The ability to infer this distribution is limited by physical and
temporal constraints of information gathering. In turn, this limits the
capacity for adaptation, and thus increases risk for the system. There are
several probabilistic inference frameworks that explicitly or implicitly
account for the limitations of information gathering processes guiding model
choice [15, 16, 17, 18, 19, 20, 21, 22]. Here we argue that a key idea arising
from these frameworks is that “complex” models, in a sense defined in
information geometry [23, 22], are difficult to infer robustly from limited
data, and that a small inaccuracy in inference of such models can lead to
large deviations from the optimal strategy.

Because of this, if information-gathering is limited, it can pay to bias
inference towards simple models. Thus, we propose a principle of “playing it
safe”, which dictates that given finite information, biological systems should
effectively bias themselves towards simpler models. We derive this principle
for a large class of probabilistic models using Bayesian probability theory
and information geometry, and show how the bias towards simplicity can be
tuned by the choice of prior. We then illustrate the efficacy of the principle
in the classic example of Kelly betting [24, 25], and demonstrate a biological
realization in the phenomenon of stochastic phenotypic switching by bacteria.
We conclude with a discussion of the implications of our finding for the types
of environments that are “learnable”, and in which living systems can prosper.

## Results

Figure 1: Landscape of model complexity dictates the required inference
precision. Inferring a model close to a boundary of parameter space requires
higher inference precision. (A) Jeffreys prior (black curve) of the Bernoulli
model with event probability $\theta$ measures model complexity and diverges
at the boundaries of the one-dimensional parameter space. Three example models
are indicated by the orange, red, and purple markers. Inset: precision with
which the models $\phi=\theta+\delta\theta$ need to be determined to keep the
Kullback-Leibler divergence, below some value $\epsilon$. (B) Jeffreys prior
for the multinomial model in $d=3$ dimensions. The parameter space is two-
dimensional, showing the typical divergences of probabilistic models on
boundaries of the space.

### Optimising betting success is hard for complex models

Betting games are unpredictable by nature. However, bet placement can be
optimised by knowing the probabilities controlling events in the game, which
often requires estimation of a probability distribution determined by $d$ real
parameters, $P(\mathbf{x};{\bm{\theta}}=\\{\theta_{1},\ldots,\theta_{d}\\})$,
denoted here by ${\bm{\theta}}$. In real-world settings, the information
required to infer parameters is typically limited. For example, we may have
$N$ independent observations $\mathbf{x}_{i}$,
$\mathbf{n}=\\{\mathbf{x}_{1},\ldots,\mathbf{x}_{N}\\}$ of a random event
which can be used for the inference of parameters. However, finite data leads
to statistical fluctuations in the inference process and we need to understand
how these affect betting success.

Success in betting games can be quantified through a loss function. In many
games the loss is a function of the game’s true probability distribution
$P_{\bm{\theta}^{*}}$ and an inferred distribution $Q_{\bm{\phi}}$, adopted by
the bettor. The true distribution can lie outside the manifold parameterized
by ${\bm{\theta}}$ within which the inference is performed, but we will
nevertheless use the notation ${\bm{\theta}}^{*}$ for the data generating
distribution, because we will often consider cases where this distribution
lies within the considered model family. We expect the loss to have a global
minimum when ${\bm{\phi}}={\bm{\theta}}^{*}$ if the truth lies within the
model family under consideration, since exact knowledge of the game’s odds
allows the bettor to place bets optimally and thus minimise the loss. If the
truth lies outside the model family, the loss should be minimized when
${\bm{\phi}}$ comes “as close as possible” in terms of the loss function.
Here, we adopt the Kullback-Leibler (KL) divergence

$\displaystyle\mathcal{L}_{\bm{\theta}^{*}}(\bm{\phi})\equiv
D_{\mathrm{KL}}(P_{\bm{\theta}^{*}}||Q_{\bm{\phi}})=\sum_{\mathbf{x}}P_{\bm{\theta}^{*}}(\mathbf{x})\log\frac{P_{\bm{\theta}^{*}}(\mathbf{x})}{Q_{\bm{\phi}}(\mathbf{x})}$
(1)

to measure the difference between the true probabilities and the adopted
model. The closer the adopted model is to the truth, the smaller the loss. The
KL divergence is the loss function of the well-known Kelly betting setup. An
important property of the KL divergence is that for many natural
parametrizations of distributions it diverges on the boundaries of the
parameter space. Changing the parameter by a small amount in such regions
leads to large changes in the KL divergence.

A way to understand the origin of this behaviour is through model complexity.
Consider the Jeffreys prior

$\displaystyle
w({\bm{\theta}})=\frac{\sqrt{\det{\mathcal{I}({\bm{\theta}})}}}{\int\mathrm{d}^{d}\bm{\psi}\sqrt{\det{\mathcal{I}(\bm{\psi})}}}$
(2)

where $\mathcal{I}({\bm{\theta}})$ is the Fisher information metric (FIM),
obtained by Taylor expansion from the KL divergence between models in the
parametric family [26, 27] (Materials and Methods). In information geometry,
this quantity is understood as the density of distinguishable models, and is
regarded as a local measure of complexity in a parametric model family [27,
22].

In regions where Jeffreys prior is small, i.e., where the model family is
simple, changing a parameter by a small amount does not move the distribution
much in model space [22]. However, Jeffreys prior can become large and even
diverge. In such regions, i.e., where the model family is complex, even models
which are close in parameter space are far from each other in model space. In
many natural parametrizations, divergences occur at the boundaries of
parameter space (see Figure 1 and Eq. (25) in Materials and Methods for the
multinomial distribution, and Appendix, Section A for the Poisson and Gaussian
distributions).

The multinomial distribution in $d=2$ variables is described by a single
parameter, typically the event probability $\theta$. The corresponding
Jeffreys prior in Figure 1A takes the form
$w(\theta)=1/\pi\sqrt{\theta(1-\theta)}$. Equivalently
$w(\theta)\,\delta\theta^{2}$ measures the KL divergence between multinomial
models parametrized by $\theta$ and $\theta+\delta\theta$. We see that the
density of distinguishable models has a minimum at $\theta=1/2$ and diverges
at the boundaries in these particular coordinates. The inset in Figure 1A
shows that to keep the KL divergence $D_{\mathrm{KL}}(P_{\theta}||P_{\phi})$
between two nearby multinomial models below some value $\epsilon$, the
parameter $\phi\equiv\theta+\delta\theta$ must lie closer to $\theta$, when
$\theta$ itself lies closer to the boundary of the one-dimensional parameter
space at $0$. By symmetry, the same holds near the other boundary at $1$.
Thus, we see that in this case when the true model is close to a boundary of
parameter space, the parameter $\phi$ must be inferred with higher precision
than when the true model lies near a value around $1/2$, in order to achieve
the same degree of statistical proximity.

It is hard to optimise the KL divergence in regions of parameter space that
are complex in the sense described above. For example, we may decide to rely
on the maximum-likelihood estimate of parameters
$\hat{\bm{\theta}}_{\mathrm{ML}}$ and define the expected loss as

$\displaystyle\big{<}\mathcal{L}_{\bm{\theta}^{*}}(\hat{\bm{\theta}}_{\mathrm{ML}})\big{>}_{\hat{\bm{\theta}}_{\mathrm{ML}}}\equiv\int\mathrm{d}^{d}\hat{\bm{\theta}}_{\mathrm{ML}}\,p(\hat{\bm{\theta}}_{\mathrm{ML}})\,\mathcal{L}_{\bm{\theta}^{*}}(\hat{\bm{\theta}}_{\mathrm{ML}})\,,$
(3)

where $p$ is the distribution of maximum-likelihood values under the true
distribution $\bm{\theta}^{*}$ given $N$ observations. In the limit that
$N\rightarrow\infty$, this distribution becomes sharply peaked around
$\bm{\theta}^{*}$. However, for finite $N$, $p$ has finite width and the
maximum-likelihood estimate $\hat{\bm{\theta}}_{\mathrm{ML}}$ is subject to
statistical fluctuations due to finite sampling, thus limiting the precision
with which the parameters can be determined. In Appendix, Section B we show
for the Bernoulli and Poisson models that the KL divergence between the true
model and the expected maximum-likelihood model plus/minus the standard error
diverges as the truth approaches the boundary, even if the standard error on
the maximum-likelihood estimate decreases to zero with increasing amounts of
data. In addition, apart from fluctuations due to finite sampling, in real-
world data the observations may additionally be corrupted, effectively
introducing an error floor in the empirical estimate.

The question then is what probabilistic model should be selected to control
the effects of fluctuations while optimising the placement of bets with finite
data?

### The playing it safe principle

We seek a model with the following properties. In the large data limit (large
$N$), the optimal model has to converge to the maximum-likelihood estimate
$\hat{\bm{\theta}}_{\mathrm{ML}}$ which itself converges to the true model
$\bm{\theta}^{*}$. In the opposite limit of small data (small $N$), we should
select a model that is less affected by statistical fluctuations. We will
parametrize this in terms of a deterministic bias
$\bm{\theta}_{\mathrm{bias}}$ and show below how to select it. For finite data
size, the optimal model should represent a compromise between the maximum-
likelihood model and the bias. Overall, we make the following Ansatz for the
model

$\displaystyle\bm{\theta}=\kappa\,\bm{\theta}_{\mathrm{bias}}+(1-\kappa)\hat{\bm{\theta}}_{\mathrm{ML}}\,,$
(4)

where the parameter $\kappa\in[0,1]$ is used to interpolate between the bias
and maximum-likelihood term. In the next section we show that this Ansatz
appears naturally in the inference of a large class of probability
distributions (exponential families with a conjugate prior). A priori, the
values of the parameters $\kappa$ and $\bm{\theta}_{\mathrm{bias}}$ can be
chosen freely, as long as $\bm{\theta}$ lies in the parameter space of the
probabilistic model. We determine the optimal choice of these parameters
through minimisation of the expected loss,
$\big{<}\mathcal{L}_{\bm{\theta}^{*}}(\bm{\theta})\big{>}_{\hat{\bm{\theta}}_{\mathrm{ML}}}$.
There are different ways to optimise with respect to these parameters. For
instance, we may keep $\kappa$ fixed and optimise for the bias, or vice versa.
These lead to different optimal models and which way we pick depends on the
setup of the inference we want to solve. We focus on the former approach to
optimisation and later discuss the second option in an example. Thus, we
optimise the expected loss with respect to the bias:

$\displaystyle\bm{\theta}^{\mathrm{opt}}_{\mathrm{bias}}=\underset{\bm{\theta}_{\mathrm{bias}}}{\arg\min}\big{<}\mathcal{L}_{\bm{\theta}^{*}}(\bm{\theta})\big{>}_{\hat{\bm{\theta}}_{\mathrm{ML}}}\,.$
(5)

To evaluate the condition, it is convenient to write it as a gradient of the
expected loss

$\displaystyle\frac{1}{\kappa}\nabla_{\bm{\theta}_{\mathrm{bias}}}\big{<}\mathcal{L}_{\bm{\theta}^{*}}(\bm{\theta})\big{>}_{\hat{\bm{\theta}}_{\mathrm{ML}}}=\nabla_{\bm{\theta}}\big{<}\mathcal{L}_{\bm{\theta}^{*}}(\bm{\theta})\big{>}_{\hat{\bm{\theta}}_{\mathrm{ML}}}=0\,.$
(6)

The gradient field is defined over parameter space and we are looking for
optimal points where the gradient vanishes. As the minima of gradient fields
are invariant under coordinate changes, the solution of the optimization
problem is invariant to changes in parameterization of the underlying
probabilistic model. However, next we also show that for a large class of
probability distributions, there exists a particular parameterization that
allows us to understand in detail how the optimum arises from an interplay of
the amount of data available about the true model and the geometry of model
space. This interplay constrains the regions of parameter space that yield
safe betting models, and we thus refer to it as the “playing it safe”
principle.

To make analytic progress, we specialise to exponential families which form a
large class of probability distributions, which include many well-known
distributions such as the Multinomial, Poisson, and Gaussian distributions
etc. [28, 29]. An exponential family is a collection of probability densities
that can be written in the canonical form

$\displaystyle
P_{\bm{\eta}}(x)=\exp\left(\bm{\eta}^{T}\cdot\mathbf{t}(x)-F(\bm{\eta})+k(x)\right)\,,$
(7)

where $\bm{\eta}$ is the _canonical parameter_ associated with the sufficient
statistic $\mathbf{t}(x)$, $F$ is the log-partition function, and $k(x)$ is
referred to as the auxiliary carrier term. Alternatively, we can use the mean
parameters $\bm{\theta}$ which are related to the canonical parameters by a
one-to-one map (see Materials and Methods). The mean parameters often also
have a simple relation to commonly used parameterizations of well-known
distributions. Examples are the event probabilities in the Multinomial
distribution, the average number of counts in the Poisson distribution, and
the variance of the Gaussian distribution with zero mean. We give examples of
these relations below and in Appendix, Sections D and E.

The KL divergence between two densities $P_{\bm{\eta}^{*}}$ and
$P_{\bm{\eta}}$ of the same exponential family, can be written as a Bregman
divergence [30, 31] between the associated mean parameters:

$\displaystyle
D_{\mathrm{KL}}(P_{\bm{\eta}^{*}(\bm{\theta}^{*})}||P_{\bm{\eta}(\bm{\theta})})=B_{F^{*}}({\bm{\theta}^{*}}:\bm{\theta})\,,$
(8)

where $F^{*}$ is the dual of the log-partition function $F$ (see Materials and
Methods). Bregman divergences form a large class of statistical distances,
including many commonly known measures, such as the squared Euclidean
distance, the Mahalanobis distance, or the Itakura–Saito distance. We note
that in our context, the particular choice of the probabilistic model
$P_{\bm{\eta}}$, implies a particular Bregman divergence through the dual
function $F^{*}$. In Appendix, Section C1 we show that for exponential
families in mean parameterisation, the optimality condition Eq. (6) takes the
simple form

$\displaystyle\big{<}\mathcal{I}(\bm{\theta})\delta\bm{\theta}\big{>}_{\hat{\bm{\theta}}_{\mathrm{ML}}}=0\,,$
(9)

which holds even for differences
$\delta\bm{\theta}\equiv\bm{\theta}^{*}-\bm{\theta}$ of finite size. In
general, the condition provides $d$ non-linear constraints on the parameters
$\kappa$ and $\bm{\theta}_{\mathrm{bias}}$ such that the expected loss is
optimised. We give an example illustrating the derivation of this condition
for the Bernoulli model in Appendix, Section C2.

We now illustrate that this condition, which holds for the mean
parameterisation, implies a bias towards models of lower complexity in the
sense described in the previous section. For this we focus on the case when
the exponential family depends on a single parameter and write out the
condition in Eq. (9) explicitly using the definitions of $\delta\theta$ and
$\theta$ in Eq. (4):

$\displaystyle\Big{<}\mathcal{I}(\theta)\left[\theta^{*}-\theta\right]\Big{>}_{\hat{{\theta}}_{\mathrm{ML}}}$
$\displaystyle=$
$\displaystyle\Big{<}\mathcal{I}(\theta)\left[\theta^{*}-\hat{{\theta}}_{\mathrm{ML}}-\kappa(\theta_{\mathrm{bias}}-\hat{{\theta}}_{\mathrm{ML}})\right]\Big{>}_{\hat{{\theta}}_{\mathrm{ML}}}$
$\displaystyle=0\,.$ (10)

Figure 2: Diverging Fisher information at the boundary leads to an asymmetric
complexity weight $\mathcal{I}(\theta)$ around the true model $\theta^{*}$. We
illustrate this for a one-dimensional parameter space with a single boundary.
The black and grey curves are the distributions of the empirically observed
maximum-likelihood value $\hat{\theta}_{\mathrm{ML}}$ for small and large $N$,
respectively. For small $N$, the distribution is wide and models to the left
of the true model $\theta^{*}$ (dashed red line), receive a higher complexity
weight when evaluating
$\langle\mathcal{I}(\theta)(\theta^{*}-\theta)\rangle_{\hat{\theta}_{\mathrm{ML}}}$
(left-hand side of Eq. (9)), than models to the right of the true model.

To illustrate the implications of this condition, we assume an infinite
parameter space with a single boundary. For example, the Poisson distribution
with parameter $\lambda\in(0,\infty)$ has such a parameter space. We consider
the scenario shown in Figure 2, where the Fisher information metric
$\mathcal{I}$ diverges to one side of the true model and falls off to small
values to the other side. We then analyse the large and small data limits as
illustrated in Figure 2, and optimise $\theta_{\mathrm{bias}}$ while $\kappa$
is fixed to a finite value:

Large $N$: The distribution of the maximum-likelihood value becomes sharply
peaked around $\theta^{*}$ and we see that in order to make the bracket inside
the expectation value vanish, we require
$\theta_{\mathrm{bias}}\rightarrow\theta^{*}$. We note that the local
complexity $\mathcal{I}$ can be arbitrarily large in this limit as it is
multiplied by the bracket with value zero. Empirical parameters
$\hat{{\theta}}_{\mathrm{ML}}$ for which
$\theta^{*}-\hat{{\theta}}_{\mathrm{ML}}$ is not zero are improbable and
multiplied by very small probability when taking the expectation value. Thus
the entire expectation value will vanish as required.

Small $N$: The distribution of the ML parameter has finite spread around the
true model $\theta^{*}$. For simplicity of argument suppose that the
distribution is symmetric around $\theta^{*}$. To understand the optimization,
we go through the possible cases for choosing the bias:

  * •

$\theta_{\mathrm{bias}}=\theta^{*}$: The expectation value becomes non-zero
(either positive or negative), because the ML values on one side of the true
model $\theta^{*}$ receive a higher complexity weight $\mathcal{I}(\theta)$
than the ML values on the opposite side.

  * •

$\theta_{\mathrm{bias}}$ on the side of $\theta^{*}$ with diverging
complexity: The expectation value assumes even larger (negative or positive)
values.

  * •

$\theta_{\mathrm{bias}}$ on the side of $\theta^{*}$ with diminishing
complexity: For the right magnitude of $\theta_{\mathrm{bias}}$ the
expectation value is adjusted to zero, because the models
$\theta(\theta_{\mathrm{bias}},\hat{\theta}_{\mathrm{ML}})$ are biased to
receive a smaller complexity weight. This demonstrates the “playing it safe”
principle for this case.

This argument in the one-dimensional case immediately generalises to the case
of higher-dimensional parameter spaces. Equation (9) provides an independent
condition for each of the $d$ dimensions of parameter space and the components
of the bias $\theta_{\mathrm{bias}}$ can thus be adjusted individually to
satisfy all conditions.

### Optimal Bayesian inference for exponential families

We now show how the “playing it safe” principle introduced in the previous
section, and in particular the Ansatz Eq. (4), arises for Bayesian inference
of exponential families with conjugate priors, for which we can provide an
analytical treatment. In short, the information contained in the prior
determines the function $\kappa$ and the bias $\bm{\theta}_{\mathrm{bias}}$ in
Eq. (4). We will show that by tuning the conjugate prior through its
hyperparameters we can control the strength of the bias relative to the
maximum-likelihood term.

We now specialise to exponential families Eq. (7) in canonical coordinates
$\bm{\eta}$. If the true model $\bm{\eta}^{*}$ is unknown, we can at best use
the $N$ independent observations $\mathbf{x}_{i}$,
$\mathbf{n}=\\{\mathbf{x}_{1},\ldots,\mathbf{x}_{N}\\}$, to determine a
Bayesian posterior distribution for ${\bm{\eta}}$,

$\displaystyle
P({\bm{\eta}}|\mathbf{n})=\frac{P(\mathbf{n}|{\bm{\eta}})P({\bm{\eta}})}{\int\mathrm{d}^{d}\bm{\psi}\,P(\mathbf{n}|\bm{\psi})P(\bm{\psi})}\,,$
(11)

where $P({\bm{\eta}})$ is a prior distribution, and
$P(\mathbf{n}|{\bm{\eta}})$ is the likelihood of seeing data ${\mathbf{n}}$
given the model parameterized by ${\bm{\eta}}$. Given the posterior, the goal
is to minimise the posterior expected loss

$\displaystyle\big{<}\mathcal{L}\big{>}_{\mathrm{post}}$
$\displaystyle=\int\mathrm{d}\tilde{\bm{\eta}}P(\tilde{\bm{\eta}}|\mathbf{n})D_{\mathrm{KL}}(Q_{\tilde{\bm{\eta}}}||Q_{\bm{\eta}})$
(12) $\displaystyle=\int
d\tilde{\bm{\theta}}\,P(\tilde{\bm{\theta}}|\mathbf{n})B_{F^{*}}(\tilde{\bm{\theta}}:\bm{\theta})\,,$
(13)

where in the second line we have used relation Eq. (8) and transformed
variables from canonical to mean parameters also in the posterior density and
integral. We want to determine the minimum value of this expression with
respect to $\bm{\theta}$

$\displaystyle\big{<}\mathcal{L}\big{>}_{\mathrm{post}}(\bm{\theta})=\underset{\bm{\theta}}{\min}\big{<}\mathcal{L}\big{>}_{\mathrm{post}}(\bm{\theta})\,.$
(14)

This optimal value is also known as the Bregman information. In [32] it was
shown that the conditional optimal parameter $\bm{\theta}$ is given by the
posterior expectation value

$\displaystyle\bm{\theta}=\int
d\tilde{\bm{\theta}}P(\tilde{\bm{\theta}}|\mathbf{n})\tilde{\bm{\theta}}\,.$
(15)

Thus, the optimal parameter depends both on the data and the choice of prior.
Our next goal is to understand how the prior choice in relation to the
available data size affects the optimal parameter.

Conjugate priors for exponential families allow us to make general analytic
statements in this respect. We denote the conjugate prior of the exponential
family by $P(\bm{\eta};\bm{\alpha})$, where $\bm{\alpha}$ is a set of
hyperparameters. In [33] it was shown for exponential families, that under the
conjugate prior the posterior mean is given by

$\displaystyle\bm{\theta}(\hat{\bm{\theta}}_{\mathrm{ML}};\bm{\alpha})=\frac{\bm{\alpha}+N\hat{\bm{\theta}}_{\mathrm{ML}}}{N_{0}+N}=\kappa\frac{\bm{\alpha}}{N_{0}}+(1-\kappa)\hat{\bm{\theta}}_{\mathrm{ML}}\,,$
(16)

where $N_{0}$ is a function of the hyperparameters $\bm{\alpha}$, and
$\kappa\equiv N_{0}/(N_{0}+N)$. From this expression we can understand the
role played by the hyperparameters. $N_{0}(\bm{\alpha})$ is an threshold
number of samples determined by the hyperparameters that sets the scale for
whether the actual sample size $N$ is large enough for the evidence to
overcome the prior. Likewise, $\bm{\alpha}/N_{0}$ functions as effective
additional parameters that appear on equal footing as the maximum-likelihood
parameter. In fact, the quotient $\bm{\alpha}/N_{0}$ is the prior expectation
for the parameter $\bm{\theta}$, and we denote it by
$\bm{\theta}_{\mathrm{prior}}\equiv\mathbb{E}(\bm{\theta}|\bm{\alpha})=\bm{\alpha}/N_{0}$.
Overall, this leaves us with the following form of the conditional optimal
parameter:

$\displaystyle\bm{\theta}=\kappa\bm{\theta}_{\mathrm{prior}}+(1-\kappa)\hat{\bm{\theta}}_{\mathrm{ML}}\,,$
(17)

with

$\displaystyle\kappa=\frac{N_{0}/N}{N_{0}/N+1}=\frac{q}{q+1}\,,$ (18)

where $q\equiv N_{0}/N$. We find $\kappa\rightarrow 0$ for $N_{0}/N\rightarrow
0$, and $\kappa\rightarrow 1$ for $N_{0}/N\rightarrow\infty$. Finally,
comparison with Eq. (4) shows that $\bm{\theta}_{\mathrm{prior}}$ takes the
role of the bias $\bm{\theta}_{\mathrm{bias}}$, such that overall this
confirms our initial Ansatz.

In the above derivation, the conjugate prior $P(\bm{\eta};\bm{\alpha})$ was
defined on the canonical parameters and lead to the linear convex form of Eq.
(17). However, in the next section and in (Appendix, Sections D and E), we
show for the multinomial, Poisson, and Gaussian model, that defining the
conjugate prior on the event probabilities, the mean number of counts, and the
variance, respectively, equally results in the optimal parameter having convex
linear form. We emphasise that for the development of our theory it is not
important on which set of parameters the conjugate prior is defined. What
counts for our argument, is that the final expression for the conditional
optimal parameter has the linear convex form.

The conditional optimal parameter in convex linear form depends on the data
and crucially also on the shape of the prior, set by its hyperparameters.
While for a given prior, the loss is minimised by the conditional optimal
model, we may further optimise the expected loss with respect to the prior
choice which we do by optimising the hyperparameters. Given the convex linear
form in Eq. (17), we see that there are two natural ways to optimise the
hyperparameters. The first, is to make a choice for
$\bm{\theta}_{\mathrm{prior}}$, and then optimise the value of $\kappa$ (or,
equivalently $q$). The second way is to set the value of $\kappa$ and then
optimise $\bm{\theta}_{\mathrm{prior}}$. We now illustrate both paths of
optimisation for the multinomial distribution.

### Optimization for the multinomial model

Figure 3: Optimization of the expected loss with respect to the prior shows a
bias to low complexity models. The true model lies at $\theta^{*}=0.2$. (A)
The bias $\theta_{\mathrm{bias}}$ is optimised while $\kappa=0.2$ ($q=0.25$)
is kept fixed. The plot shows the landscape of the expected loss as a function
of the bias and the number of observations. We are using a non-linear colormap
based on the cumulative distribution of the loss. (B) We show the ratio of the
local complexity as measured by Jeffreys prior of the expected optimal to that
of the true model for three different true models $\theta^{*}$. (C) The
parameter $\kappa$ is optimised while the bias $\theta_{\mathrm{bias}}=1/2$
(imposing isotropy on hyperparameter space). The plot shows the landscape of
the expected loss as a function of the parameter $\kappa$ and the number of
observations. (D) We show the ratio of the local complexity as measured by
Jeffreys prior of the expected optimal to that of the true model for three
different true models $\theta^{*}$.

The multinomial model is the probability distribution of $N$ events, each
assuming one of $d$ possible outcomes (Materials and Methods). The true model
is parameterised by the $d$-dimensional set of probabilities
${\bm{\theta}}^{*}$. We note that these probabilities are related to the mean
parameters of the multinomial distribution by rescaling the mean parameters
with a factor of $1/N$. For simplicity of notation and to emphasize their
close relation, we use the same symbol $\bm{\theta}$ to denote the
probabilities that we used above for the mean parameters. It is natural to
collect information on the model by counting event outcomes and computing
empirical probabilities $\hat{\bm{\theta}}$. The data $\mathbf{n}$ consists of
the counts of occurrence $x_{i}$ for each of the $d$ possible outcomes over
$N$ observations ($\sum_{i=1}^{d}x_{i}=N$). Following our earlier discussion,
we use Bayesian inference with a conjugate prior and optimise the conditional
model to minimise the expected loss. The likelihood function is the
multinomial distribution and its conjugate prior is a Dirichlet distribution
given by

$\displaystyle
P(\mathbf{n}|{\bm{\theta}})=\frac{1}{B(\mathbf{n})}\prod_{i=1}^{d}\theta_{i}^{x_{i}}\,,\quad
P({\bm{\theta}};\bm{\alpha})=\frac{1}{B(\bm{\alpha})}\prod_{i=1}^{d}\theta_{i}^{\alpha_{i}-1}\,.$
(19)

The set of $d$ hyperparameters $\alpha_{i}>0$ determines the prior shape.
Direct calculation yields, as expected, the linear convex form Eq. (17) of the
conditional optimal model

$\displaystyle\theta_{i}(\hat{\bm{\theta}};\bm{\alpha})=\frac{\frac{\alpha_{i}}{N}+\hat{\theta}_{i}}{\frac{1}{N}\sum_{i}\alpha_{i}+1}\,,$
(20)

where $\hat{\theta}_{i}=x_{i}/N$, are the empirical probabilities. The
expression in Eq. (20) is also known as a pseudo-Laplacian law [15, 34, 35]
with the $\alpha_{i}$ being “pseudocounts”. We now specialise to $d=2$, i.e.
two possible event outcomes. The conditional optimal model depends on two
parameters, $\alpha_{1}$ and $\alpha_{2}$ from which we define

$\displaystyle
N_{0}\equiv\alpha_{1}+\alpha_{2}\,,\;\bm{\theta}_{\mathrm{bias}}\equiv\frac{\bm{\alpha}}{N_{0}}\,.$
(21)

The next step is to optimize these hyperparameters in the two possible ways
outlined in the previous section.

First, we consider optimization of the bias $\bm{\theta}_{\mathrm{bias}}$,
while $\kappa$ is fixed. Since $\theta_{2}=1-\theta_{1}$, we are dealing with
optimization in one dimension. For simplicity of notation, we suppress
component indices in the following. The landscape of the KL loss is shown in
Figure 3A for $\theta^{*}=0.2$ and $\kappa=0.2$. The curve of the optimal bias
choice and the optimized parameter $\langle\theta_{\mathrm{opt}}\rangle$
(averaged over empirical observations), as a function of $N$ is shown in
white. Asymptotically, the curve approaches $\theta^{*}$ for both parameters.
In Figure 3B, we show how the local complexity of the model with optimal
parameter compares to the complexity of the true model as a function of $N$
and for different values of the true model. Asymptotically, the complexity
ratio approaches unity, but the closer the true model lies to the boundary,
the more data is required to reach the full complexity of the true model. We
also observe in Figs. 3A and B that for small $N$ there is an initial dip in
the curves. At first sight this may seem to be in tension with the “playing it
safe” principle. However, this is explained and resolved by the existence of
the second boundary of the parameter space of the $d=2$ multinomial model at
$\theta=1$. For small $N$, the distribution of maximum-likelihood values is so
broad that diverging complexity at the second boundary influences the
optimization, pushing the parameter away from the boundary.

Next, we consider optimization of the parameter $\kappa$, while
$\theta_{\mathrm{bias}}$ is held fixed. We make the choice
$\alpha_{1}=\alpha_{2}$ which, by Eq. (21) implies
$\theta_{\mathrm{bias}}=1/2$. This is an appealing choice as it implies equal
prior knowledge in all directions of parameter space. The landscape of the KL
loss is shown in Figure 3C with $\theta^{*}=0.2$. The curve of the optimal
choice for $\kappa$ and the optimized parameter
$\langle\theta_{\mathrm{opt}}\rangle$ as a function of $N$ is shown in white.
Asymptotically, $\kappa$ approaches zero and the optimal model approaches
$\theta^{*}$. In Figure 3D, we show for different true models, that the
complexity ratio approaches unity in the asymptotic limit as expected, but the
closer the true model lies to the boundary, the more data is required.

We thus see that either way of optimising implies “playing it safe” as
strategy. In (Appendix, Sections D and E), we demonstrate the optimization of
the bias for the Poisson model with unknown average number of counts, and the
Gaussian model with unknown variance parameter.

### Playing it safe in bacterial phenotypic switching

Figure 4: Bacteria that place bets according to the “playing it safe”
principle grow more safely. (A) Long-term growth rate of the bettor’s capital
for Kelly betting on independent Bernoulli events as a function of the
bettor’s assumed model $\theta$ for three values of the true model
$\theta^{*}$ (orange, red, and purple). The maximum is reached when
$\theta=\theta^{*}$. (B) Three bacterial populations, each of initial size
$10$, are exposed to a sequence of binary environmental states (growth or
antibiotic) and the development of population sizes is followed. The
individuals of each population stochastically switch between a growth and
resistant phenotype independently from one and another. All populations keep a
memory of the last $N=5$ past environmental states, but each population uses a
different learning strategy, set by the hyperparameters $N_{0}=2\alpha$ and
$\theta_{\mathrm{bias}}=1/2$. The environmental state sequence contains a
single antibiotic state (skull), while all other states are growth states
(apple). For each population we show the mean population size and the standard
deviation. The population with the least safe strategy (dark blue) grows the
fastest, but is extinguished by the antibiotic environment. The population
with the safest strategy (yellow) grows more slowly, but shows the smallest
loss due to antibiotic treatment. The population with intermediate safety
(teal), experiences larger loss after antibiotic treatment, but outgrows the
population with highest safety. The learning strategy of intermediate safety
represents a compromise between too much risk and too much safety. Population
curves are computed from fifty independent trials and data points have a small
horizontal offset to increase legibility. (C) The magnitude of fluctuations of
the Bayesian conditional optimal model
$\theta(\hat{\theta};N_{0},\theta_{\mathrm{bias}}=1/2)$ is shown as the plus
and minus standard error of the conditional optimal model,
$\theta_{\pm}=\theta(\theta^{*}\pm\sigma_{N};N_{0},\theta_{\mathrm{bias}}=1/2)$,
as a function of the number of observations $N$ and the prior choice. Three
prior choices from the family of the conjugate prior, set by different values
of the hyperparameter $N_{0}$ (dark blue, teal and yellow) and the fixed bias
$\theta_{\mathrm{bias}}=1/2$. The choice $N_{0}=0$ corresponds to maximum-
likelihood estimation. For large number $N$ of observations, the fluctuations
diminish ($q\rightarrow 0$ and $\kappa\rightarrow 0$) and the conditional
optimal model $\theta$ converges to the true model $\theta^{*}=0.2$ (dotted
line). The safe region (gray) constitutes models $\theta$ for which the growth
rate function $G_{\theta^{*}}(\theta)$ is non-negative. For larger values of
$N_{0}$, the standard error of the conditional optimal model $\theta$ is
pulled into the safe region, such that a smaller loss in the growth rate is
expected.

Betting strategies in biology are often found in collective systems where bets
can be distributed (or hedged), across the possible outcomes of random events.
There are many cases of bet hedging systems in biology [36, 10, 37]; a
particularly well-known example is stochastic phenotypic switching in
bacterial populations [38, 2]. In the simplest setup, the population lives in
an environment which can assume one of two different states: a nutrient-rich
state, allowing the bacterial population to grow, or an antibiotic, hazardous
state that kills off bacteria. Each bacterium switches stochastically and
independently from the others between a growth and a resistant phenotype. The
growth phenotype prospers in the growth state of the environment, but is
instantly killed if the environment is in the antibiotic state. Conversely,
the resistant phenotype does not grow in the growth environment, but survives
antibiotic treatment. It has been shown that in order to maximise population
growth, the bacteria should perform Kelly betting (also known as proportional
betting) [24, 2]. According to Kelly’s theory, the bacteria should
stochastically select each phenotype with a probability that is proportional
to the probability of occurrence of the corresponding environmental state. We
consider the case where the bacteria only have imprecise information about
environmental probabilities, show how the “playing it safe” principle implies
safer betting strategies in this case, and how these choices translate into
optimisation of growth of the bacterial population.

In our toy model of stochastic phenotypic switching, successive discrete
environmental states are chosen independently according to a Bernoulli model
with a probability $\theta^{*}$ for the nutrient state, and probability
$1-\theta^{*}$ for the antibiotic states. The bacteria can maximise population
growth by betting according to Kelly’s theory, namely by matching their
probabilistic model to the environment. In Kelly betting any mismatch is
penalised and the loss is captured by the KL divergence between the true
probabilistic model $P_{\theta^{*}}$ and the model $P_{\theta}$ assumed by the
bacterial population (the bettor). Instead of the loss, it is also instructive
to consider the closely related growth rate function derived by Kelly in [24,
25]:

$\displaystyle
G_{\theta^{*}}(\theta)=[H(P_{1/2})-H(P_{\theta^{*}})]-D_{\mathrm{KL}}(P_{\theta^{*}}||P_{\theta})\,,$
(22)

where $H(\cdot)$ is the entropy. This function captures the long-term
population growth rate, i.e., the growth rate after a long succession of
environmental states. The term in brackets is the maximal achievable growth
rate from which the KL loss is subtracted. In Figure 4A we show the long-term
growth rate $G_{\theta^{*}}$ of Kelly betting. The growth rate has a single
maximum when the assumed model is equal to the truth, $\theta=\theta^{*}$ for
which the loss is zero. The further away the adopted model lies from the
truth, the smaller the growth rate; it can even become negative if the assumed
model lies too far away from the truth. We will see (Figures 4A,B,C and S7)
that there is a region of safe model choices. This region is given by models
$\theta$, for which the long-term population growth rate is non-negative.
Models in this region are safe, because while they might not be optimal, they
guarantee that there are no long-term losses.

We now assess, for Kelly betting, how the amount of information available
about the true model and the prior choice affect the magnitude of fluctuations
in the assumed model. The empirical observations are subject to statistical
fluctuations and as discussed above translate to the conditional optimal
model. To quantify the magnitude of fluctuations, we consider the standard
error of the empirical estimate $\hat{\theta}$. For $N$ observations the
standard error is given by $\sigma_{N}=\sqrt{\theta^{*}(1-\theta^{*})/N}$. The
fluctuations of the empirical estimate, determine the fluctuations of the
conditional optimal model. The magnitude of fluctuations of the conditional
optimal strategy Eq. (17) within one standard deviation, is given by the
bounds
$\theta_{\pm}\equiv\theta(\hat{\theta}^{\pm}_{\mathrm{ML}};N_{0},\theta_{\mathrm{bias}})$,
where $\hat{\theta}^{\pm}_{\mathrm{ML}}=\theta^{*}\pm\sigma_{N}$.

We will examine a toy model of bacterial phenotypic switching, in which
environmental states are chosen and phenotypes expressed simultaneously at
discrete time steps. The bacterial population then performs collective
betting. We know that with perfect information about the environmental
probabilities the optimal strategy would be to perform proportional betting.
This means matching the probability of expressing the growth phenotype with
that of the occurrence of the growth state of the environment. However, the
bacteria do not know the probability of environmental states. We will suppose
that instead they keep an implicit memory of the previous $N$ environmental
states, which they use to compute empirical probabilities.

To illustrate the effect of safe learning, we compare the growth of three
types of bacterial colonies distinguished only by their learning strategy, as
they respond to a sequence of environmental states which include a single
antibiotic episode. In Figure 4B we show how the size of each population
develops over time. All bacteria have a memory over the past $N=5$
environmental states. In the extreme case of maximum-likelihood estimation
($N_{0}=0$, dark blue curve), the population is killed completely by the
antibiotic state. In contrast, the safer the learning strategy (larger
$N_{0}$), the smaller the loss a bacterial population suffers due to
antibiotic treatment. After the antibiotic treatment, bacteria with the safest
strategy recover their numbers most quickly (yellow curve), but they are
eventually overtaken by bacteria with less safe learning strategies (teal
curve). The observed trade-off between the ability to survive antibiotic
treatment and achieve high growth rates is the reason for the existence of an
optimally safe learning strategy. The precise choice of optimal strategy
depends on the statistics of the environment, as well as potentially other
factors such as a future horizon over which safe growth needs to be
guaranteed.

In environments which fluctuate more the optimal learning strategy will need
to be safer, because the probabilities inferred from a finite memory of past
events will have a larger variance. However, with a maximally safe learning
strategy ($\theta=0.5$), the bacterial population will not be able to achieve
any long-term growth, but only keep its population size constant as shown in
Figure 4A. In turn, this means that long-term physical growth is only possible
in environments which are not too random. In fact, the more environments are
suitable for life to prosper, the closer they are to being deterministic.
Near-deterministic environments are thus livable, but as we have shown, rare,
yet potentially catastrophic events require organisms to “play it safe” to
survive in even in such worlds.

Finally, in Figure 4C we contrast different prior choices to show how these
affect model safety and plot the standard error bounds of the conditional
optimal model $\theta_{\pm}$ as function of the observations $N$. In this
example we fix $\theta_{\mathrm{bias}}=1/2$ which, by Eq. (21), implies
$N_{0}=2\alpha$. We observe that a large value of $N_{0}$ reduces the
magnitude of statistical fluctuations and biases the conditional optimal model
$\theta_{\mathrm{opt}}$ towards a probability of $1/2$, i.e. minimal
complexity. Overall, the inferred model is pulled into the safe growth region
and larger values of $N_{0}$ thus result in a safer learning strategy. From
(Appendix, Section F and Figure S6) we see that this effect is the more
pronounced the closer the true model lies to the boundary of parameter space.
Notably, maximum-likelihood estimation ($N_{0}=2\alpha=0$) is the least safe
strategy. For a large number of observations $N$, the fluctuations diminish
and the inferred model converges to the true model. However, we observe that
safer learning strategies (larger $N_{0}$) converge to the true model more
slowly than less safe strategies. Furthermore, in the limit of maximum safety
$q=N_{0}/N\rightarrow\infty$, the optimal model
$\theta_{\mathrm{opt}}\rightarrow 1/2$ for which the growth rate vanishes.
This illustrates the trade-off between safe learning and growth maximisation.

## Discussion

We have proposed that biological systems acting with bounded information
should “play it safe” by biasing themselves towards models of the environment
that are less complex, in the sense of being easier to infer despite
statistical fluctuations in limited observed data. We derived this result
analytically for exponential families, a very general class of probability
distributions, including, e.g., the widely employed multinomial, Poisson, and
Gaussian distributions. In the context of collective betting, our results
imply that optimal adaptation involves a balance between risk created by
uncertainty in inference of the latent distribution, and making full use of
the potential for growth or reward.

We illustrated these ideas in the example of bacteria in an environment which
mostly supports growth, but occasionally manifests an antibiotic. In this
situation, it pays to maintain diversity between high-growth and low-growth
but resistant phenotypes, as a bet-hedging mechanism [38, 2]. Indeed, wild
type Escherichia coli show a resistant fraction of $\sim 10^{-6}-10^{-5}$
[39]. We showed how this sort of strategy arises in terms of “playing it safe”
in an uncertain world. It will be interesting in the future to apply our
analysis to other situations where populations of organisms must make
decisions with limited information, for example in bacterial chemotaxis [40]

We also cast our analysis in the language of Bayesian inference, and noted
that the bias towards simpler models can be implemented by an appropriate
choice of prior. The optimal prior depends on knowledge of the true model, or
rather its distribution in cases where the generative process in the world
itself changes stochastically. It would be interesting to determine the sorts
of evolutionary dynamics that would allow systems to adjust their priors to
“play it safe” across a given distribution of generative processes.

Another interesting question is to understand how our considerations operate
if we use the information maximizing discrete priors discussed in [41, 42]. We
instead used the continuous conjugate priors because it seemed natural that a
biological agent might have ways of tuning over a continous family of choices,
e.g., by adjusting concentration levels of a protein, or activity levels of a
circuit. By contrast, it is harder to imagine biological methods for selecting
over a discrete set of values. Neverthless, it is important and interesting to
consider how our “playing it safe” principle might interact with the discrete
priors of [41, 42] because the authors of these articles are precisely
concerned with the limitations imposed on inference by having a finite amount
of data. Indeed, their discrete priors depend on the number of data points
$N$, and are designed to maximize the amount of information gained about the
parameters from each observation. Perhaps there is some way for an organism to
incorporate an understanding of the amount of data it is able to accumulate
from the past, to effectively discretize its prior on the parameters of the
world. It would be useful to understand whether and how some variant of our
“playing it safe” principle could be realized in that context. Models of
biological systems often depend on a very large number of parameters which
have been shown to be close to impossible to determine with high precision
[43]. This important observation may suggest that we should define a “playing
it safe” principle on effective parameters.

Finally, in a time-varying world, inferences are only useful up to some future
time horizon. This restricted future utility should limit the value and
precision required of the optimal model, while also increasing the imperative
for “safety” in light of possible future changes in the world. We can expect
trade-offs associated with keeping past memory to play a role in such
situations [44, 45]. Perhaps this observation also has a bearing on the
widespread observation that Occam’s Razor, understood as a preference for
simplicity or parsimony, seems to be an organizing principle for mental
functions [46, 47, 48, 49, 50, 51], and that biases towards simple models are
helpful in many statistical settings that require inference and prediction
[52, 22, 53, 54, 41, 55].

## Materials and Methods

### Information geometry basics

We define the FIM by Taylor expansion of the KL divergence between two
parametric distributions $P_{\bm{\theta}}$ and $P_{\bm{\phi}}$. The Taylor
expansion between two models which are close in parameter space
${\bm{\phi}}\equiv{\bm{\theta}}^{*}-\delta{\bm{\theta}}$ is given by

$\displaystyle
D_{\text{KL}}(P_{\bm{\theta}^{*}}||P_{{\bm{\theta}^{*}}-\delta{\bm{\theta}}})=\frac{1}{2}\mathcal{I}_{ij}\delta\theta^{i}\delta\theta^{j}+\ldots\,,$
(23)

where summation over repeated indices is implied and the FIM $\mathcal{I}$ is
given by the coefficient of the quadratic term

$\displaystyle\mathcal{I}_{ij}({\bm{\theta}})=\partial_{\delta\theta^{i}}\partial_{\delta\theta^{j}}D_{\text{KL}}(P_{\bm{\theta}}||P_{{\bm{\theta}}-\delta{\bm{\theta}}})\Big{|}_{\delta{\bm{\theta}}=\mathbf{0}}\,.$
(24)

As an important example we consider the multinomial distribution

$\displaystyle
P(\mathbf{x};{\bm{\theta}})=\frac{1}{B(\mathbf{x})}\prod_{i=1}^{d}\theta_{i}^{x_{i}}\,,\text{
with }\;\sum_{i=1}^{d}\theta_{i}=1\,,\;x_{i}\in\\{0,...,N\\}\,,$ (25)
$\displaystyle\sum_{i=1}^{d}x_{i}=N\,,\text{and multivariate Beta
}B(\mathbf{x})\equiv\frac{\prod_{i=1}^{d}\Gamma(x_{i})}{\sum_{i=1}^{d}\Gamma(x_{i})}$

where the parameters $\theta_{i}$ are probabilities and the mean is
$N{\bm{\theta}}$. For the multinomial distribution Jeffreys prior takes the
form of a Dirichlet distribution

$\displaystyle w({\bm{\theta}})$
$\displaystyle=\frac{\Gamma(d/2)}{\sqrt{\pi}^{d}}\prod_{i=1}^{d}\frac{1}{\sqrt{\theta_{i}}}$
(26)

with $\theta_{d}=1-\sum_{i=1}^{d-1}\theta_{i}$. From this form we see that
Jeffreys prior has a minimum at $\theta_{i}=1/d$, but diverges when any of the
parameters tend to zero, $\theta_{i}\rightarrow 0$. Analogously, in (Appendix,
Sections A) we show for the Poisson and Gaussian distribution, respectively,
that Jeffreys prior diverges on the boundaries of the standard
parameterization of these distributions.

### Exponential families

An exponential family consists of a collection of densities

$\displaystyle
P_{\bm{\eta}}(x)=\exp\left(\bm{\eta}^{T}\cdot\mathbf{t}(x)-F(\bm{\eta})+k(x)\right)\,,$
(27)

where $\bm{\eta}=[\eta_{1}\;\eta_{2}\ldots\;\eta_{d}]^{T}$ is the _canonical
parameter_ associated with the sufficient statistic $\mathbf{t}(x)$, $k(x)$ is
the auxiliary carrier term, and $F$ is the log-partition function:

$\displaystyle
F(\bm{\eta})=\log\left(\int\exp(\bm{\eta}^{T}\cdot\mathbf{t}(x)+k(x))\mathrm{d}\nu(x)\right)\,,$
(28)

where $\nu(x)$ the Lebesgue or counting measure. The log-partition function is
convex. If the sufficient statistic is a collection of $d$ functions, the
canonical parameter $\bm{\eta}$ takes values in the convex set

$\displaystyle\Omega=\\{\bm{\eta}\in\mathbb{R}^{d}\,|\,F(\bm{\eta})<+\infty\\}\,.$
(29)

We consider exponential families where $\Omega$ is an open set, in which case
the exponential family is called ${\it regular}$. Since the log-partition
function is a convex function, the set $\Omega$ is also convex. We restrict
ourselves to minimal exponential families, which are defined by not having
linear constraints amongst the parameters and also not amongst the components
of the sufficient statistic. Three examples appear in the table below.

Distribution | $\mathbf{t}(x)$ | $F(\bm{\eta})$ | $k(x)$
---|---|---|---
Multinomial | $[x_{1}\;x_{2}\;\ldots\;x_{d-1}]^{T}$ | $N\log\left(1+\sum_{i=1}^{d-1}e^{\eta_{i}}\right)$ | $\log\frac{N!}{\prod_{i=1}^{d}x_{i}!}$
Poisson | $x$ | $e^{\eta}$ | $-\log x!$
Gaussian | $[x\;x^{2}]^{T}$ | $-\frac{\eta_{1}^{2}}{4\eta_{2}}+\frac{1}{2}\log(-\pi/\eta_{2})$ | $0$

An alternative set of parameters is given by the mean parameters $\bm{\theta}$
defined as the expectation of the sufficient statistic

$\displaystyle\bm{\theta}=\mathbb{E}[\mathbf{t}(x)]=\nabla_{\bm{\eta}}F(\bm{\eta})\,,$
(30)

where the last relation can be verified by direct calculation. The
distribution in mean parameterisation is obtained by inverting this defining
equation of the mean parameter and substituting for $\bm{\eta}$ in the
canonical form of the exponential family. The analogue of the log-partition
function $F$ for the mean parameters is the dual function of $F$ defined by
the Legendre-Fenchel transform

$\displaystyle
F^{*}(\bm{\theta})=\sup_{\bm{\eta}\in\Omega}\left\\{\bm{\eta}^{T}\cdot\bm{\theta}-F(\bm{\eta})\right\\}\,.$
(31)

The dual function $F^{*}$ can thus be obtained by explicit substitution for
$\bm{\eta}$ in terms of $\bm{\theta}$, after inverting the defining equation
of the mean parameter. Finally, we note a relation between $F^{*}$ and
negative entropy $-H(P_{\bm{\eta}(\bm{\theta})})$ of the distribution

$\displaystyle-H(P_{\bm{\eta}(\bm{\theta})})$ $\displaystyle=\int
P_{\bm{\eta}(\bm{\theta})}\log P_{\bm{\eta}(\bm{\theta})}\mathrm{d}\nu(x)$
(32)
$\displaystyle=\mathbb{E}[\bm{\eta}^{T}\cdot\mathbf{t}(x)-F(\bm{\eta})]+\mathbb{E}[k(x)]$
(33) $\displaystyle=F^{*}(\bm{\theta})+\mathbb{E}[k(x)]\,,$ (34)

where we have used the definition of the mean parameters and the defining
relation of the dual function. For distributions which afford a vanishing
auxiliary carrier term $k(x)$, the dual function is thus given by the negative
entropy. This is for example the case for the Gaussian distribution.

###### Acknowledgements.

This work was supported in part by the Simons Foundation MMLS Grant 400425 and
by the NIH grant R01EB026945. VB thanks the Galileo Galileo Institute where he
was a Simons Visiting Scientist, and the Aspen Center for Physics (supported
by NSF grant PHY-1607611), for hospitality as this work was completed.

##  Appendix A Jeffreys prior for discrete and continuous distributions

We compute Jeffreys prior for common examples of discrete and continuous
distributions. To this end, we define the Fisher information metric (FIM) as
the coefficient of the second order term in the Taylor expansion of the
Kullback-Leibler (KL) divergence Eq. (1). If two parameterised distributions,
are defined by discrete sets of parameters
${\bm{\theta}}=\\{\theta_{1}\,\ldots,\theta_{d}\\}$ and
${\bm{\phi}}=\\{\phi_{1}\,\ldots,\phi_{d}\\}$, respectively, we consider the
Taylor expansion of the KL divergence between two models which are close in
parameter space ${\bm{\phi}}\equiv{\bm{\theta}}-\delta{\bm{\theta}}$:

$\displaystyle
D_{\text{KL}}(P_{\bm{\theta}}||P_{\bm{\theta}-\delta\bm{\theta}})=\frac{1}{2}{\mathcal{I}}_{ij}\delta\theta^{i}\delta\theta^{j}+\ldots\,,$
(35)

where summation over repeated indices is implied and $\mathcal{I}$ is the FIM,
given by

$\displaystyle{\mathcal{I}}_{ij}({\bm{\theta}})=\partial_{\delta\theta^{i}}\partial_{\delta\theta^{j}}D_{\text{KL}}(P_{\bm{\theta}}||P_{\bm{\theta}-\delta\bm{\theta}})\Big{|}_{\delta{\bm{\theta}}=\mathbf{0}}\,.$
(36)

The normalised Jeffreys prior is computed from the FIM through the relation

$\displaystyle
w({\bm{\theta}})=\frac{\sqrt{\det\mathcal{I}({\bm{\theta}})}}{\int\mathrm{d}^{d}{\bm{\psi}}\sqrt{\det\mathcal{I}({\bm{\psi}})}}\,.$
(37)

###  A.1 Multinomial distribution

We define the multinomial distributions with parameters ${\bm{\theta}}$ and
${\bm{\theta}}^{\prime}$ as:

$\displaystyle
P_{\bm{\theta}}(\mathbf{x})=\frac{1}{B(\mathbf{x})}\prod_{i=1}^{d}\theta_{i}^{x_{i}}\,,\;Q_{\bm{\theta}^{\prime}}(\mathbf{x})=\frac{1}{B(\mathbf{x})}\prod_{i=1}^{d}{\theta^{\prime}}_{i}^{x_{i}}$
(38)

with

$\displaystyle\sum_{i=1}^{d}\theta_{i}=\sum_{i=1}^{d}{\theta^{\prime}}_{i}=1\,,\;x_{i}\in\\{0,...,N\\}\text{
and }\sum_{i=1}^{d}x_{i}=N\,.$ (39)

The KL divergence between the two distributions is given by

$\displaystyle D_{\mathrm{KL}}(P_{\bm{\theta}}||Q_{\bm{\theta}^{\prime}})$
$\displaystyle=\sum_{\mathbf{x}}\frac{1}{B(\mathbf{x})}\prod^{d}_{i=1}\theta_{i}^{x_{i}}\log\left(\frac{\prod_{j}\theta_{j}^{x_{j}}}{\prod_{j}{\theta^{\prime}}_{j}^{x_{j}}}\right)$
(40)
$\displaystyle=\sum_{j}(\log\theta_{j}-\log\theta_{j}^{\prime})\left(\sum_{\mathbf{x}}\frac{x_{j}}{B(\mathbf{x})}\prod_{i=1}^{d}\theta_{i}^{x_{i}}\right)$
(41)
$\displaystyle=\sum_{j}(\log\theta_{j}-\log\theta_{j}^{\prime})N\theta_{j}\,,$
(42)

where in the last step we have evaluated the mean of the multinomial
distribution. The next step is to Taylor expand the KL divergence around
${\bm{\theta}}$, by setting
${\bm{\theta}}^{\prime}={\bm{\theta}}-\delta{\bm{\theta}}$. For the Taylor
expansion we find

$\displaystyle D_{\mathrm{KL}}(P_{\bm{\theta}}||Q_{\bm{\theta}^{\prime}})$
$\displaystyle=\sum_{j}\left(\frac{\delta\theta_{j}}{\theta_{j}}+\frac{\delta\theta_{j}^{2}}{2\theta_{j}^{2}}+O(\delta\theta_{j}^{3})\right)N\theta_{j}=\sum_{j}\left(\frac{\delta\theta_{j}^{2}}{2\theta_{j}^{2}}+O(\delta\theta_{j}^{3})\right)N\theta_{j}\,,$
(43)

where in the last step the linear term vanishes since
$\sum_{i}\theta_{i}=\sum_{i}\theta_{i}^{\prime}=1$ and
$\theta_{i}^{\prime}=\theta_{i}-\delta\theta_{i}$. From the quadratic term we
read off the form of the FIM

$\displaystyle{\mathcal{I}}_{ij}({\bm{\theta}})=N\frac{\delta_{ij}}{\theta_{i}}\,,$
(44)

where $\delta_{ij}$ is the Kronecker delta. Taking the square root of the
determinant and computing the normalisation factor finally gives us Jeffreys
prior:

$\displaystyle
w({\bm{\theta}})=Dir({\bm{\theta}};1/2,\ldots,1/2)=\frac{\Gamma(d/2)}{\sqrt{\pi}^{d}}\frac{1}{\sqrt{1-\sum_{i=1}^{d-1}\theta_{i}}}\prod_{i=1}^{d-1}\frac{1}{\sqrt{\theta_{i}}}\,,$
(45)

where $\mathrm{Dir}$ is the Dirichlet distribution.

###  A.2 Poisson distribution

We consider two Poisson distributions with parameters ${\lambda^{\prime}}$ and
$\lambda$:

$\displaystyle
P_{\lambda^{\prime}}(n)=\frac{{\lambda^{\prime}}^{n}e^{-{\lambda^{\prime}}}}{n!}\,,\;Q_{\lambda}(n)=\frac{\lambda^{n}e^{-\lambda}}{n!}\,.$
(46)

The KL divergence between the two distributions is given by

$\displaystyle D_{\mathrm{KL}}(P_{\lambda^{\prime}}||Q_{\lambda})$
$\displaystyle=\sum_{n}\frac{{\lambda^{\prime}}^{n}e^{-{\lambda^{\prime}}}}{n!}\log\left(\frac{{\lambda^{\prime}}^{n}e^{-{\lambda^{\prime}}}}{n!}\frac{n!}{\lambda^{n}e^{-\lambda}}\right)={\lambda^{\prime}}\log\frac{{\lambda^{\prime}}}{\lambda}-{\lambda^{\prime}}+\lambda\,.$
(47)

The next step is to Taylor expand the KL divergence around
${\lambda^{\prime}}$, by setting $\lambda={\lambda^{\prime}}+\delta\lambda$.
For the Taylor expansion we find

$\displaystyle D_{\mathrm{KL}}(P_{\lambda^{\prime}}||Q_{\lambda})$
$\displaystyle=\frac{\delta\lambda^{2}}{2{\lambda^{\prime}}}+O(\delta\lambda^{3})\,.$
(48)

From the quadratic term we read off the form of the FIM

$\displaystyle{\mathcal{I}}({\lambda^{\prime}})=\frac{1}{{\lambda^{\prime}}}\,.$
(49)

Taking the square root gives us Jeffreys prior:

$\displaystyle w({\lambda^{\prime}})=\frac{1}{\sqrt{{\lambda^{\prime}}}}\,,$
(50)

Jeffreys prior is shown in Figure S1. It diverges on the boundary of parameter
space for ${\lambda^{\prime}}\rightarrow 0$ and approach zero as for
${\lambda^{\prime}}\rightarrow\infty$.

Figure S1: Jeffreys prior for the Poisson distribution with rate $\lambda$,
diverges at the boundary of parameter space for $\lambda\rightarrow 0$.

###  A.3 Gauss distribution

The KL divergence for continuous distributions is defined as

$\displaystyle
D_{\mathrm{KL}}(P||Q)=\int\mathrm{d}x\,P(x)\log\frac{P(x)}{Q(x)}\,.$ (51)

The KL divergence between two Gaussian distributions with mean $\mu$ and
variances $\sigma^{2}$ and $\sigma^{\prime 2}$

$\displaystyle
P_{\mu,\,\sigma^{2}}(x)=\frac{1}{\sqrt{2\pi\sigma^{2}}}e^{-\frac{(x-\mu)^{2}}{2\sigma^{2}}}\,,\;P_{\mu,\,\sigma^{\prime
2}}(x)=\frac{1}{\sqrt{2\pi\sigma^{\prime
2}}}e^{-\frac{(x-\mu)^{2}}{2\sigma^{\prime 2}}}\,,$ (52)

is given by

$\displaystyle D_{\mathrm{KL}}(P_{\mu,\,\sigma^{2}}||P_{\mu,\,\sigma^{\prime
2}})=-\frac{1}{2}\left[\log\sigma^{2}-\log\sigma^{\prime
2}+1-\frac{\sigma^{2}}{\sigma^{\prime 2}}\right]\,.$ (53)

The next step is to Taylor expand the KL divergence around $\sigma$, by
setting $\sigma^{\prime}=\sigma+\delta\sigma$. For the Taylor expansion we
find

$\displaystyle D_{\mathrm{KL}}(P_{\mu,\,\sigma^{2}}||P_{\mu,\,\sigma^{\prime
2}})$
$\displaystyle=\frac{1}{2}\left[\partial_{\sigma^{\prime}}^{2}D_{\mathrm{KL}}(P_{\mu,\,\sigma^{2}}||P_{\mu,\,\sigma^{\prime
2}})\right]_{\delta\sigma=0}\delta\sigma^{2}+O(\delta\sigma^{3})=-\frac{1}{4}\left[\frac{2}{\sigma^{\prime
2}}-6\frac{\sigma^{2}}{\sigma^{\prime
4}}\right]_{\delta\sigma=0}\delta\sigma^{2}+O(\delta\sigma^{3})$ (54)

From this we find the FIM:

$\displaystyle{\mathcal{I}}(\sigma)=\frac{1}{\sigma^{2}}$ (55)

and Jeffreys prior

$\displaystyle w(\sigma)=\frac{1}{\sigma}\,.$ (56)

Jeffreys prior for the Gaussian is shown in Figure S2. It diverges on the
boundary of parameter space as $\sigma\rightarrow 0$ and approaches zero for
$\sigma\rightarrow\infty$.

Figure S2: Jeffreys prior for the Gaussian model with variance $\sigma^{2}$,
diverges at the boundary of parameter space for $\sigma\rightarrow 0$.

##  Appendix B Standard error and diverging KL divergence

Consider the KL divergence for the Bernoulli model

$\displaystyle
D_{\mathrm{KL}}(P_{\theta^{*}}||P_{\theta})=H_{\theta^{*}}(\theta)-H_{\theta^{*}}(\theta^{*})\,,$
(57)

where $H_{x}(y)$ is the cross-entropy function

$\displaystyle H_{x}(y)=-x\log y-(1-x)\log(1-y)\,.$ (58)

For $\theta\rightarrow 0,1$, the KL divergence diverges to plus infinity. What
is the behaviour of the KL divergence when we try to estimate the value of
$\theta$ using data from many Bernoulli trials generated by the true model
$\theta^{*}$? In particular, how does the KL divergence behave when
$\theta^{*}$ lies very close to a boundary of parameter space?

To answer this question we consider the empirically estimated parameter
$\hat{\theta}$. This parameter is subject to statistical fluctuations due to
finite data size. In general, the standard error for a parameter estimated
from $N$ independent trials is given by

$\displaystyle\sigma_{N}=\frac{\sigma}{\sqrt{N}}\,,$ (59)

where $\sigma$ is the standard deviation of the distribution we are sampling
from. Within the standard error the estimated parameter lies within the bounds

$\displaystyle\hat{\theta}_{\pm}=\theta^{*}\pm\sigma_{N}\,.$ (60)

For the Bernoulli model, the empirical estimate of the probability parameter
is obtained as the mean number of events observed in $N$ independent Bernoulli
trials and $\sigma$ is the standard deviation of the Bernoulli distribution

$\displaystyle\sigma=\sqrt{\theta^{*}(1-\theta^{*})}=\theta^{*}\sqrt{1/\theta^{*}-1}\,.$
(61)

The statistically expected bounds are

$\displaystyle\hat{\theta}_{\pm}=\theta^{*}\left(1-\sqrt{\frac{1/\theta^{*}-1}{N}}\right)\,.$
(62)

Let us now focus on the boundary at $\theta=0$ and recall that the KL
divergence diverges at this boundary. Therefore, one way to check whether we
can expect the KL divergence to diverge for a given model $\theta^{*}$ and
number of observations $N$, is to check whether the boundary $\theta=0$ lies
within the standard error bound $\hat{\theta}_{-}$. From the last equation we
see that the boundary is included if the condition

$\displaystyle\sqrt{\frac{1/\theta^{*}-1}{N}}>1$ (63)

is satisfied. Rearranging this inequality gives

$\displaystyle\theta^{*}<\frac{1}{N+1}\,.$ (64)

Thus true models which are within $1/(N+1)$ of the boundary, the boundary
where the KL divergence diverges lies within the standard error and
statistically we can expect the KL divergence to diverge. See Figure S3 for
examples plots.

For the Poisson model, the empirical estimate of the parameter lambda is
obtained as the mean number of events observed in $N$ independent Poisson
trials and $\sigma$ is the standard deviation of the Poisson distribution

$\displaystyle\sigma_{N}=\sqrt{\frac{\lambda^{*}}{N}}=\lambda^{*}\frac{1}{\sqrt{\lambda^{*}N}}\,.$
(65)

The statistically expected bounds are

$\displaystyle\hat{\lambda}_{\pm}=\lambda^{*}\pm\sigma_{N}=\lambda^{*}\left(1\pm\frac{1}{\sqrt{\lambda^{*}N}}\right)\,.$
(66)

The boundary of parameter space lies at $\lambda=0$. Therefore we consider
$\hat{\lambda}_{-}$ as the important bound and the condition for the boundary
of parameter space to be included is

$\displaystyle\frac{1}{\sqrt{\lambda^{*}N}}>1\,,$ (67)

or

$\displaystyle\lambda^{*}<\frac{1}{N}\,.$ (68)

Therefore, for any given $N$ we can find a models $\lambda^{*}$ which lie
close enough to the boundary according to the above condition, and for which
we should expect the KL divergence to diverge.

Figure S3: Diverging behaviour of
$D_{\mathrm{KL}}(P_{\theta^{*}}||P_{\hat{\theta}_{-}=\theta^{*}-\sigma_{N}})$
for three values of $N$. Vertical gray lines show the position of diverging
behaviour at $1/(N+1)$.

##  Appendix C Optimality condition

###  C.1 The optimality condition for exponential families in mean
parameterisation

This derivation makes use of the fact that the Kullback-Leibler divergence
between two distributions from the same exponential family in canonical form,
can be expressed as the Bregman divergence (defined with respect to the dual
of the log partition function, $F^{*}$) between the mean parameters
$\bm{\theta}$:

$\displaystyle
D_{\mathrm{KL}}(Q_{\bm{\eta}^{*}}||Q_{\bm{\eta}})=B_{F^{*}}(\bm{\theta}^{*}:\bm{\theta})\,.$
(69)

Writing out the definition of the Bregman divergence, the right-hand side is
of the form

$\displaystyle
B_{F^{*}}(\bm{\theta}^{*}:\bm{\theta})=F^{*}(\bm{\theta}^{*})-F^{*}(\bm{\theta})-[\bm{\theta}^{*}-\bm{\theta}]^{T}\cdot\nabla_{\bm{\theta}}F^{*}(\bm{\theta})\,,$
(70)

where

$\displaystyle\bm{\theta}=\kappa\bm{\theta}_{\mathrm{bias}}+(1-\kappa)\hat{\bm{\theta}}_{\mathrm{ML}}\,.$
(71)

To compute the condition for the optimal model, we take the derivative of the
Bregman divergence with respect to the bias:

$\displaystyle\frac{1}{\kappa}\nabla_{\theta_{\mathrm{bias}\,i}}B_{F^{*}}(\bm{\theta}^{*}:\bm{\theta})=\nabla_{\theta_{i}}B_{F^{*}}(\bm{\theta}^{*}:\bm{\theta})$
$\displaystyle=-\nabla_{\theta_{i}}F^{*}(\bm{\theta})-\bm{\theta}^{*T}\cdot\nabla_{\theta_{i}}\nabla_{\bm{\theta}}F^{*}(\bm{\theta})+\nabla_{\theta_{i}}\\{\bm{\theta}^{T}\cdot\nabla_{\bm{\theta}}F^{*}(\bm{\theta})\\}$
(72)
$\displaystyle=-\nabla_{\theta_{i}}F^{*}(\bm{\theta})-\sum_{j=1}^{d}\theta^{*}_{j}\nabla_{\theta_{i}}\nabla_{\theta_{j}}F^{*}(\bm{\theta})+\nabla_{\theta_{i}}F^{*}(\bm{\theta})+\sum_{j=1}^{d}\theta_{j}\nabla_{\theta_{i}}\nabla_{\theta_{j}}F^{*}(\bm{\theta})$
(73)
$\displaystyle=-\sum_{j=1}^{d}[\theta^{*}_{j}-\theta_{j}]\nabla_{\theta_{i}}\nabla_{\theta_{j}}F^{*}(\bm{\theta})\,.$
(74)

Next, we use the relation

$\displaystyle\mathcal{I}_{ij}(\bm{\theta})=\nabla_{\theta_{i}}\nabla_{\theta_{j}}F^{*}(\bm{\theta})\,,$
(75)

relating the Fisher information metric in mean parameterization to the dual
function. We derive this relation starting with Eq. (8)

$\displaystyle
D_{\mathrm{KL}}(Q_{\bm{\eta}^{*}(\bm{\theta}^{*})}||Q_{\bm{\eta}(\bm{\theta})})=B_{F^{*}}({\bm{\theta}^{*}}:\bm{\theta})\,,$
(76)

replacing $\bm{\theta}=\bm{\theta}^{*}-\delta\bm{\theta}$, and Taylor
expanding in $\delta\bm{\theta}$. Using the definition of the Bregman
divergence Eq. (70) and noting that
$\nabla_{\delta\theta_{i}}=-\nabla_{\theta_{i}}$, the coefficients of the
first three orders in $\delta\bm{\theta}$ are given by

$\displaystyle
B_{F^{*}}({\bm{\theta}^{*}}:\bm{\theta}^{*}-\delta\bm{\theta})\big{|}_{\delta\bm{\theta}=0}$
$\displaystyle=0\,,$ (77)
$\displaystyle\nabla_{\delta\theta_{i}}B_{F^{*}}({\bm{\theta}^{*}}:\bm{\theta}^{*}-\delta\bm{\theta})\big{|}_{\delta\bm{\theta}=0}$
$\displaystyle=\left[\nabla_{\theta_{i}}F^{*}(\bm{\theta})-\nabla_{\theta_{i}}F^{*}(\bm{\theta})+\sum_{j}\delta\theta_{j}\nabla_{\theta_{i}}\nabla_{\theta)j}F^{*}(\bm{\theta})\right]_{\delta\bm{\theta}=0}=0\,,$
(78)
$\displaystyle\nabla_{\delta\theta_{k}}\nabla_{\delta\theta_{i}}B_{F^{*}}({\bm{\theta}^{*}}:\bm{\theta}^{*}-\delta\bm{\theta})\big{|}_{\delta\bm{\theta}=0}$
$\displaystyle=\Big{[}-\nabla_{\delta\theta_{k}}\nabla_{\theta_{i}}F^{*}(\bm{\theta})+\nabla_{\delta\theta_{k}}\nabla_{\theta_{i}}F^{*}(\bm{\theta})+\delta_{kj}\nabla_{\theta_{i}}\nabla_{\theta_{j}}F^{*}(\bm{\theta})$
(79)
$\displaystyle-\sum_{j}\delta\theta_{j}\nabla_{\delta\theta_{k}}\nabla_{\theta_{i}}\nabla_{\theta_{j}}F^{*}(\bm{\theta})\Big{]}_{\delta\bm{\theta}=0}=\nabla_{\theta_{i}}\nabla_{\theta_{k}}F^{*}(\bm{\theta}^{*})\,.$
(80)

Using this relation, we rewrite the above condition as

$\displaystyle\frac{1}{\kappa}\nabla_{\bm{\theta}_{\mathrm{bias}}}B_{F^{*}}(\bm{\theta}^{*}:\bm{\theta})=-\mathcal{I}(\bm{\theta})\left[\bm{\theta}^{*}-\bm{\theta}\right]=\mathcal{I}(\bm{\theta})\delta\bm{\theta}\,,$
(81)

where $-\delta\bm{\theta}\equiv\bm{\theta}^{*}-\bm{\theta}$. Thus, in mean
parameterisation, the optimality constraint in Eq. (6) of the main text, takes
the simple form

$\displaystyle\frac{1}{\kappa}\nabla_{\bm{\theta}_{\mathrm{bias}}}\big{<}\mathcal{L}_{\bm{\theta}^{*}}(\bm{\theta})\big{>}_{\hat{\bm{\theta}}_{\mathrm{ML}}}=\big{<}\mathcal{I}(\bm{\theta})\delta\bm{\theta}\big{>}_{\hat{\bm{\theta}}_{\mathrm{ML}}}=0\,.$
(82)

###  C.2 The optimality condition for the Bernoulli model in mean
parameterization derived from Taylor expansion

As an example, we evaluate the optimality condition

$\displaystyle\nabla_{\bm{\theta}_{\mathrm{bias}}}\big{<}D_{\mathrm{KL}}(P_{\bm{\theta}^{*}}||P_{\bm{\theta}})\big{>}_{\hat{\bm{\theta}}_{\mathrm{ML}}}=0$
(83)

for the Bernoulli model written in mean parameters. In this derivation we will
use Taylor expansion of the KL divergence and find that all but a single term
of the infinite series cancel to yield a simplified form of the condition in
Eq. (9) of main text.

We consider KL divergence between two Bernoulli distributions, denoted by
$P_{\bm{\theta}^{*}}$ and $P_{\bm{\theta}}$, respectively:

$\displaystyle D_{\mathrm{KL}}(P_{{\bm{\theta}}^{*}}||P_{{\bm{\theta}}})$
$\displaystyle=D_{\mathrm{KL}}(P_{{\bm{\theta}}+\delta{\bm{\theta}}}||P_{{\bm{\theta}}})=\sum_{x}P_{{\bm{\theta}}+\delta{\bm{\theta}}}\log\left(\frac{P_{{\bm{\theta}}+\delta{\bm{\theta}}}}{P_{{\bm{\theta}}}}\right)\,,$
(84)

where ${\bm{\theta}}^{*}\equiv{\bm{\theta}}+\delta{\bm{\theta}}$ and expand
the right-hand side in $\delta{\bm{\theta}}$. For the first few orders we find

$\displaystyle
D_{\mathrm{KL}}(P_{{\bm{\theta}}^{*}}||P_{{\bm{\theta}}})\Big{|}_{\delta{\bm{\theta}}=0}$
$\displaystyle=0\,,$ (85) $\displaystyle\partial
D_{\mathrm{KL}}(P_{{\bm{\theta}}^{*}}||P_{{\bm{\theta}}})\Big{|}_{\delta{\bm{\theta}}=0}$
$\displaystyle=1\,,$ (86)
$\displaystyle\partial^{2}D_{\mathrm{KL}}(P_{{\bm{\theta}}^{*}}||P_{{\bm{\theta}}})\Big{|}_{\delta{\bm{\theta}}=0}$
$\displaystyle=1/{\bm{\theta}}\,,$ (87)
$\displaystyle\partial^{3}D_{\mathrm{KL}}(P_{{\bm{\theta}}^{*}}||P_{{\bm{\theta}}})\Big{|}_{\delta{\bm{\theta}}=0}$
$\displaystyle=-1/\bm{\theta}^{2}\,,$ (88)
$\displaystyle\partial^{4}D_{\mathrm{KL}}(P_{{\bm{\theta}}^{*}}||P_{{\bm{\theta}}})\Big{|}_{\delta{\bm{\theta}}=0}$
$\displaystyle=2/\bm{\theta}^{3}\,,$ (89)
$\displaystyle\partial^{5}D_{\mathrm{KL}}(P_{{\bm{\theta}}^{*}}||P_{{\bm{\theta}}})\Big{|}_{\delta{\bm{\theta}}=0}$
$\displaystyle=-6/\bm{\theta}^{4}\,,$ (90) $\displaystyle\vdots$

where $\partial\equiv\partial_{\delta{\bm{\theta}}}$ and powers on the right-
hand side act element-wise. The Taylor expansion is thus

$\displaystyle D_{\mathrm{KL}}(P_{\bm{\theta}^{*}}||P_{\bm{\theta}})$
$\displaystyle=\sum_{i}\left\\{\delta\theta_{i}+\frac{1}{2!}\frac{\delta\theta_{i}^{2}}{\theta_{i}}-\frac{1}{3!}\frac{\delta\theta_{i}^{3}}{\theta_{i}^{2}}+\frac{1}{4!}\frac{2\delta\theta_{i}^{4}}{\theta_{i}^{3}}-\frac{1}{5!}\frac{6\delta\theta_{i}^{5}}{\theta_{i}^{4}}+O(\delta_{i}^{6})\right\\}$
(91)
$\displaystyle=1+\sum_{i}\left\\{\frac{1}{2!}\frac{\delta\theta_{i}^{2}}{\theta_{i}}-\frac{1}{3!}\frac{\delta\theta_{i}^{3}}{\theta_{i}^{2}}+\frac{1}{4!}\frac{2\delta\theta_{i}^{4}}{\theta_{i}^{3}}-\frac{1}{5!}\frac{6\delta\theta_{i}^{5}}{\theta_{i}^{4}}+O(\delta_{i}^{6})\right\\}$
(92)

Next, we act with the derivative with respect to the bias:

$\displaystyle\partial_{\theta_{\mathrm{bias}\,i}}D_{\mathrm{KL}}(P_{\bm{\theta}^{*}}||P_{\bm{\theta}})=(\partial\theta_{i}/\partial{\theta_{\mathrm{bias}\,i}})\partial_{\theta_{i}}D_{\mathrm{KL}}(P_{\bm{\theta}^{*}}||P_{\bm{\theta}})=\kappa\partial_{\theta_{i}}D_{\mathrm{KL}}(P_{\bm{\theta}^{*}}||P_{\bm{\theta}})\,.$
(93)

Remembering that $\delta\bm{\theta}$ is a function of $\bm{\theta}$,
$\partial\delta\theta/\partial\theta=-1$, we obtain:

$\displaystyle\partial_{\theta_{i}}D_{\mathrm{KL}}(P_{\bm{\theta}^{*}}||P_{\bm{\theta}})$
$\displaystyle=-\frac{1}{\theta_{i}}\delta\theta_{i}+\left[-\frac{1}{2!}\frac{\delta\theta_{i}^{2}}{\theta_{i}^{2}}+\frac{3}{3!}\frac{\delta\theta_{i}^{2}}{\theta_{i}^{2}}\right]+\left[\frac{2}{3!}\frac{\delta\theta_{i}^{3}}{\theta_{i}^{3}}-\frac{4}{4!}\frac{2\delta\theta_{i}^{3}}{\theta_{i}^{3}}\right]+...$
(94) $\displaystyle=-\frac{1}{\theta_{i}}\delta\theta_{i}\,,$ (95)

where on in the first line on the right-hand side we have grouped terms such
that expressions in the brackets vanish to show how terms of successive order
in $\delta\theta$ cancel. Using $\theta_{2}=1-\theta_{1}$ , we finally see
that in matrix notation the expression takes the form

$\displaystyle\partial_{\bm{\theta}}D_{\mathrm{KL}}(P_{\bm{\theta}^{*}}||P_{\bm{\theta}})=-\begin{pmatrix}1/\theta_{1}&0\\\
0&1/(1-\theta_{1})\end{pmatrix}\delta\bm{\theta}=-\mathcal{I}(\bm{\theta})\delta\bm{\theta}\,,$
(96)

where we recognised the $2\times 2$ matrix to be the Fisher information
matrix. Thus, we obtain the optimality condition

$\displaystyle-\big{<}\partial_{\bm{\theta}}D_{\mathrm{KL}}(P_{\bm{\theta}^{*}}||P_{\bm{\theta}})\big{>}_{\hat{\bm{\theta}}_{\mathrm{ML}}}=\big{<}\mathcal{I}(\bm{\theta})\delta\bm{\theta}\big{>}_{\hat{\bm{\theta}}_{\mathrm{ML}}}=0\,.$
(97)

##  Appendix D Optimal parameter calculation for Poisson distribution with
unknown rate parameter

Figure S4: Optimization of the expected loss for the Poisson distribution with
respect to the prior shows a bias to low complexity models. The true model
lies at $\lambda^{*}=3$. (A) The bias $\lambda_{\mathrm{bias}}$ is optimised
while $\kappa=0.2$ ($q=0.25$) is kept fixed. The plot shows the landscape of
the expected loss as a function of the bias and the number of observations. We
are using a non-linear colormap based on the cumulative distribution of the
loss. (B) The plot shows for different true models $\lambda^{*}$ the ratio of
the local complexity measure by Jeffreys prior of the expected optimal to the
true model.

We consider two Poisson distributions with parameters $\lambda^{*}$ and
$\lambda$:

$\displaystyle
P_{\lambda^{*}}(n)=\frac{\lambda^{*\,n}e^{-\lambda^{*}}}{n!}\,,\;P_{\lambda}(n)=\frac{\lambda^{n}e^{-\lambda}}{n!}\,.$
(98)

The KL divergence between the two distributions is given by

$\displaystyle D_{\mathrm{KL}}(P_{\lambda^{*}}||Q_{\lambda})$
$\displaystyle=\sum_{n=0}^{\infty}\frac{\lambda^{*\,n}e^{-\lambda^{*}}}{n!}\log\left(\frac{\lambda^{*\,n}e^{-\lambda^{*}}}{n!}\frac{n!}{\lambda^{n}e^{-\lambda}}\right)=\lambda^{*}\log\frac{\lambda^{*}}{\lambda}+(\lambda-\lambda^{*})\,.$
(99)

We want to minimise the posterior mean loss

$\displaystyle\big{<}\mathcal{L}_{\lambda^{*}}\big{>}_{\mathrm{post}}(\lambda)=\int\mathrm{d}\tilde{\lambda}P(\tilde{\lambda}|\mathbf{n})D_{\mathrm{KL}}(P_{\tilde{\lambda}}||Q_{\lambda})$
(100)

with respect to $\lambda$. Taking the derivative with respect to $\lambda$ we
obtain the condition

$\displaystyle\partial_{\lambda}\big{<}\mathcal{L}_{\lambda^{*}}\big{>}_{\mathrm{post}}=\partial_{\lambda}\int\mathrm{d}\tilde{\lambda}P(\tilde{\lambda}|\mathbf{n})D_{\mathrm{KL}}(P_{\tilde{\lambda}}||Q_{\lambda})=\int\mathrm{d}\tilde{\lambda}P(\tilde{\lambda}|\mathbf{n})\left(-\tilde{\lambda}/\lambda+1\right)=0\,.$
(101)

The condition is satisfied by

$\displaystyle\lambda=\int\mathrm{d}\tilde{\lambda}P(\tilde{\lambda}|\mathbf{n})\tilde{\lambda}\,.$
(102)

The conjugate prior for the mean number of counts $\lambda$ is the Gamma
distribution and the posterior takes the form

$\displaystyle
P(\lambda|\mathbf{n})=\mathrm{Gamma}(\alpha^{\prime},\beta^{\prime})=\frac{{\beta^{\prime}}^{\alpha^{\prime}}}{\Gamma(\alpha^{\prime})}\lambda^{\alpha^{\prime}-1}e^{-\beta^{\prime}\lambda}$
(103)

with

$\displaystyle\alpha^{\prime}$
$\displaystyle=\alpha+N\bar{X}=\alpha+N\frac{1}{N}\sum_{i=1}^{N}x_{i}\,,$
(104) $\displaystyle\beta^{\prime}$ $\displaystyle=\beta+N\,,$ (105)

and the hyperparameters $\alpha,\beta>0$. $\bar{X}$ is an empirical estimate
of the rate $\lambda$ from $N$ observations. The posterior mean which gives
the optimal parameter choice is given by

$\displaystyle\lambda(\hat{\lambda}_{\mathrm{ML}};\alpha,\beta)=\frac{\alpha^{\prime}}{\beta^{\prime}}=\frac{\alpha+N\hat{\lambda}_{\mathrm{ML}}}{\beta+N}=\frac{\frac{\beta}{N}\frac{\alpha}{\beta}+\hat{\lambda}_{\mathrm{ML}}}{\frac{\beta}{N}+1}=\frac{q\lambda_{\mathrm{bias}}+\hat{\lambda}_{\mathrm{ML}}}{q+1}$
(106)

and we define as before

$\displaystyle
N_{0}\equiv\beta\,,\quad\lambda_{\mathrm{bias}}\equiv\frac{\alpha}{\beta}\,,\quad\text{and
}\quad q\equiv\frac{\beta}{N}\,.$ (107)

The hyperparameter $\beta$ can thus be interpreted as virtual observations and
$\alpha$ as the sum of counts across $\beta$ observations. If we do not trust
in the empirical estimate, then we should take $\alpha$ large compared to
$\beta$. In this limit, the optimal $\lambda$ values becomes large and lies
further away from the boundary of parameter space situated at $\lambda=0$.

Up to a scaling by $1/N$, the maximum-likelihood value of the average count
$\hat{\lambda}_{\mathrm{ML}}$ computed from $N$ observations, is distributed
according to the true Poisson distribution with the rate rescaled by $N$. This
can be seen from the characteristic function of the Poisson distribution given
by $\exp(\lambda(e^{it}-1))$ and the fact that the characteristic function of
the sum of $N$ i.i.d. random variables is given by the $N$’th power of the
characteristic function of a single random variable. The empirically averaged
expected loss is then given by

$\displaystyle\big{<}\mathcal{L}_{\lambda^{*}}(\lambda)\big{>}_{\hat{\lambda}_{\mathrm{ML}}}=\sum_{n=0}^{\infty}\,\mathrm{Poi}(n;N\lambda^{*})D_{\mathrm{KL}}(P_{\lambda^{*}}||Q_{{\lambda}(n/N;\alpha,\beta)})\,.$
(108)

We minimise this loss with respect to the hyperparameters, for instance by
keeping $\beta$ fixed an minimising with respect to $\alpha$. In Figure S4 we
show the optimal choice of the hyperparameter $\alpha$ as a function of the
data size, confirming that for small $N$, $\alpha$ is large, while for
increasing data size $\alpha\rightarrow 0$.

##  Appendix E Optimal parameter calculation for Gaussian with known mean and
unknown variance

Consider two Gaussian distributions with the same known mean $\mu$ and unknown
variances ${\sigma^{*}}^{2}$ and $\sigma^{2}$

$\displaystyle
P_{\sigma^{*2}}(x)=\frac{1}{\sqrt{2\pi{\sigma^{*}}^{2}}}e^{-\frac{(x-\mu)^{2}}{2{\sigma^{*}}^{2}}}\,,\;P_{\sigma^{2}}(x)=\frac{1}{\sqrt{2\pi\sigma^{2}}}e^{-\frac{(x-\mu)^{2}}{2\sigma^{2}}}\,.$
(109)

The KL divergence two between these distributions is given by

$\displaystyle
D_{\mathrm{KL}}(P_{\sigma^{*2}}||P_{\sigma^{2}})=-\frac{1}{2}\left[\log{\sigma^{*}}^{2}-\log\sigma^{2}+1-\frac{{\sigma^{*}}^{2}}{\sigma^{2}}\right]\,.$
(110)

We want to minimise the posterior mean loss

$\displaystyle\big{<}\mathcal{L}_{\sigma^{*}}\big{>}_{\mathrm{post}}(\sigma)=\int\mathrm{d}{\tilde{\sigma}}^{2}P({\tilde{\sigma}}^{2}|\mathbf{n})D_{\mathrm{KL}}(P_{\tilde{\sigma}^{2}}||P_{\sigma^{2}})$
(111)

with respect to $\sigma^{2}$, where $P({\tilde{\sigma}}^{2}|\mathbf{n})$ is
the posterior. Taking the derivative with respect to $\sigma$ we obtain the
condition

$\displaystyle\partial_{\sigma}\int\mathrm{d}{\tilde{\sigma}}^{2}P({\tilde{\sigma}}^{2}|\mathbf{n})D_{\mathrm{KL}}(P_{\tilde{\sigma}^{2}}||P_{\sigma^{2}})=\int\mathrm{d}{\tilde{\sigma}}^{2}P({\tilde{\sigma}}^{2}|\mathbf{n})\left(\frac{1}{\sigma}-\frac{{\tilde{\sigma}}^{2}}{\sigma^{3}}\right)=0$
(112)

satisfied by

$\displaystyle\sigma^{2}=\int\mathrm{d}{\tilde{\sigma}}^{2}{\tilde{\sigma}}^{2}P({\tilde{\sigma}}^{2}|\mathbf{n})\,.$
(113)

Figure S5: Optimization of the expected loss for the Gaussian distribution
with respect to the prior shows a bias to low complexity models. The true
model lies at $\sigma^{*2}=2$ and we set $\mu=0$ for the known mean value. (A)
The bias $\sigma_{\mathrm{bias}}^{2}$ is optimised while $\kappa=0.2$
($q=0.25$) is kept fixed. The plot shows the landscape of the expected loss as
a function of the bias and the number of observations. We are using a non-
linear colormap based on the cumulative distribution of the loss. (B) The plot
shows for different true models $\sigma^{*2}$ the ratio of the local
complexity measure by Jeffreys prior of the expected optimal to the true
model.

The conjugate prior for the variance is the scaled inverse-chi-squared
distribution $\chi^{-2}$. In [56], the posterior distribution for the
variance, is given as

$\displaystyle
P(\sigma^{2}|\nu_{N},\sigma_{N}^{2})=\chi^{-2}(\sigma^{2}|\nu_{N},\sigma_{N}^{2})=\frac{1}{\Gamma(\nu_{N}/2)}\left(\frac{\nu_{N}{\sigma}^{2}_{N}}{2}\right)^{\nu_{N}/2}x^{-\nu_{N}/2-1}e^{-\frac{\nu_{N}\sigma_{N}^{2}}{2x}}\,,$
(114)

where

$\displaystyle\nu_{N}$ $\displaystyle=\nu+N\,,\quad\text{with}\quad\nu>2\,,$
(115) $\displaystyle\sigma_{N}^{2}$
$\displaystyle=\frac{1}{\nu_{N}}\left(\nu\sigma_{0}^{2}+\sum_{i=1}^{N}(x_{i}-\mu)^{2}\right)\,.$
(116)

The loss-minimizing value for $\sigma^{2}$ is given by the Bayesian mean which
comes out as

$\displaystyle\sigma^{2}(\hat{\sigma}_{\mathrm{ML}}^{2};\nu,\sigma_{0}^{2})=\int\mathrm{d}{\tilde{\sigma}}^{2}{\tilde{\sigma}}^{2}\chi^{-2}(\tilde{\sigma}^{2}|\nu_{N},\sigma_{N}^{2})$
$\displaystyle=\frac{\nu_{N}}{\nu_{N}-2}\sigma_{N}^{2}=\frac{\nu\sigma_{0}^{2}+\sum_{i}(x_{i}-\mu)^{2}}{\nu-2+N}=\frac{\frac{\nu-2}{N}\frac{\nu\sigma_{0}^{2}}{\nu-2}+\frac{1}{N}\sum_{i}(x_{i}-\mu)^{2}}{\frac{\nu-2}{N}+1}$
(117)

and we define

$\displaystyle N_{0}$ $\displaystyle\equiv\nu-2\,,$ (118)
$\displaystyle\sigma_{\mathrm{bias}}^{2}$
$\displaystyle\equiv\frac{\nu\sigma_{0}^{2}}{\nu-2}\,,$ (119)
$\displaystyle\hat{\sigma}_{\mathrm{ML}}^{2}$
$\displaystyle\equiv\frac{1}{N}\sum_{i}(x_{i}-\mu)^{2}\,.$ (120)

Up to a scaling of $\sigma^{*2}/N$, the empirical variance is distributed
according to the chi-squared distribution of order $N$ (recall, the variance
is given by the sum of $N$ i.i.d. squared Gaussian random variables divided by
$N$) given by

$\displaystyle\chi^{2}_{N}(x)=\frac{1}{2^{N/2}\Gamma(N/2)}\;x^{N/2-1}e^{-x/2}\,.$
(121)

Finally, we optimise the expected loss

$\displaystyle\big{<}\mathcal{L}_{\sigma^{*}}(\sigma)\big{>}_{\hat{\sigma}_{\mathrm{ML}}}=\int\mathrm{d}x\,\chi^{2}_{N}(x)D_{\mathrm{KL}}(P_{{\sigma^{*}}^{2}}||P_{{\sigma}^{2}(\sigma^{*2}x/N;\nu,\sigma_{0}^{2})})$
(122)

with respect to the bias term, while $\kappa$ is fixed. In this equation we
are using the continuous version of the Kullback-Leibler divergence. By Eq.
(18) of the main text, $\kappa$ fixes $q$ which in turn provides a value for
$N_{0}$ given $N$:

$\displaystyle q=\frac{N_{0}}{N}=\frac{\kappa}{1-\kappa}\,,\quad
N_{0}=Nq=\nu-2\,,$ (123)

where we have also used the definition of $N_{0}$ in terms of the
hyperparameter $\nu$. For the hyperparameters $\nu$ and the bias we can thus
write:

$\displaystyle\nu$ $\displaystyle=N_{0}+2\,,$ (124)
$\displaystyle\sigma_{\mathrm{bias}}^{2}$
$\displaystyle=\frac{N_{0}+2}{N_{0}}\sigma_{0}^{2}\,.$ (125)

Minimising the expected loss given $N$, yields the optimal choice for
$\sigma_{\mathrm{bias}}^{2}$. In Figure S5 we show the optimal bias value as a
function of the data size.

##  Appendix F Statistical fluctuations of different learning strategies

Figure S6: The prior choice in Bayesian inference controls the safety of the
learning strategy. Left column: long-term growth rate as a function of the
adopted model $\phi$ for three different values of the true model
$\theta^{*}=0.4$, $0.2$ and $0.05$. Shaded boxes indicate the regions of non-
negative growth rate. Right column: Magnitude of fluctuations of the
conditional optimal model as a function of the number of observations $N$ for
two different learning strategies: maximum-likelihood estimation ($\alpha=0$;
dashed line) and Jeffreys prior ($\alpha=1/2$; solid line).

In Figure S6, we show the magnitude of statistical fluctuations of the
inferred model in relation to the region of non-negative long-term growth
rates. Subfigures A, B and C show different true models $\theta^{*}=0.4$,
$0.2$ and $0.05$, respectively. In the left column of subfigures we show the
growth rates $G_{\theta^{*}}(\theta)$ as a function of the model choice
$\theta$. Regions where the growth rate is non-negative are indicated by color
shaded boxes. In the right column of subfigures we show the statistical
fluctuations of the inferred model as function of $N$ for two different
learning strategies. The first learning strategy is maximum-likelihood
estimation $\alpha=0$ and the second strategy uses $\alpha=1/2$. The later
strategy ($\alpha=1/2$) is safer, because the magnitude of fluctuations is
shrunk and more contained in the interval of non-negative growth rates in
comparison to the maximum-likelihood strategy.

## References

  * Burnet [1976] F. M. Burnet, A modification of Jerne’s theory of antibody production using the concept of clonal selection, CA: A Cancer Journal for Clinicians 26, 119 (1976).
  * Kussell and Leibler [2005] E. Kussell and S. Leibler, Phenotypic diversity, population growth, and information in fluctuating environments, Science 309, 2075 (2005).
  * Xue and Leibler [2018] B. Xue and S. Leibler, Benefits of phenotypic plasticity for population growth in varying environments, Proceedings of the National Academy of Sciences 115, 12745 (2018).
  * Hopfield [1974] J. J. Hopfield, Kinetic proofreading: A new mechanism for reducing errors in biosynthetic processes requiring high specificity, Proceedings of the National Academy of Sciences 71, 4135 (1974).
  * Martincorena  _et al._ [2012] I. Martincorena, A. S. Seshasayee, and N. M. Luscombe, Evidence of non-random mutation rates suggests an evolutionary risk management strategy, Nature 485, 95 (2012).
  * Harder and Real [1987] L. D. Harder and L. A. Real, Why are bumble bees risk averse?, Ecology 68, 1104 (1987).
  * Bhat  _et al._ [2020] U. Bhat, C. P. Kempes, and J. D. Yeakel, Scaling the risk landscape drives optimal life-history strategies and the evolution of grazing, Proceedings of the National Academy of Sciences 117, 1580 (2020).
  * Sutton and Barto [2018] R. S. Sutton and A. G. Barto,  _Reinforcement learning: An introduction_ (MIT press, 2018).
  * Vergassola  _et al._ [2007] M. Vergassola, E. Villermaux, and B. I. Shraiman, ‘infotaxis’ as a strategy for searching without gradients, Nature 445, 406 (2007).
  * Cohen [1966] D. Cohen, Optimizing reproduction in a randomly varying environment, Journal of Theoretical Biology 12, 119 (1966).
  * Rubenstein [2011] D. R. Rubenstein, Spatiotemporal environmental variation, risk aversion, and the evolution of cooperative breeding as a bet-hedging strategy, Proceedings of the National Academy of Sciences 108, 10816 (2011).
  * Levins [1962] R. Levins, Theory of fitness in a heterogeneous environment. i. the fitness set and adaptive function, The American Naturalist 96, 361 (1962).
  * Nemenman [2012] I. Nemenman, Information theory and adaptation, Quantitative biology: from molecular to cellular systems 4, 73 (2012).
  * Perkins and Swain [2009] T. J. Perkins and P. S. Swain, Strategies for cellular decision-making, Molecular systems biology 5, 326 (2009).
  * Laplace, marquis de [1814] P. S. Laplace, marquis de,  _Essai philosophique sur les probabilités_ (Courcier, 1814).
  * Keynes [1921] J. M. Keynes,  _A Treatise on Probability_ (Dover Publications, 1921).
  * Kullback [1997] S. Kullback,  _Information theory and statistics_ (Courier Corporation, 1997).
  * Jaynes [2003] E. T. Jaynes,  _Probability Theory: The Logic of Science_, edited by G. L. Bretthorst (Cambridge University Press, 2003).
  * Jaynes [1957a] E. T. Jaynes, Information theory and statistical mechanics, Phys. Rev. 106, 620 (1957a).
  * Jaynes [1957b] E. T. Jaynes, Information theory and statistical mechanics. ii, Phys. Rev. 108, 171 (1957b).
  * Rissanen [1978] J. Rissanen, Modeling by shortest data description, Automatica 14, 465 (1978).
  * Balasubramanian [1997] V. Balasubramanian, Statistical inference, occam’s razor, and statistical mechanics on the space of probability distributions, Neural Comput. 9, 349–368 (1997).
  * Jeffreys [1946] H. Jeffreys, An invariant form for the prior probability in estimation problems, Proceedings of the Royal Society of London. Series A. Mathematical and Physical Sciences 186, 453 (1946).
  * Kelly Jr. [1956] J. L. Kelly Jr., A new interpretation of information rate, Bell System Technical Journal 35, 917 (1956).
  * Rivoire and Leibler [2011] O. Rivoire and S. Leibler, The value of information for populations in varying environments, Journal of Statistical Physics 142, 1124 (2011).
  * Fisher and Russell [1922] R. A. Fisher and E. J. Russell, On the mathematical foundations of theoretical statistics, Philosophical Transactions of the Royal Society of London. Series A, Containing Papers of a Mathematical or Physical Character 222, 309 (1922).
  * Amari [1983] S.-I. Amari, A foundation of information geometry, Electronics and Communications in Japan (Part I: Communications) 66, 1 (1983).
  * Nielsen and Garcia [2011] F. Nielsen and V. Garcia, Statistical exponential families: A digest with flash cards (2011), arXiv:0911.4863 [cs.LG] .
  * Wainwright  _et al._ [2008] M. J. Wainwright, M. I. Jordan, _et al._ , Graphical models, exponential families, and variational inference, Foundations and Trends in Machine Learning 1, 1 (2008).
  * Bregman [1967] L. M. Bregman, The relaxation method of finding the common point of convex sets and its application to the solution of problems in convex programming, USSR computational mathematics and mathematical physics 7, 200 (1967).
  * Nielsen [2020] F. Nielsen, An elementary introduction to information geometry, Entropy 22, 1100 (2020).
  * Banerjee  _et al._ [2005] A. Banerjee, S. Merugu, I. S. Dhillon, J. Ghosh, and J. Lafferty, Clustering with Bregman Divergences., Journal of machine learning research 6 (2005).
  * Diaconis and Ylvisaker [1979] P. Diaconis and D. Ylvisaker, Conjugate priors for exponential families, The Annals of Statistics 7, 269 (1979).
  * Lidstone [1920] G. J. Lidstone, Note on the general case of the bayes-laplace formula for inductive or a posteriori probabilities, Transactions of the Faculty of Actuaries 8, 13 (1920).
  * Nemenman  _et al._ [2001] I. Nemenman, F. Shafee, and W. Bialek, Entropy and inference, revisited, in  _Proceedings of the 14th International Conference on Neural Information Processing Systems: Natural and Synthetic_ , NIPS’01 (MIT Press, Cambridge, MA, USA, 2001) p. 471–478.
  * Mayer  _et al._ [2015] A. Mayer, V. Balasubramanian, T. Mora, and A. M. Walczak, How a well-adapted immune system is organized, Proceedings of the National Academy of Sciences 112, 5950 (2015).
  * Gremer and Venable [2014] J. R. Gremer and D. L. Venable, Bet hedging in desert winter annual plants: optimal germination strategies in a variable environment., Ecology letters 17 3, 380 (2014).
  * Balaban  _et al._ [2004] N. Balaban, J. Merrin, R. Chait, L. Kowalik, and S. Leibler, Bacterial persistence as a phenotypic switch, Science 305, 1622 (2004).
  * Moyed and Bertrand [1983] H. S. Moyed and K. P. Bertrand, hipa, a newly recognized gene of escherichia coli k-12 that affects frequency of persistence after inhibition of murein synthesis., Journal of Bacteriology 155, 768 (1983).
  * Mattingly  _et al._ [2021] H. Mattingly, K. Kamino, B. Machta, and T. Emonet, E. coli chemotaxis is information-limited, bioRxiv 10.1101/2021.02.22.432091 (2021).
  * Mattingly  _et al._ [2018] H. H. Mattingly, M. K. Transtrum, M. C. Abbott, and B. B. Machta, Maximizing the information learned from finite data selects a simple model, Proceedings of the National Academy of Sciences 115, 1760 (2018).
  * Quinn  _et al._ [2021] K. N. Quinn, M. C. Abbott, M. K. Transtrum, B. B. Machta, and J. P. Sethna, Information geometry for multiparameter models: New perspectives on the origin of simplicity (2021), arXiv:2111.07176 [cond-mat.stat-mech] .
  * Gutenkunst  _et al._ [2007] R. N. Gutenkunst, J. J. Waterfall, F. P. Casey, K. S. Brown, C. R. Myers, and J. P. Sethna, Universally sloppy parameter sensitivities in systems biology models, PLoS computational biology 3, e189 (2007).
  * Schnaack  _et al._ [2021] O. H. Schnaack, L. Peliti, and A. Nourmohammad, Risk-utility tradeoff shapes memory strategies for evolving patterns, arXiv preprint arXiv:2110.15008 (2021).
  * Lee  _et al._ [2022] E. D. Lee, J. C. Flack, and D. C. Krakauer, Outsourcing memory through niche construction, bioRxiv , 2022 (2022).
  * Feldman [2016] J. Feldman, The simplicity principle in perception and cognition, Wiley Interdisciplinary Reviews: Cognitive Science 7, 330 (2016).
  * Koffka [2013] K. Koffka,  _Principles of Gestalt psychology_ , Vol. 44 (routledge, 2013).
  * Chater and Vitányi [2003] N. Chater and P. Vitányi, Simplicity: a unifying principle in cognitive science?, Trends in cognitive sciences 7, 19 (2003).
  * Gershman and Niv [2013] S. J. Gershman and Y. Niv, Perceptual estimation obeys occam’s razor, Frontiers in psychology 4, 623 (2013).
  * Tavoni  _et al._ [2022] G. Tavoni, T. Doi, C. Pizzica, V. Balasubramanian, and J. I. Gold, Human inference reflects a normative balance of complexity and accuracy, Nature Human Behaviour 6, 1153 (2022).
  * Piasini  _et al._ [2023] E. Piasini, S. Liu, P. Chaudhari, V. Balasubramanian, and J. I. Gold, How occam’s razor guides human decision-making, bioRxiv , 2023 (2023).
  * Rissanen [1996] J. J. Rissanen, Fisher information and stochastic complexity, IEEE transactions on information theory 42, 40 (1996).
  * Grünwald [2007] P. D. Grünwald,  _The minimum description length principle_ (MIT press, 2007).
  * Bialek  _et al._ [2001] W. Bialek, I. Nemenman, and N. Tishby, Predictability, complexity, and learning, Neural computation 13, 2409 (2001).
  * Sachdeva  _et al._ [2021] V. Sachdeva, T. Mora, A. M. Walczak, and S. E. Palmer, Optimal prediction with resource constraints using the information bottleneck, PLOS Computational Biology 17, e1008743 (2021).
  * Murphy [2007] K. P. Murphy, Conjugate bayesian analysis of the gaussian distribution, https://www.cs.ubc.ca/~murphyk/Papers/bayesGauss.pdf (2007), (version: 2023-04-13).