Learning in Volatile Environments With the Bayes Factor Surprise

Animals, humans, and similarly reinforcement learning agents may safely assume that the world is stochastic and stationary during some intervals of time interrupted by change points. The position of leaves on a tree, a stock market index, or the time it takes to travel from A to B in a crowded city is often well captured by stationary stochastic processes for extended periods of time. Then sudden changes may happen, such that the distribution of leaf positions becomes different due to a storm, the stock market index is affected by the enforcement of a new law, or a blocked road causes additional traffic jams. The violation of an agent's expectation caused by such sudden changes is perceived by the agent as surprise, which can be seen as a measure of how much the agent's current belief differs from reality.

Surprise, with its physiological manifestations in pupil dilation (Nassar et al., 2012; Preuschoff, t Hart, & Einhauser, 2011) and EEG signals (Mars et al., 2008; Modirshanechi, Kiani, & Aghajan, 2019; Ostwald et al., 2012), is believed to modulate learning, potentially through the release of specific neurotransmitters (Gerstner, Lehmann, Liakoni, Corneil, & Brea, 2018; Yu & Dayan, 2005) so as to allow animals and humans to adapt quickly to sudden changes. The quick adaptation to novel situations has been demonstrated in a variety of learning experiments (Behrens, Woolrich, Walton, & Rushworth, 2007; Glaze, Kable, & Gold, 2015; Heilbron & Meyniel, 2019; Nassar et al., 2012; Nassar, Wilson, Heasly, & Gold, 2010; Yu & Dayan, 2005). The bulk of computational work on surprise-based learning can be separated into two groups. Studies in the field of computational neuroscience have focused on biological plausibility with little emphasis on the accuracy of learning (Behrens et al., 2007; Bogacz, 2017; Faraji, Preuschoff, & Gerstner, 2018; Friston, 2010; Friston, FitzGerald, Rigoli, Schwartenbeck, & Pezzulo, 2017; Nassar et al., 2012; Nassar, Wilson, Heasly, & Gold, 2010; Ryali, Reddy, & Yu, 2018; Schwartenbeck, FitzGerald, Dolan, & Friston, 2013; Yu & Dayan, 2005), whereas exact and approximate Bayesian online methods (Adams & MacKay, 2007; Fearnhead & Liu, 2007) for change point detection and parameter estimation have been developed without any focus on biological plausibility (Aminikhanghahi & Cook, 2017; Cummings, Krehbiel, Mei, Tuo, & Zhang, 2018; Lin, Sharpnack, Rinaldo, & Tibshirani, 2017; Masegosa et al., 2017; Wilson, Nassar, & Gold, 2010).

In this work, we take a top-down approach to surprise-based learning. We start with a generative model of change points similar to the one that has been the starting point of multiple experiments (Behrens et al., 2007; Findling, Chopin, & Koechlin, 2019; Glaze et al., 2015; Heilbron & Meyniel, 2019; Nassar et al., 2012, 2010; Yu & Dayan, 2005). We demonstrate that Bayesian inference on such a generative model can be interpreted as modulation of learning by surprise; we show that this modulation leads to a natural definition of surprise that is different from but closely related to the Shannon Surprise (Shannon, 1948). Moreover, we derive three novel approximate online algorithms with update rules that inherit the surprise-modulated adaptation rate of exact Bayesian inference. The overall goal of this article is to give a Bayesian interpretation for surprise-based learning in the brain and to find approximate methods that are computationally efficient and biologically plausible while maintaining a high level of learning accuracy. As a by-product, our approach provides theoretical insights on commonalities and differences among existing surprise-based and approximate Bayesian approaches. Importantly, our approach makes specific experimental predictions.

