Error-driven learning rules have received considerable attention because of their close relationships to both optimal theory and neurobiological mechanisms. However, basic forms of these rules are effective under only a restricted set of conditions in which the environment is stable. Recent studies have defined optimal solutions to learning problems in more general, potentially unstable, environments, but the relevance of these complex mathematical solutions to how the brain solves these problems remains unclear. Here, we show that one such Bayesian solution can be approximated by a computationally straightforward mixture of simple error-driven ‘Delta’ rules. This simpler model can make effective inferences in a dynamic environment and matches human performance on a predictive-inference task using a mixture of a small number of Delta rules. This model represents an important conceptual advance in our understanding of how the brain can use relatively simple computations to make nearly optimal inferences in a dynamic world.
Decisions are often guided by beliefs about the probability and utility of potential outcomes. These beliefs are learned through past experiences that, in stable environments, can be used to generate accurate predictions. However, in dynamic environments, changes can occur that render past experiences irrelevant for predicting future outcomes. For example, after a change in government, historical tax rates may no longer be a reliable predictor of future tax rates. Thus, an important challenge faced by a decision-maker is to identify and respond to environmental change-points, corresponding to when previous beliefs should be abandoned and new beliefs should be formed. A toy example of such a situation is shown in figure 1A, where we plot the price of a fictional stock over time. In this example, the stock price on a given day (red dots) is generated by sampling from a Gaussian distribution with variance $1 and a mean (dashed black line) that starts at $10 before changing abruptly to $20 at a change-point, perhaps caused by the favorable resolution of a court case. A trader only sees the stock price and not the underlying mean but has to make predictions about the stock price on the next day. One common strategy for computing this prediction is based on the Delta rule: (1) According to this rule, an observation, , is used to update an existing prediction, , based on the learning rate, and the prediction error, . Despite its simplicity, this learning rule can provide effective solutions to a wide range of machine-learning problems [1], [2]. In certain forms, it can also account for numerous behavioral findings that are thought to depend on prediction-error signals represented in brainstem dopaminergic neurons, their inputs from the lateral habenula, and their targets in the basal ganglia and the anterior cingulate cortex [3]–[15]. Unfortunately, this rule does not perform particularly well in the presence of change-points. We illustrate this problem with a toy example in figure 1B and C. In panel B, we plot the predictions of this model for the toy data set when is set to 0. 2. In this case, the algorithm does an excellent job of computing the mean stock value before the change-point. However, it takes a long time to adjust its predictions after the change-point, undervaluing the stock for several days. In figure 1C, we plot the predictions of the model when. In this case, the model responds rapidly to the change-point but has larger errors during periods of stability. One way around this problem is to dynamically update the learning rate on a trial-by-trial basis between zero, indicating that no weight is given to the last observed outcome, and one, indicating that the prediction is equal to the last outcome [16], [17]. During periods of stability, a decreasing learning rate can match the current belief to the average outcome. After change-points, a high learning rate shifts beliefs away from historical data and towards more recent, and more relevant, outcomes. These adaptive dynamics are captured by Bayesian ideal-observer models that determine the rate of learning based on the statistics of change-points and the observed data [18]–[20]. An example of the behavior of the Bayesian model is shown in figure 1D. In this case, the model uses a low learning rate in periods of stability to make predictions that are very close to the mean, then changes to a high learning rate after a change-point to adapt more quickly to the new circumstances. Recent experimental work has shown that human subjects adaptively adjust learning rates in dynamic environments in a manner that is qualitatively consistent with these algorithms [16], [17], [21]. However, it is unlikely that subjects are basing these adjustments on a direct neural implementation of the Bayesian algorithms, which are complex and computationally demanding. Thus, in this paper we ask two questions: 1) Is there a simpler, general algorithm capable of adaptively adjusting its learning rate in the presence of change-points? And 2) Does the new model better explain human behavioral data than either the full Bayesian model or a simple Delta rule? We address these questions by developing a simple approximation to the full Bayesian model. In contrast to earlier work that used a single Delta rule with an adaptive learning rate [17], [21], our model uses a mixture of biologically plausible Delta rules, each with its own, fixed learning rate, to adapt its behavior in the presence of change-points. We show that the model provides a better match to human performance than the other models. We conclude with a discussion of the biological plausibility of our model, which we propose as a general model of human learning.
Human subject protocols were approved by the University of Pennsylvania internal review board. Informed consent was given by all participants prior to taking part in the study. To familiarize readers with change-point processes and the Bayesian model, we first review these topics in some detail and then turn our attention to the reduced model. In this paper we are concerned with data generated from change-point processes. An example of such a process generating Gaussian data is given in figure 2. We start by defining a hazard rate, , that in the general case can be variable over time but for our purposes is assumed to be constant. Change-point locations are then generated by sampling from a Bernoulli distribution with this hazard rate, such that the probability of a change-point occurring at time is (figure 2A). In between change-points, in periods we term ‘epochs, ’ the generative parameters of the data are constant. Within each epoch, the values of the generative parameters, , are sampled from a prior distribution, for some hyper-parameters and that will be described in more detail in the following sections. For the Gaussian example, is simply the mean of the Gaussian at each time point. We generate this mean for each epoch (figure 2B) by sampling from the prior distribution shown in figure 2C. Finally, we sample the data points at each time, from the generative distribution (figure 2D and E). The goal of the full Bayesian model [18], [19] is to make accurate predictions in the presence of change-points. This model infers the predictive distribution, , over the next data point, , given the data observed up to time, . In the case where the change-point locations are known, computing the predictive distribution is straightforward. In particular, because the parameters of the generative distribution are resampled independently at a change-point (more technically, the change-points separate the data into product partitions [22]) only data seen since the last change-point are relevant for predicting the future. Therefore, if we define the run-length at time, , as the number of time steps since the last change-point, we can write (2) where we have introduced the shorthand to denote the predictive distribution given the last time points. Assuming that our generative distribution is parameterized by parameters, then is straightforward to write down (at least formally) as the marginal over (3) where is the inferred distribution over given the last time points, and is the likelihood of the data given the generative parameters. When the change-point locations are unknown the situation is more complex. In particular we need to compute a probability distribution over all possible values for the run-length given the observed data. This distribution is called the run-length distribution. Once we have the run-length distribution, we can compute the predictive distribution in the following way. First we compute the expected run-length on the next trial, ; i. e. , (4) where the sum is over all possible values of the run-length at time and is the change-point prior that describes the dynamics of the run-length over time. In particular, because the run-length either increases by one, with probability in between change-points, or decreases to zero, with probability at a change-point, the change-point prior, , takes the following form (5) Given the distribution, we can then compute the predictive distribution of the data on the next trial, in the following manner, (6) where the sum is over all possible values of the run-length at time. All that then remains is to compute the run-length distribution itself, which can be done recursively using Bayes' rule (7) Substituting in the form of the change-point prior for we get (8) Thus for each value of the run-length, all but two of the of the terms in equation 7 vanish and the algorithm has complexity of computations per timestep. Unfortunately, although this is a substantial improvement compared to complexity of a more naïve change-point model, this computation is still quite demanding. In principle, the total number of run-lengths we must consider is infinite, because we must allow for the possibility that a change-point occurred at any time in the past. In practice, however, it is usual to introduce a maximum run-length, , and define the change-point prior here to be (9) With this procedure, the complexity of the computation is bounded but still can remain dauntingly high. Despite the elegance of the full Bayesian algorithm, it is complex, requiring a memory of a large number () of different run-lengths, which, in the worst case, is equivalent to keeping track of all the past data. Thus, it seems an unlikely model of human cognition, and a key question is whether comparable predictive performance can be achieved with a simpler, more biologically plausible algorithm. Here we introduce an approximation to the full model that addresses these issues. First we reduce the model' s complexity by removing nodes from the update graph (Figure 3). Then we transform the update equation for into a Delta-rule update equation in which the sufficient statistic on each node updates independently of the other nodes. The resulting algorithm is a biologically plausible mixture of Delta-rules that is able to flexibly adapt its overall learning rate in the presence of change-points and whose performance is comparable with that of the full Bayesian model at a fraction of the computational cost. Below we derive new update equations for the sufficient statistics and the weights of each new node for this reduced model. To more easily distinguish the full and reduced models, we use to denote run-length in the reduced model and to denote run-length in the full model. Thus, the reduced model has nodes, where node has run-length. The set of run-lengths, , are ordered such that. Unlike the full model, the run-lengths in the reduced model can take on non-integer values, which allows greater flexibility. The first step in our approximation is to remove nodes from the update graph. This step reduces the memory demands of the algorithm but also requires us to change the update rule for the sufficient statistic and the form of the change-point prior. Consider a node with run-length. In the full Bayesian model, the sufficient statistic for this node would be (20) Note that this form of the update relies on having computed, which is the sufficient statistic at run length. In the full Bayesian model, this procedure is straightforward because all possible run-lengths are represented. In contrast, the reduced model includes only a subset of possible run-lengths, and thus a node with run-length will not exist for some values of. Therefore, the reduced model must include a new method for updating the sufficient statistic and a new form of the change-point prior. We first note that another way of writing the update for is as (21) This sliding-window update equation depends only on information available at node and thus does not rely on knowing the sufficient statistic at node. However, this update also has a high memory demand because, to update the sliding window, we have to subtract, which we can only do if we keep track of the previous data points on each node. In our model, we remove the dependence on, and hence the additional memory demands, by taking the average of equation 21. This procedure leads to a memoryless (yet approximate) form of the update equation for each node. In particular, if we take the average of equation 21 with respect to, we have (22) where we have introduced as the Delta-rule' s approximation to the mean sufficient statistic and (23) as the mean of the node. Dividing equation 21 by gives us the following form of the update for the mean (24) Note that this equation for the update of is a Delta rule, just like equation 1, with a fixed learning rate, . Thus, the reduced model simply has to keep track of for each node and update it using only the most recent data point. This form of update rule also allows us to interpret non-integer values of the run-length, , in terms of changes in the learning rate of the Delta rule on a continuum. In figure 4 we show the effect of this approximation on the extent to which past data points are used to compute the mean of each node. The sliding window rule computes the average across the last data points, ignoring all previous data. In contrast, the Delta rule computes a weighted average using an exponential that decays over time, which tends to slightly under-emphasize the contributions of recent data and over-emphasize the contributions of distant data relative to the sliding window. Reducing the number of nodes in the model also requires us to change how we update the weights of each node. In particular the update for the weights, , is given as (25) This equation is similar to equation 7 but differs in the number of run-lengths available. Crucially, this difference requires an adjustment to the change-point prior. The adjusted prior should approximate the full change-point prior (Eq. 5) as closely as possible. Recall that the full prior captures the fact that the run-length either decreases to zero if there is a change-point (with prior probability) or increases by one if there is no change-point (with prior probability). To see how to compute this adjusted prior in the reduced model, we first decompose the change-point prior into two terms corresponding to the possibility that a change-point will occur or not; i. e. , (26) where is the probability that the run-length is given that there was a change-point and that the previous run-length was. Similarly is the probability that the run-length is given that the previous run-length was and there was not a change-point. The change-point case is straightforward, because a change-point always results in a transition to the shortest run-length; i. e. , is zero, except when when it takes value 1. The no change-point case, however, is more difficult. In the full model the run-length increases by 1 when there is no change-point, thus we would like to have (27) However, because the nodes have variable spacing in the reduced model, this form is not possible as there may be no node with a run-length. We thus seek an approximation such that the prior defines an average increase in run-length of 1 if there is not a change-point. That is, we require (28) For we can match this expectation exactly by setting (29) For we approximate using (30) In this case we do not match the expected increase in run-length. For the final node, , it is impossible to transition to a longer run-length and so we simply have a self transition with probability 1; i. e. , (31) Taken together with equation 26, equations 29,30 and 31 define the change-point prior in the reduced model. Like the full Bayesian model, our reduced model also has a graphical interpretation. Again each node, , keeps track of two quantities: 1) the mean, computed according to equation 24, and 2) the weight. As in the full model, the weights are computed by passing messages along the edge of the graph. However, the structure of the graph is slightly different, with no increasing message being sent by node and an extra ‘self’ message from to itself. The increasing message has weight (32) the self message has weight (33) and the change-point message has weight (34) Finally the new weight for each node is computed by summing all of the incoming messages to implement equation 25.
First we consider the simplest cases of one and two nodes with Gaussian data. These cases have particularly simple update rules, and their output is easy to understand. We then consider the more general case of many nodes to show how the reduced model retains many of the useful properties of the full model, such as keeping track of an approximate run-length distribution and being able to handle different kinds of data. Here we derive an approximate, but analytic, expression for the average discrepancy between the predictions made by the reduced model and the ground truth generative parameters. We then use this result to compute approximately optimal node arrangements for a variety of conditions and investigate how the error varies as a function of the parameters in the model. In this section, we ask how well our model describes human behavior by fitting versions of the model to behavioral data from a predictive-inference task [24]. Briefly, in this task, 30 human subjects (19 female, 11 male) were shown a sequence of numbers between 0 and 300 that were generated by a Gaussian change-point process. This process had a mean that was randomly sampled at every change-point and a standard deviation that was constant (set to either 5 or 10) for blocks of 200 trials. Samples were constrained to be between 0 and 300 by keeping the generative means away from these bounds (the generative means were sampled from uniform distribution [from 40 to 260]) and resampling the small fraction of samples outside of this range until they lay within the range. The hazard rate was set at 0. 1 except for the first three trials following a change-point, in which case the hazard rate was zero. The subjects were required to predict the next number in the sequence and obtained more reward the closer their predictions were to the actual outcome. In particular, subjects were required to minimize the mean absolute error between prediction and outcome, which we denote. Because prediction errors depended substantially on the specific sequence of numbers generated for the given session, the exact conversion between error and monetary reward was computed by comparing performance with two benchmarks: a lower benchmark (LB) and an higher benchmark (HB). The LB was computed as the mean absolute difference between sequential generated numbers. The HB was the mean difference between mean of the generative distribution on the previous trial and the generated number. Payout was then computed as follows: (50) A benefit of this task design is that the effective learning rates used by subjects on a trial-by-trial basis can be computed in terms of their predictions following each observed outcome, using the relationships in equation 1. Our previous studies indicated that these learning rates varied systematically as a function of properties of the generative process, including its standard deviation and the occurrence of change-points [17], [24]. To better understand the computational basis for these behavioral findings, we compared five different inference models: the full Bayesian model (‘full’), the reduced model with 1 to 3 nodes and the approximately Bayesian model of Nassar et al [17]. The Nassar et al model instantiates an alternative hypothesis to the mixture of fixed Delta rules by using a single Delta rule with a single, adaptive learning rate to approximate Bayesian inference. On each trial, each of these models, , produces a prediction about the location of the next data point. To simulate the effects of decision noise, we assume that the subjects' reported predictions, , are subject to noise, such that (51) where is sampled from a Gaussian distribution with mean 0 and standard deviation that we fit as a free parameter for all models. In addition to this noise parameter, we fit the following free parameters for each model: The full model and the model of Nassar et al. have a hazard rate as their only other parameter, the one-node model has a single learning rate and the remaining models with nodes () have a hazard rate as well as the learning rates. Our fits identified the model parameters that maximized the log likelihood of the observed human predictions, , given each of the models, , which is given by (52) We used the maximum likelihood value to approximate the log Bayesian evidence, for each model using the standard Bayesian information criterion (BIC) approximation [25], which takes into account the different numbers of parameters in the different models; i. e. , (53) where is the number of free parameters in model. Models were then compared at the group level using the Bayesian method of Stephan et al. [26]. Briefly, this method aggregates the evidence from each of the models for each of the subjects to estimate two measures of model fit. The first, which we refer to as the ‘model probability’, is an estimate of how likely it is that a given model generated the data from a randomly chosen subject. The second, termed the ‘exceedance probability’, is the probability that one model is more likely than any of the others to have generated the behavior of all of the subjects. An important question when interpreting the model fits is the extent to which the different models are identifiable using these analyses. In particular we are interested in the extent to which different models can be separated on the basis of their behavior and the accuracy with which the parameters of each model can be fit. The question of model identifiability is addressed in figure 10, where we plot two confusion matrices showing the model probability (A) and the exceedance probability (B) for simulated data. These matrices were generated using simulations that matched the human-subjects experiments, with the same values of the observed stimuli, the same number of trials per experiment and the same parameter settings as found by fitting the human data. Ideally, both confusion matrices should be the identity matrix, indicating that data fit to model is always generated by model and never by any other model (e. g. , [27]). However, because of noise in the data and the limited number of trials in the experiment, it is often the case that not all of the models are completely separable. In the present case, there is good separation for the Nassar et al. , full, 1-node, and 2-node models and reasonable separation between the 3-node model and others. When we extended this analysis to include 4- and 5-node models, we found that they were indistinguishable from the 3-node model. Thus, these models are not included in our analyses, and we consider the ‘3-node model’ to represent a model with 3 or more nodes. Note that the confusion matrix showing the exceedance probability (figure 10B) is closer to diagonal than the model probability confusion matrix (figure 10A). This result reflects the fact that exceedance probability is computed at the group level (i. e. , that all the simulated data sets were generated by model M), whereas model probability computes the chance that any given simulation is best by model. To address the question of parameter estimability, we computed correlations between the simulated parameters and the parameter values recovered by the fitting procedure for each of the models. There was strong correspondence between the simulated and fit parameter values for all of the models and all correlations were significant (see supplementary table S1). The 3-node model most effectively describes the human data (Figure 11), producing slightly better fits than the model of Nassar et al. at the group level. Figure 11A shows model probability, the estimated probability that any given subject is best fit by each of the models. This measure showed a slight preference for the 3-node model over the model of Nassar et al. Figure 11B shows the exceedance probability for each of the models, the probability that each of the models best fits the data at the group level. Because this measure aggregates across the group it magnifies the differences between the models and showed a clearer preference for the 3-node model. Table 1 reports the means of the corresponding fit parameters for each of the models (see also supplementary figure S1 for plots of the full distributions of the fit parameters). Consistent with the optimal parameters derived in the previous section (figure 9E), for the 2- and 3-node models, the learning rate of the 1st node is close to one (mean ∼0. 95).
To address this question, we derived an approximation to the Bayesian model based on a mixture of Delta rules, each implemented in a separate ‘node’ of a connected graph. In this reduced model, each Delta rule has its own, fixed learning rate. The overall prediction is generated by computing a weighted sum of the predictions from each node. Because only a small number of nodes are required, the model is substantially less complex than the full Bayesian model. Qualitatively, the outputs of the reduced and full Bayesian models share many features, including the ability to quickly increase the learning rate following a change-point and reduce it during periods of stability. These features were apparent for the reduced model even with a small number of (2 or 3) nodes. Thus, effective solutions to change-point problems can be achieved with minimal computational cost. For future work, it would be interesting to consider other generative distributions, such as a Gaussian with unknown mean and variance or multidimensional data (e. g. , multidimensional Gaussians) to better assess the generality of this solution. In principle, these extensions should be straightforward to deal with in the current model, which would simply require the sufficient statistic to be a vector instead of a scalar. Another obvious extension would be to consider generative parameters that drift over time (perhaps in addition to abrupt changes at change-points) or a hazard rate that changes as a function of run-length and/or time. To address this question, we used a model-based analysis of human behavior on a prediction task with change-points. The reduced model fit the behavioral data better than either the full Bayesian model or a single learning-rate Delta rule. Our fits also suggest that a three-node model can in many cases be sufficient to explain human performance on the task. However, our experiment did not have the power to distinguish models with more that three nodes. Thus, although the results imply that the three-node model is better than the other models we tested, we cannot rule out the possibility that humans use significantly more that three learning rates. Despite this qualification, it is an intriguing idea that the brain might use just a handful of learning rates. Our theoretical analysis suggests that this scheme would yield only a small cost in performance for the variety of different problems considered here. In this regard, our model can be seen as complementary to recent work showing that in many probabilistic-inference problems faced by humans [28] and pigeons [29], as few as just one sample from the posterior can be enough to generate good solutions. It is also interesting to note that, for models with more than one node, the fastest learning rate was always close to one. Such a high learning rate corresponds to a Delta rule that does not integrate any information over time and simply uses the last outcome to form a prediction. This qualitative difference in the behavior of the fastest node could indicate a very different underlying process such as working memory for the last trial as is proposed in [30], [31]. One situation in which many nodes would be advantageous is the case in which the hazard rate changes as a function of run-length. In this case, only having a few run-lengths available would be problematic, because the changing hazard rate would be difficult to represent. Experiments designed to measure the effects of variable hazard rates on the ability to make predictions might therefore be able to distinguish whether multiple Delta rules are indeed present. The question of biological plausibility is always difficult to answer in computational neuroscience. This difficulty is especially true when the focus of the model is at the algorithmic level and is not directly tied to a specific neural architecture, like in this study. Nevertheless, one useful approach to help guide an answer to this question is to associate key components of the algorithm to known neurobiological mechanisms. Here we support the biological plausibility of our reduced model by showing that signatures of all the elements necessary to implement it have been observed in neural data. In the reduced model, the update of each node uses a simple Delta rule with a fixed learning rate. The ‘Delta’ of such an update rule corresponds to a prediction error, correlates of which have been found throughout the brain, including notably brainstem dopaminergic neurons and their targets, and have been used extensively to model behavioral data [3]–[15]. More recently, several studies have also shown evidence for representations of different learning rates, as required by the model. Human subjects performing a statistical-learning task used a pair of learning rates, one fast and one slow, that were associated with BOLD activity in two different brain areas, with the hippocampus responsible for slow learning and the striatum for fast learning [32]. A related fMRI study showed different temporal integration in one network of brain areas including the amygdala versus another, more sensory network [33]. Complementary work at the neural level found a reservoir of many different learning rates in three brain regions (anterior cingulate cortex, dorsolateral prefrontal cortex, and the lateral intraparietal area) of monkeys performing a competitive game [34]. Likewise, neural correlates of different learning rates have been identified in each of the ventral tegmental area and habenula [35]. Finally, outside of the reward system, other fMRI studies using scrambled movies have found evidence for temporal receptive fields of increasingly long time scales (equivalent to decreasingly small learning rates) up the sensory processing hierarchy [36]. Applied to our model, these results suggest that each node is implemented in a distinct, although not necessarily anatomically separated, population of neurons. For our task and the above-referenced studies, in which trials last on the order of seconds, we speculate that the mean of a node is encoded in persistent firing of neurons. Alternatively, for tasks requiring learning over longer timescales, other mechanisms such as changes in synaptic weights might play key roles in these computations. Our model also depends on the run-length distribution, . Functionally, this distribution serves as a weighting function, determining how each of the different nodes (corresponding to different run lengths) contributes to the final prediction. In this regard, the run-length distribution can be thought of as an attentional filter, similar to mechanisms of spatial or feature-based attention, evident in multiple brain regions that enhance the output of certain signals and suppress others. For longer timescales, this kind of weighting process might have analogies to certain mechanisms of perceptual decision-making that involve the readout of appropriate sensory neurons [37]. Intriguingly, these readout mechanisms are thought to be shaped by experience – governed by a Delta-rule learning process – to ultimately enhance the most reliable sensory outputs and suppress the others [38], [39]. We speculate that a similar process might help select, from a reservoir of nodes with different learning rates, those that can most effectively solve a particular task. The brain must also solve another challenge to directly implement the run-length distribution in our model. In particular, the update equation for the weights (Eq. 25) includes a constant of proportionality that serves to normalize the probability distribution. On a computer, ensuring that the run-length distribution is normalized is relatively straightforward: after the update we just divide by the sum of the node weights. In the brain, this procedure requires some kind of global divisive normalization among all areas coding different nodes. While such divisive normalization is thought to occur in the brain [40], it may be more difficult to implement over different brain regions that are far apart.
The ability to make accurate predictions is important to thrive in a dynamic world. Many predictions, like those made by a stock picker, are based, at least in part, on historical data thought also to reflect future trends. However, when unexpected changes occur, like an abrupt change in the value of a company that affects its stock price, the past can become irrelevant and we must rapidly update our beliefs. Previous research has shown that, under certain conditions, human predictions are similar to those of mathematical, ideal-observer models that make accurate predictions in the presence of change-points. Despite this progress, these models require superhuman feats of memory and computation and thus are unlikely to be implemented directly in the brain. In this work, we address this conundrum by developing an approximation to the ideal-observer model that drastically reduces the computational load with only a minimal cost in performance. We show that this model better explains human behavior than other models, including the optimal model, and suggest it as a biologically plausible model for learning and prediction.
learning
cognitive neuroscience
psychology
social and behavioral sciences
cognitive psychology
decision making
biology
neuroscience
learning and memory
A Mixture of Delta-Rules Approximation to Bayesian Inference in Change-Point Problems