In section 2, we introduce the generative model and we present our surprise-based interpretation of Bayesian inference and our three approximate algorithms. Next, we use simulations to compare our algorithms with existing ones on two different tasks inspired by and closely related to real experiments (Behrens et al., 2007; Mars et al., 2008; Nassar et al., 2012, 2010; Ostwald et al., 2012). At the end of section 2, we formalize two experimentally testable predictions of our theory and illustrate them with simulations. A brief review of related studies as well as a few directions for further work are supplied in section 3.

Given a sequence of observations $y1:t$⁠, the agent's belief $π(t)(θ)$ about the parameter $Θt$ at time $t$ is defined as the posterior probability distribution $P(Θt=θ|y1:t)$⁠. In the online learning setting studied here, the agent's goal is to update the belief $π(t)(θ)$ to the new belief $π(t+1)(θ)$⁠, or an approximation thereof, upon observing $yt+1$⁠.

A simplified real-world example of such an environment is illustrated in Figure 1B. Imagine that every day, a friend meets you at the coffee shop, starting after work from her office (see Figure 1B, left). To do so, she needs to cross a river via a bridge. The time of arrival of your friend (⁠$yt$⁠) exhibits some variability due to various sources of stochasticity (e.g., traffic and your friend's daily workload), but it has a stable average over time (⁠$θt$⁠). However, if a new bridge is opened, your friend arrives earlier, since she no longer has to take a detour (see Figure 1B, right). The moment of opening the new bridge is indicated by $ct+1=1$ in our framework and the sudden change in the average arrival time of your friend by a sudden change from $θt$ to $θt+1$⁠. Even without any explicit discussion with your friend about this situation and only by observing her actual arrival time, you can notice the abrupt change and hence adapt your schedule to the new situation.

Note that if $π$ is the exact Bayesian belief defined as above in equation 2.5, then $P(y;π(t'))=P(Yt'+1=y|y1:t',ct'+1=0)$⁠. In section 2.2 we will also use $P(y;π^(t'))$ for an arbitrary $π^(t')$⁠. Two particularly interesting cases of equation 2.6 are $P(yt+1;π(t))$⁠, the probability of a new observation $yt+1$ with the current belief $π(t)$⁠, and $P(yt+1;π(0))$⁠, the probability of a new observation $yt+1$ with the prior belief $π(0)$⁠.

This definition of surprise measures how much more probable the current observation is under the naive prior $π(0)$ relative to the current belief $π(t)$ (see section 3 for further interpretation). This probability ratio is the Bayes factor (Efron & Hastie, 2016; Kass & Raftery, 1995) that tests the prior belief $π(0)$ against the current belief $π(t)$⁠. We emphasize that our definition of surprise is not arbitrary, but essential in order to write the exact inference in equation 2.5 on the generative model in the compact recursive form indicated in the proposition that follows. Moreover, as we show later, this term can be identified in multiple learning algorithms (among them, Nassar et al., 2010, 2012), but it has never been interpreted as a surprise measure. In the following sections we establish the generality of this computational mechanism and identify it as a common feature of many learning algorithms.

$πB(t+1)(θ)$ describes the incorporation of the new information into the current belief via Bayesian updating.

Using the above definitions and equation 2.5, we have for the generative model of Figure 1A and equations 2.1 to 2.4 the following proposition:

As a conclusion, our first main result is that a split as in equation 2.12 with a weighting factor (“adaptation rate” $γ$⁠) as in equation 2.9 is exact and always possible for the class of environments defined by our hierarchical generative model. This surprise modulation gives rise to specific testable experimental predictions discussed later.

Despite the simplicity of the recursive formula in equation 2.10, the updated belief $π(t+1)$ is generally not in the same family of distributions as the previous belief $π(t)$⁠, for example, the result of averaging two normal distributions is not a normal distribution. Hence it is in general impossible to find a simple and exact update rule for, say, some sufficient statistic. As a consequence, the memory demands for $π(t+1)$ scale linearly in time, and updating $π(t+1)$ using $π(t)$ needs $O(t)$ operations. In the following sections, we investigate three approximations, algorithms 1 to 3, that have simple update rules and finite memory demands, so that the updated belief remains tractable over a long sequence of observations.

As our second main result, we show that all our three approximate algorithms inherit the surprise-modulated adaptation rate from the exact Bayesian approach, equations 2.9 and 2.12. The first algorithm adapts an earlier algorithm of surprise minimization learning (SMiLe, Faraji et al., 2018) to variational learning. We refer to our novel algorithm as Variational SMiLe and abbreviate it as VarSMiLe (see algorithm 1). The second algorithm is based on message passing (Adams & MacKay, 2007) restricted to a finite number of messages $N$⁠. We refer to this algorithm as MP$N$ (see algorithm 2). The third algorithm uses the ideas of particle filtering (Gordon, Salmond, & Smith, 1993) for an efficient approximation for our hierarchical generative model. We refer to our approximate algorithm as Particle Filtering with $N$ particles and abbreviate it by pf$N$ (see algorithm 3). All algorithms are computationally efficient, have finite memory demands, and are biologically plausible; Particle Filtering has possible neuronal implementations (Huang & Rao, 2014; Kutschireiter, Surace, Sprekeler, & Pfister, 2017; Legenstein & Maass, 2014; Shi & Griffiths, 2009), MP$N$ can be seen as a greedy version of pf$N$ without sampling, and Variational SMiLe may be implemented by neo-Hebbian (Gerstner et al., 2018; Lisman, Grace, & Duzel, 2011) update rules. Simulation results show that the performance of the three approximate algorithms is comparable to and more robust across environments than other state-of-the-art approximations.

where $γt+1=γSBF(yt+1;π^(t)),m$ is given by equation 2.9 with a free parameter $m>0$⁠, which can be tuned to each environment. By doing so, we still have the explicit trade-off between two terms as in equation 2.10, but in the logarithms, yet an advantageous consequence of averaging over logarithms is that if the likelihood function $PY$ is in the exponential family and the initial belief $π(0)$ is its conjugate prior, then $π^(t+1)$ and $π(0)$ are members of the same family. In this particular case, we arrive at a simple update rule for the parameters of $π^(t+1)$ (see algorithm 1 for pseudocode and section 4 for details). As it is common in variational approaches (Beal, 2003), the price of this simplicity is that except for the trivial cases of $pc=0$ and $pc=1$⁠, there is no evidence other than simulations that the update rule of equation 2.13 will end up at an approximate belief close to the exact Bayesian belief.

where the bound $Bt+1∈0,DKL[π^B(t+1)(θ)||P(θ|yt+1)]$ is a decreasing function of the Bayes Factor Surprise $SBF(yt+1;π^(t))$ (see section 4 for proof), and $P(θ|yt+1)$ is given by equation 2.11.

Because of the similarity of the constraint optimization problem in equations 2.14 and 2.15 to the Surprise Minimization Learning rule, (SMiLe) (Faraji et al., 2018), we call this algorithm the Variational Surprise Minimization Learning rule, or for short, the Variational SMiLe rule. The differences between SMiLe and variational SMiLe are discussed in section 4.

Our variational method, and particularly its surprise-modulated adaptation rate, is complementary to earlier studies (Masegosa et al., 2017; Özkan, Šmídl, Saha, Lundquist, & Gustafsson, 2013) in machine learning that assumed different generative models and used additional assumptions and different approaches for deriving the learning rule.

For a hierarchical generative model similar to ours, a message-passing algorithm has been used to perform exact Bayesian inference (Adams & MacKay, 2007), where the algorithm's memory demands and computational complexity scale linearly in time $t$⁠. In this section, we first explain the idea of the message-passing algorithm of Adams and MacKay (2007) and its relation to our proposition. We then present our approximate version of this algorithm, which has a constant (in time) computational complexity and memory demands.

The computational complexity and memory requirements of the complete message-passing algorithm increase linearly with time $t$⁠. To deal with this issue and to have a constant computation and memory demands over time, we implemented a message-passing algorithm of the form of equations 2.17 to 2.19, but with a fixed number $N$ of particles, chosen as those with the highest weights $wt(k)$⁠. Therefore, our second algorithm adds a new approximation step to the full message-passing algorithm of Adams and MacKay (2007): Whenever $t>N$⁠, after adding the new particle with the weight as in equation 2.18 and updating the previous weights as in equation 2.19, we discard the particle with the smallest weight (i.e., set its weight equal to 0) and renormalize the weights. By doing so, we always keep the number of particles with nonzero weights equal to $N$⁠. Note that for $t≤N$⁠, our algorithm is exact and identical to the message-passing algorithm of Adams and MacKay (2007). We call our modification of the message-passing algorithm of Adams and MacKay (2007) “Message Passing $N$ (MP$N$⁠).”

To deal with the computational complexity and memory requirements, one may alternatively keep only the particles with weights greater than a cut-off threshold (Adams & MacKay, 2007). However, such a constant cut-off leads to a varying number (smaller or equal to $t$⁠) of particles in time. Our approximation MP$N$ can therefore be seen as a variation of the thresholding algorithm in Adams and MacKay (2007) with fixed number of particles $N$⁠, and hence a variable cut-off threshold. The work of Fearnhead and Liu (2007) follows the same principle as Adams and MacKay (2007), but employs stratified resampling to eliminate particles with negligible weights in order to reduce the total number of particles. Their resampling algorithm involves solving a complicated nonlinear equation at each time step, which makes it unsuitable for a biological implementation. In addition, we experienced that in some cases, the small errors introduced in the resampling step of the algorithm of Fearnhead and Liu (2007) accumulated and led to worse performance than our MP$N$ algorithm, which simply keeps the $N$ particles with the highest weight at each time step.

For the case where the likelihood function $PY(y|θ)$ is in the exponential family and $π(0)$ is its conjugate prior, the resulting algorithm of MP$N$ has a simple update rule for the belief parameters (see algorithm 2 and section 4 for details). For comparison, we also implemented in our simulations the full message-passing algorithm of Adams and MacKay (2007) with an almost zero cut-off (machine precision), which we consider as our benchmark “Exact Bayes,” as well as the stratified optimal resampling algorithm of Fearnhead and Liu, which we call ``SORN.''

Interestingly, the above formulas entail the same surprise modulation and the same trade-off as proposed by the proposition's equation 2.10. For the weight update, there is a trade-off between an exact Bayesian update and keeping the value of the previous time step, controlled by a adaptation rate modulated exactly in the same way as in equation 2.10. Note that in contrast to equation 2.10, the trade-off for the particles' weights is not between forgetting and integrating, but between maintaining the previous knowledge and integrating. However, the change probability (see equation 2.23) for sampling is equal to the adaptation rate and is an increasing function of surprise. As a result, although the weights are updated less for surprising events, a higher surprise causes a higher probability for change, indicated by $ct+1(i)=1$⁠, which implies forgetting, because for a particle $i$ with $ct+1(i)=1$⁠, the associated belief $π^i(t+1)=P(Θt+1=θ|ct+1(i)=1,c1:t(i),y1:t+1)$ is equal to $P(Θt+1=θ|ct+1(i)=1,yt+1)=P(θ|yt+1)$ (see Figure 1A), which is equivalent to a reset of the belief as in equation 2.12. In other words, while in MP$N$ and the exact Bayesian inference in the proposition's equation 2.10, the trade-off between integration and reset is accomplished by adding at each time step a new particle with weight $γt+1$⁠, in Particle Filtering, it is accomplished via sampling. As a conclusion, the above formulas are essentially the same as the update rules of MP$N$ (see equations 2.19 and 2.18) and have the same spirit as the recursive update of the proposition's equation 2.10.

Equations 2.21 and 2.23 can be applied to the case where the likelihood function $PY(y|θ)$ is in the exponential family and $π(0)$ is its conjugate prior. The resulting algorithm (algorithm 3) has a particularly simple update rule for the belief parameters (see section 4 for details).

The theory of particle filter methods is well established (Doucet et al., 2000; Gordon et al., 1993; Särkkä, 2013). Particle filters in simpler (Brown & Steyvers, 2009) or more complex (Findling et al., 2019) forms have also been employed to explain human behavior. Here we derived a simple particle filter for the general case of generative models of equations 2.2 to 2.4. Our main contribution is to show that the use of the optimal proposal distribution in this particle filter leads to a surprise-based update scheme.

Other existing algorithms (Adams & MacKay, 2007; Faraji et al., 2018; Fearnhead & Liu, 2007; Nassar et al., 2012, 2010) can also be formulated in the surprise-modulation framework of equations 2.9 and 2.12 (see section 4). Moreover, in order to allow for a transparent discussion and for fair comparisons in simulations, we extended the algorithms of Nassar et al. (2012, 2010) to a more general setting. Here we give a brief summary of the algorithms we considered. A detailed analysis is provided in section 4.5.

The algorithms of Nassar et al. (2012, 2010) were originally designed for a gaussian estimation task (see section 2.3 for details of the task) with a broad uniform prior. We extended them to the more general case of gaussian tasks with gaussian priors, and we call our extended versions Nas10$*$ and Nas12$*$ for Nassar et al. (2010) and Nassar et al. (2012), respectively (for a performance comparison between our extended algorithms and their original versions see Figures 12 and 13). Both algorithms have the same surprise-modulation as in our proposition (see equation 2.10). There are multiple interpretations of the approaches of Nas10$*$ and Nas12$*$ and links to other algorithms. One such link we identify is in relation to Particle Filtering with a single particle (pf1). More specifically, one can show that pf1 behaves in expectation similar to Nas10$*$ and Nas12$*$ (see section 4 and the appendix).

With the goal of gaining a better understanding of different approximate algorithms, we evaluated the departure of their performance from the exact Bayesian algorithm in terms of mean squared error (MSE) of estimating $Θt$ (see section 4), on two tasks inspired by and closely related to real experiments (Behrens et al., 2007; Maheu, Dehaene, & Meyniel, 2019; Mars et al., 2008; Nassar et al., 2012, 2010; Ostwald et al., 2012): a gaussian and a categorical estimation task.

We compared our three novel algorithms VarSMiLe, Particle Filtering (pf$N$⁠, where $N$ is the number of particles), and Message Passing with finite number of particles $N$ (MP$N$⁠) to the online exact Bayesian Message Passing algorithm (Adams & MacKay, 2007) (Exact Bayes), which yields the optimal solution with $Θ^t=Θ^tOpt$⁠. Furthermore, we included in the comparison the stratified optimal resampling algorithm (Fearnhead & Liu, 2007) (SOR$N$⁠, where $N$ is the number of particles), our variant of Nassar et al. (2010) (Nas10$*$⁠) and Nassar et al. (2012) (Nas12$*$⁠), the Surprise-Minimization Learning algorithm of Faraji et al. (2018) (SMiLe), as well as a simple Leaky Integrator (Leaky; see section 4).

The algorithms Exact Bayes and SORN come from the field of change point detection, and whereas the former has high memory demands, the latter has the same memory demands as our algorithms pfN and MPN. The algorithms Nas10$*$⁠, Nas12$*$⁠, and SMiLe, on the other hand, come from the human learning literature and are more biologically oriented.

The task is a generalized version of the experiment of Nassar et al. (2012, 2010). The goal of the agent is to estimate the mean $θt=μt$ of observed samples, which are drawn from a gaussian distribution with known variance $σ2$⁠: $yt+1|μt+1∼N(μt+1,σ2)$⁠. The mean $μt+1$ is itself drawn from a gaussian distribution $μt+1∼N(0,1)$ whenever the environment changes. In other words, the task is a special case of the generative model of equations 2.2 to 2.4, with $π(0)(μt)=N(μt;0,1)$ and $PY(yt|μt)=N(yt;μt,σ2)$⁠. An example of the task is shown in Figure 2A.

We simulated the task for all combinations of $σ∈{0.1,0.5,1,2,5}$ and $pc∈{0.1,0.05,0.01,0.005,0.001,0.0001}$⁠. For each combination of $σ$ and $pc$⁠, we first tuned the free parameter of each algorithm—$m$ for SMiLe and Variational SMiLe, the leak parameter for the Leaky Integrator, and the $pc$ of Nas10$*$ and Nas12$*$—by minimizing the MSE on three random initializations of the task. For the Particle Filter (pf$N$⁠), the Exact Bayes, the MP$N$⁠, and the SOR$N$⁠, we empirically checked that the true $pc$ of the environment was indeed the value that gave the best performance, and we used this value for the simulations. We evaluated the performance of the algorithms on 10 different random task instances of $105$ steps each for $pc∈{0.1,0.05,0.01,0.005}$ and $106$ steps each for $pc∈{0.001,0.0001}$ (in order to sample more change points). Note that the parameter $σ$ is not estimated, and its actual value is used by all algorithms except the Leaky Integrator.

In Figure 2B we show the $MSE[Θ^t|Rt=n]$ in estimating the parameter after $n$ steps since the last change point, for each algorithm, computed over multiple changes, for two exemplar task settings. The Particle Filter with 20 particles (pf20), the VarSMiLe, and the Nas12$*$ have an overall performance very close to that of the Exact Bayes algorithm (⁠$MSE[Θ^tOpt|Rt=n]$⁠), with much lower memory requirements. VarSMiLe sometimes slightly outperforms the other two early after an environmental change (see Figure 2B, right), but shows slightly higher error values at later phases. The MP$N$ algorithm is the closest one to the optimal solution (i.e., Exact Bayes) for low $σ$ (see Figure 2B, left), but its performance is much worse for the case of high $σ$ and low $pc$ (see Figure 2B, right). For SOR$N$ we observe a counter intuitive behavior in the regime of low $σ$⁠: the inclusion of more particles leads to worse performance (see Figure 2B, left). At higher $σ$ levels, the performance of SOR20 is close to optimal and better than the MP20 in later time steps. This may be due to the fact that the MP$N$ discards particles in a deterministic and greedy way (i.e., the one with the lowest weight), whereas for the SOR$N$⁠, there is a component of randomness in the process of particle elimination, which may be important for environments with higher stochasticity.

For the Leaky Integrator, we observe a trade-off between good performance in the transient phase and the stationary phase; a fixed leak value cannot fulfill both requirements. The SMiLe rule, by construction, never narrows its belief $π^(θ)$ below some minimal value, which allows it to have a low error immediately after a change but leads later to high errors. Its performance deteriorates for higher $σ$ (see Figure 2B, right). The Nas10$*$ performs well for low but not for higher values of $σ$⁠. Despite the fact that a Particle Filter with one particle (pf1) is in expectation similar to Nas10$*$ and Nas12$*$ (see section 4), it performs worse than these two algorithms on trial-by-trial measures. Still, it performs better than the MP1 and identically to the SOR1.

In Figure 3A, we have plotted the average of $MSE[Θ^tOpt]$ of the Exact Bayes algorithm over the whole simulation time for each of the considered $σ$ and $pc$ levels. The difference between the other algorithms and this benchmark is called $ΔMSE[Θ^t]$ (see section 4) and is plotted in Figures 3C to 3F. All algorithms except for the SOR20 have lower average error values for low $σ$ and low $pc$ than high $σ$ and high $pc$⁠. The Particle Filter pf20 and the Message Passing MP20 have the smallest difference from the optimal solution. The average error of MP20 is higher than that of pf20 for high $σ$ and low $pc$⁠, whereas pf20 is more robust across levels of environmental parameters. The worst-case performance for pf20 is $ΔMSE[Θ^t]=0.033$ for $σ=5$ and $pc=0.0001$⁠, and for SOR20, it is $ΔMSE[Θ^t]=0.061$ for $σ=0.1$ and $pc=0.1$⁠. The difference between these two worst-case scenarios is significant (⁠$p-value=2.79×10-6$⁠, two-sample $t$-test, 10 random seeds for each algorithm). Next in performance are the algorithms, Nas12$*$ and VarSMiLe. VarSMiLe exhibits its largest deviation from the optimal solution for high $σ$ and low $pc$⁠, but is still more resilient compared to the MP$N$ algorithms for this type of environment. Among the algorithms with only one unit of memory demands—pf1, MP1, SOR1, VarSMiLe, SMiLe, Leaky, Nas10$*$ and Nas12$*$—the winners are VarSMiLe and Nas12$*$⁠. The SOR20 has low error overall but unexpectedly high error for environmental settings that are presumably more relevant for biological agents (intervals of low stochasticity marked by abrupt changes). The simple Leaky Integrator performs well at low $σ$ and $pc$ but deviates more from the optimal solution as these parameters increase (see Figure 3F). The SMiLe rule performs best at lower $σ$⁠, that is, more deterministic environments.

A summary graph, where we collect the $ΔMSE[Θ^t]$ across all levels of $σ$ and $pc$⁠, is shown in Figure 4. We can see that pf20, Nas12$*$⁠, and VarSMiLe give the lowest worst case (lowest maximum value) $ΔMSE[Θ^t]$ and are statistically better than the other eight algorithms (the error bars indicate the standard error of the mean across the ten random task instances).

The task is inspired by the experiments of Behrens et al. (2007), Maheu et al. (2019), Mars et al. (2008), and Ostwald et al. (2012). The goal of the agent is to estimate the occurrence probability of five possible states. Each observation $yt+1∈{1,…,5}$ is drawn from a categorical distribution with parameters $θt+1=pt+1$⁠, that is, $yt+1|pt+1∼Cat(yt+1;pt+1)$⁠. When there is a change $Ct+1=1$ in the environment, the parameters $pt+1$ are drawn from a Dirichlet distribution $Dir(s·1)$⁠, where $s∈(0,∞)$ is the stochasticity parameter. In relation to the generative model of equations 2.2 to 2.4, we thus have $π(0)(pt)=Dir(pt;s·1)$ and $PY(yt|pt)=Cat(yt;pt)$⁠. An illustration of this task is depicted in Figure 5A.

We considered the combinations of stochasticity levels $s∈{0.01$⁠, 0.1, 0.14, 0.25, 1, 2, 5$}$ and change probability levels $pc∈{0.1$⁠, 0.05, 0.01, 0.005, 0.001, 0.0001$}$⁠. The algorithms of Nassar et al. (2012, 2010) were specifically developed for a gaussian estimation task and cannot be applied here. All other algorithms were first optimized for each combination of environmental parameters before an experiment starts, and then evaluated on 10 different random task instances, for $105$ steps each for $pc∈{0.1$⁠, 0.05, 0.01, 0.005$}$ and $106$ steps each for $pc∈{0.001$⁠, 0.0001$}$⁠. The parameter $s$ is not estimated, and its actual value is used by all algorithms except the Leaky Integrator.

In summary, our simulation results of the two tasks collectively suggest that our Particle Filtering (pf$N$⁠) and Message Passing (MP$N$⁠) algorithms achieve a high level of performance, very close to the one of biologically less plausible algorithms with higher (Exact Bayes) and same (SOR$N$⁠) memory demands. Moreover, their behavior is more consistent across tasks. Finally, among the algorithms with memory demands of one unit, VarSMiLe performs best.

In all algorithms we considered, the environment's hyperparameters are assumed to be known. We can distinguish between two types of hyperparameters in our generative model: the parameters of the likelihood function (e.g., $σ$ in the gaussian task), and the $pc$ and the parameters of the conjugate prior (e.g., $s$ in the categorical task). Hyperparameters of the first type can be added to the parameter vector $θ$ and be inferred with the same algorithm. However, learning the second type of hyperparameters is not straightforward. By assuming that these hyperparameters are learned more slowly than $θ$⁠, one can fine-tune them after each $n$ (e.g., 10) change points, while change points can be detected by looking at the particles for the Particle Filter and at the peaks of surprise values for VarSMiLe. Other approaches to hyperparameter estimation can be found in Doucet and Tadić (2003), George and Doss (2017), Liu and West (2001), and Wilson et al. (2010).

When the hyperparameters are fixed, a mismatch between the assumed values and the true values is a possible source of errors. In this section, we investigate the robustness of the algorithms to a mismatch between the assumed and the actual probability of change points. To do so, we first tuned each algorithm's parameter for an environment with a change probability $pc$ and then tested the algorithms in environments with different change probabilities while keeping the parameter fixed. For each new environment with a different change probability, we calculated the difference between the MSE of these fixed parameters and the optimal MSE, that is, the resulting MSE for the case that the Exact Bayes' parameter is tuned for the actual $pc$⁠.

More precisely, if we denote as $MSE[Θ^t;pc',pc]$ the MSE of an algorithm with parameters tuned for an environment with $pc'$⁠, applied in an environment with $pc$⁠, we calculated the mean regret, defined as $MSE[Θ^t;pc',pc]-MSE[Θ^tOpt,pc]$⁠, over time; note that the second term is equal to $MSE[Θ^t;pc,pc]$ when the algorithm Exact Bayes is used for estimation. The lower the values and the flatter the curve of the mean regret, the better the performance and the robustness of the algorithm in the face of lacking knowledge of the environment. The slope of the curve indicates the degree of deviations of the performance as we move away from the optimally tuned parameter. We ran three random (and same for all algorithms) tasks initializations for each $pc$ level.

In Figure 8 we plot the mean regret for each algorithm for the gaussian task for four pairs of $s$ and $pc'$ levels. For $σ=0.1$ and $pc'=0.04$ (see Figure 8A), the Exact Bayes and the MP20 show the highest robustness (smallest regret) and are closely followed by the pf20, VarSMiLe, and Nas12$*$ (note the regret's small range of values). The lower the actual $pc$⁠, the higher the regret, but the changes are still very small. The curves for the SMiLe and the Leaky Integrator are also relatively flat, but the mean regret is much higher. The SOR20 is the least robust algorithm.

Similar observations can be made for $σ=0.1$ and $pc'=0.004$ (see Figure 8B). In this case, the performance of all algorithms deteriorates strongly when the actual $pc$ is higher than the assumed one.

However, for $σ=5$ (see Figures 8C and 8D), the ranking of algorithms changes. The SOR20 is very robust for this level of stochasticity. The pf20 and MP20 perform similarly for $pc=0.04$⁠, but for lower $pc'$⁠, the pf20 is more robust and the MP20 exhibits high fluctuations in its performance. The Nas12$*$ is quite robust at this $σ$ level. Overall for Exact Bayes, SOR20, pf20, VarSMiLe, and Nas12$*$⁠, a mismatch of the assumed $pc$ from the actual one does not deteriorate the performance dramatically for $σ=5$⁠, $pc'=0.004$ (see Figure 8D). The SMiLe and the Leaky Integrator outperform the other algorithms for higher $pc'$ if $pc<pc'$ (see Figure 8C). A potential reason is that the optimal behavior for the Leaky Integrator (according to the tuned parameters) is to constantly integrate new observations into its belief (i.e., to act like a Perfect Integrator) regardless of the $pc'$ level. This feature makes it blind to the $pc$ and therefore very robust against the lack of knowledge of it (see Figure 8C).