Both theoretical (Bogacz and Gurney, 2007; Bogacz et al., 2010; Ratcliff and Frank, 2012; Dunovan and Verstynen, 2016; Dunovan et al., 2019; Vich et al., 2022) and experimental (Yartsev et al., 2018; Thura et al., 2022) evidence suggests that the CBGT circuits play a critical role in the evidence accumulation process (for a review, see Gupta et al., 2021). The canonical CBGT circuit (Figure 1A) includes two dissociable control pathways: the direct (facilitation) and indirect (suppression) pathways (Albin et al., 1995; Friend and Kravitz, 2014). A critical assumption of the canonical model is that the basal ganglia are organized into multiple ‘channels,’ mapped to specific action representations (Mink, 1996; Alexander et al., 1986), each containing a direct and indirect pathway. It is important to note that, for the sake of parsimony, we adopt a simple and canonical model of CBGT pathways, with action channels that are agnostic as to the location of representations (e.g., lateralization), simply assuming that actions have unique population-level representations. While a strict, segregated action channel organization may not accurately reflect the true underlying circuitry, striatal neurons have been shown to organize into task-specific spatiotemporal assemblies that qualitatively reflect independent action representations (Adler et al., 2013; Klaus et al., 2017; Barbera et al., 2016; Carrillo-Reid et al., 2011; Badreddine et al., 2022).

Figure 1 with 4 supplements see all

Download asset Open asset

Within these action channels, activation of the direct pathway, via cortical excitation of D1-expressing spiny projection neurons (SPNs) in the striatum, releases GABAergic signals that can suppress activity in the CBGT output nucleus (internal segment of the globus pallidus, GPi, in primates or substantia nigra pars reticulata, SNr, in rodents) (Kravitz et al., 2012; Gurney et al., 2001; Alexander et al., 1986; Mink, 1996). This relieves the thalamus from tonic inhibition, thereby exciting postsynaptic cortical cells and facilitating action execution. Conversely, activation of the indirect pathway via D2-expressing SPNs in the striatum controls firing in the external segment of the globus pallidus (GPe) and the subthalamic nucleus (STN), resulting in strengthened basal ganglia inhibition of the thalamus. This weakens drive to postsynaptic cortical cells and reduces the likelihood that an action is selected in cortex.

Critically, the direct and indirect pathways converge in the GPi/SNr (Kitano et al., 1998; Foster et al., 2021). This suggests that these pathways compete to control whether each specific action is selected (Dunovan et al., 2015). The apparent winner-take-all selection policy and action-channel like coding (Adler et al., 2013; Klaus et al., 2017; Barbera et al., 2016; Carrillo-Reid et al., 2011; Badreddine et al., 2022) also imply that action representations themselves compete. Altogether, this neuroanatomical evidence suggests that competition both between and within CBGT pathways control the rate of evidence accumulation during decision-making (Dunovan et al., 2019; Bariselli et al., 2019; Vich et al., 2022).

To simulate this process, we designed a spiking neural network model of the CBGT circuits, shown in Figure 1A, with dopamine-dependent plasticity occurring at the corticostriatal synapses (Vich et al., 2020; Rubin et al., 2021). Critically, although this model simulates dynamics that happen on a microscale, it can be mapped upward to infer macroscale properties, like inter-region dynamics and complex behavior, making it a useful theoretical tool for bridging across levels of analysis. The network performed a probabilistic two-arm bandit task with switching reward contingencies (see Figure 3—figure supplement 1) that followed the same general structure as our prior work (Bond et al., 2021), with the exception that block switches were deterministic for the model, happening every 10 trials, whereas in actual experiments they are generated probabilistically so as to increase the uncertainty of participant expectations of the timing of outcome switches. In brief, the network selected one of two targets, each of which returned a reward according to a specific probability distribution. The relative reward probabilities for each target were held constant at 75 and 25% and the action–outcome contingency was changed every 10 trials, on average. For the purpose of this study, we focus primarily on the neural and behavioral effects associated with a switch in the identity of the optimal target. We used four different network instances (see Supp. Methods) as a proxy for simulating individual differences over human participants.

Figure 1B shows the firing rates of dSPNs and iSPNs in the left action channel, time-locked to selection onset (when thalamic units exceed 30 Hz, t = 0), for both fast (<196 ms) and slow (>314.5 ms) decisions (see Figure 1—figure supplement 1 for node-by-node firing-rates). As expected, the dSPNs show a ramping of activity as decision onset is approached and the slope of this ramp scales with response speed. In contrast, we see that iSPN firing is sustained during slow movements and weakly ramps during fast movements. However, iSPN firing was relatively insensitive to left versus right decisions. This is consistent with our previous work showing that differences in direct pathways track primarily with choice while indirect pathway activity modulates overall response speeds (Dunovan et al., 2019; Vich et al., 2022) as supported by experimental studies (Yttri and Dudman, 2016; Forstmann et al., 2010; Maia and Frank, 2011).

We then modeled the behavior of the CBGT network using a hierarchical version of the DDM (Wiecki et al., 2013), a canonical formalism for the process of evidence accumulation during decision-making (Ratcliff, 1978; Figure 2A). This model returns four key parameters with distinct influences on evidence accumulation. The drift rate ($v$) represents the rate of evidence accumulation, the boundary height ($a$) represents the amount of evidence required to cross the decision threshold, nondecision time ($t$) is the delay in the onset of the accumulation process, and starting bias ($z$) is a bias to begin accumulating evidence for one choice over another (see ‘Materials and methods’ section).

Figure 2 with 1 supplement see all

Download asset Open asset

We tracked internal estimates of action-value and environmental change using trial-by-trial estimates of two ideal observer parameters, the belief in the value of the optimal choice ($\mathrm{\Delta }B$) and change point probability ($\mathrm{\Omega }$), respectively (see Bond et al., 2021; Nassar et al., 2010 and ‘Materials and methods’ for details). Using these estimates, we evaluated how a suspected change in the environment and the belief in optimal choice value influenced underlying decision parameters. Consistent with prior observations in humans (Bond et al., 2021), we found that both $v$ and $a$ were the most pliable parameters across experimental conditions for the network. Specifically, we found that the model mapping $\mathrm{\Delta }B$ to drift rate and $\mathrm{\Omega }$ to boundary height and the model relating $\mathrm{\Delta }B$ to drift rate provided equivocal best fits to the data over human participants ($\mathrm{\Delta }DI{C}_{\text{null}}=-29.85\pm 12.76$ and $\mathrm{\Delta }DI{C}_{\text{null}}=-22.60\pm 7.28$, respectively; see Burnham and Anderson, 1998 and ‘Materials and methods’ for guidelines on model fit interpretation). All other models failed to provide a better fit than the null model (Supplementary file 1). Consistent with prior work (Bond et al., 2021), we found that the relationship between $\mathrm{\Omega }$ and the boundary height was unreliable (mean ${\beta }_{a\sim \mathrm{\Omega }}=0.069\pm 0.152$; mean ${\displaystyle \mathrm{p}=0.232\pm 0.366}$). However, drift rate reliably increased with $\mathrm{\Delta }B$ in three of four participants (mean ${\beta }_{v\sim \mathrm{\Delta }B}=0.934\pm 0.386$; mean p<0.001 participants p<0.001; Supplementary file 2).

These effects reflect a stereotyped trajectory around a change point, whereby $v$ immediately plummets and $a$ briefly increases, with $a$ quickly recovering and $v$ slowly growing as reward feedback reinforces the new optimal target (Bond et al., 2021). Because prior work has shown that the change in $v$ is more reliable than changes in $a$ (Bond et al., 2021) and because $v$ determines the direction of choice, we focus the remainder of our analysis on the control of $v$.

To test whether these shifts in $v$ are driven by competition within and between action channels, we predicted the network’s decision on each trial using a LASSO-PCR trained on the pre-decision firing rates of the network (see ‘Measuring neural action representations’). The choice of LASSO-PCR was based on prior work building reliable classifiers from whole-brain-evoked responses that maximize inferential utility (see Meissner et al., 2011). The method is used when models are over-parameterized, as when there are more voxels than observations, relying on a combination of dimensionality reduction and sparsity constraints to find the true, effective complexity of a given model. While these are not considerations with our network model, they are with the human validation experiment that we describe next. Thus, we used the same classifier on our model as on our human participants to directly compare theoretical predictions and empirical observations. Model performance was cross-validated at the run level using a leave-one-run-out procedure, resulting in 45 folds per subject (five runs for each of the nine sessions). We then classified all trials in the hold-out set to evaluate prediction accuracy. The cross-validated accuracy for the four models, simulating individual participants, is shown in Figure 2B as ROC curves. The classifier was able to predict the chosen action with approximately 75% accuracy (72–80%) for each simulated participant, with an average area under the curve (AUC) of approximately 0.75, ranging from 0.71 to 0.77.

Examining the encoding pattern in the simulated network, we see lateralized activation over left and right action channels (Figure 1A), with opposing weights in GPi and thalamus, and, to a lesser degree, contralateral encoding in STN and in both indirect and direct SPNs in striatum. We do not observe contralateral encoding in cortex, which likely reflects the emphasis on basal ganglia structures and lumped representation of cortex in the model design.

To quantify the competition between action channels, we took the unthresholded prediction from the LASSO-PCR classifier, ${\widehat{y}}_t$, and calculated its distance from the optimal target (i.e., target with the highest reward probability) on each trial (Figure 2C). This provided an estimate of the uncertainty driven by the separability of pre-decision activity across action channels. In other words, the distance from the optimal target should increase with increased co-activation of circuits that represent opposing actions. The decision to model aggregate trial dynamics with a classifier stems from the limitations of the hemodynamic response that we will use next to vet the model predictions in humans. The low temporal resolution of the evoked BOLD signal makes finer-grained temporal analysis for the human data impossible as the signal is a low-pass filtered version of the aggregate response over the entire trial. So, we chose to represent the macroscopic network dynamics as classifier uncertainty, which cleanly links the cognitive model results to both behavior and neural dynamics at the trial-by-trial level using only two variables (drift rate and classifier uncertainty). This approach allows us to directly compare model and human results.

If the competition in action channels is also driving $v$, then there should be a negative correlation between the classifier’s uncertainty and $v$, particularly around a change point. Indeed, this is exactly what we see (Figure 2D). In fact, the uncertainty and $v$ are consistently negatively correlated across all trials in every simulated participant and in aggregate (Figure 2E). Thus, in our model of the CBGT pathways, competition between action representations drives changes in $v$ in response to environmental change.

Next, in order to rule out the possibility that these adaptive network effects emerged due to the specific parameter scheme that we used for the simulations, we re-ran our simulations using different parameter schemes. For this, we used a constrained sampling procedure (see Vich et al., 2022) to sample a range of different networks with varying degrees of speed and accuracy. This parameter search was constrained to permit regimes that result in biologically realistic firing rates (Figure 1—figure supplement 1). The simulations above arose from a parameter scheme lying in the middle of this response speed distribution (intermediate). We then chose two parameter regimes, one that produces response speeds in the upper quartile of the distribution (slow) and one that produces response speeds in the lower quartile (fast; Figure 1—figure supplement 2A and B). We repeated the simulation experiments with these new more ‘extreme’ networks. As expected, our general classifier accuracy held across the range of regimes, with comparable performance across all three model types (Figure 1—figure supplement 2C). In addition, the reciprocal relationship between classifier uncertainty and $v$ were replicated in the fast and slow networks (Figure 2—figure supplement 1A), with the fast network showing a more expansive dynamic range of drift rates than the intermediate or slow networks. When we look at the correlation between classifier uncertainty and $v$, we again see a consistent negative association across parameter regimes (Figure 2—figure supplement 1B). The variability of this effect increases when networks have faster response times, suggesting that certain parameter regimes increase overall behavioral variability. Despite this, our key simulation effects appear to be robust to variation in parameter scheme.

To test the predictions of our model, a sample of humans (N = 4) played a dynamic two-armed bandit task under experimental conditions similar to those used for the simulated CBGT network and prior behavioral work (Bond et al., 2021) as whole brain hemophysiological signals were recorded using functional magnetic resonance imaging (fMRI) (Figure 3—figure supplement 1). On each trial, participants were presented with a male and female Greeble (Gauthier and Tarr, 1997). The goal was to select the Greeble most likely to give a reward. Selections were made by pressing a button with their left or right hand to indicate the left or right Greeble on the screen.

Crucially, we designed this experiment such that each participant acted as an out-of-set replication test, having performed thousands of trials individually. Specifically, to ensure we had the statistical power to detect effects on a participant-by-participant basis, we collected an extensive data set comprising 2700 trials over 45 runs from nine separate imaging sessions for each of four participants. Consequently, we amassed a grand total of 36 hr of imaging data over all participants, which was used to evaluate the replicability of our findings at the participant-by-participant level. Therefore, our statistical analyses were able to estimate effects on a single-participant basis.

Behaviorally, human performance in the task replicated our prior work (Bond et al., 2021). Both response speed and accuracy changed across conditions in a way that matched what we observed in Experiment 2 in Bond et al., 2021. Specifically, we see a consistent effect of change point on both RT and accuracy that matches the behavior of our network (Figure 3—figure supplement 2). To address how a change in the environment shifted underlying decision dynamics, we used a hierarchical DDM modeling approach (Wiecki et al., 2013) as we did with the network behavior (see ‘Materials and methods’ for details). Given previous empirical work (Bond et al., 2021) and the results from our CBGT network model showing that only $v$ and, less reliably, $a$ respond to a shift in the environment, we focused our subsequent analysis on these two parameters. We compared models where single parameters changed in response to a switch, pairwise models where both parameters changed, and a null model that predicts no change in decision policy (Supplementary files 1 and 2). Consistent with the predictions from our CBGT model, we found equivocal fits for the model mapping both $\mathrm{\Delta }B$ to $v$ and $\mathrm{\Omega }$ to $a$ and a simpler model mapping $\mathrm{\Delta }B$ to $v$ (see Supplementary file 1 for average results). This pattern was fairly consistent at the participant level, with 3/4 participants showing $\mathrm{\Delta }B$ modulating $v$ (Supplementary file 2). These results suggest that as the belief in the value of the optimal choice approaches the reward value for the optimal choice, the rate of evidence accumulation increases.

Taken altogether, we confirm that humans rapidly shift how quickly they accumulate evidence (and, to some degree, how much evidence they need to make a decision) in response to a change in action–outcome contingencies. This mirrors the decision parameter dynamics predicted by the CBGT model. We next evaluated how this change in decision policy tracks with competition in neural action representations.

To measure competition in action representations, we first needed to determine how individual regions (i.e., voxels) contribute to single decisions. For each participant, trial-wise responses at every voxel were estimated by means of a general linear model (GLM), with trial modeled as a separate condition in the design matrix. Therefore, the ${\widehat{\beta }}_{t,v}$ estimated at voxel $v$ reflected the magnitude of the evoked response on trial $t$. As in the CBGT model analysis, these whole-brain, single-trial responses were then submitted to a LASSO-PCR classifier to predict left/right response choices (Figure 3—figure supplement 3). The performance of the classifier for each participant was evaluated with a 45-fold cross-validation, iterating through all runs so that each one corresponded to the hold-out test set for one fold.

Our classifier was able to predict single-trial responses well above chance for each of the four participants (Figure 3A and B), with mean prediction accuracy ranging from 65 to 83% (AUCs from 0.72 to 0.92). Thus, as with the CBGT network model, we were able to reliably predict trial-wise responses for each participant. Figure 3C shows the average encoding map for our model as an illustration of the influence of each voxel on our model predictions (Figure 3—figure supplement 4 displays individual participant maps). These maps effectively show voxel-tuning toward rightward (blue) or leftward (red) responses. Qualitatively, we see that cortex, striatum, and thalamus all exhibit strongly lateralized influences on contralateral response prediction. Indeed, when we average the encoding weights in terms of principal CBGT nuclei (Figure 3D), we confirm that these three regions largely predict contralateral responses. See Figure 3—figure supplement 4 for a more detailed summary of the encoding weights across multiple cortical and subcortical regions.

Figure 3 with 4 supplements see all

Download asset Open asset

These results show that we can reliably predict single-trial choices from whole-brain hemodynamic responses for individual participants. Further, key regions of the CBGT pathway contribute to these predictions. Next, we set out to determine whether competition between these representations for left and right actions correlates with changes in the drift rate, as predicted by the CBGT network model (Figure 2C).

To evaluate whether competition between action channels correlates with the magnitude of $v$ on each trial, as the CBGT network predicts (Figure 2C), we focused our analysis on trials surrounding the change point, following analytical methods identical to those described in the previous section and shown in Figure 2C.

Consistent with the CBGT network model predictions, following a change point, $v$ shows a stereotyped drop and recovery as observed in the CBGT network (Figure 2C) and prior behavioral work (Bond et al., 2021; Figure 4A). This drop in $v$ tracked with a relative increase in classifier uncertainty, and subsequent recovery, in response to a change in action–outcome contingencies (mean bootstrapped ${\displaystyle \beta :-0.021}$ to ${\displaystyle -0.001;t}$ range: $-3.996$ to ${\displaystyle -1.326;p{s}_1=0.057}$, ${\displaystyle p_{S2}<0.001}$, ${\displaystyle p_{\text{All}}<0.001}$). As with the CBGT network simulations (Figure 2D), we also observe a consistent negative correlation between $v$ and classifier uncertainty over all trials, irrespective of their position to a change point, in each participant and in aggregate (Figure 4B; Spearman’s $\rho$ range: $-0.08$ to ${\displaystyle -0.04;p}$ range: ${\displaystyle <0.001}$ to 0.043, see Figure 4—figure supplement 1 for null effect on $a$).

Figure 4 with 1 supplement see all

Download asset Open asset

These results clearly suggest that, as predicted by our CBGT network simulations and prior work (Dunovan et al., 2019; Vich et al., 2022; Rubin et al., 2021), competition between action representations drives changes in the rate of evidence accumulation during decision-making in humans.

As the world changes and certain actions become less optimal, successful behavioral adaptation requires flexibly changing how sensory evidence drives decisions. Our simulations and hemophysiological experiments in humans show how this process can occur within the CBGT circuits. Here, a shift in action–outcome contingencies induces competition between encoded action plans by modifying the relative balance of direct and indirect pathway activity in CBGT circuits, both within and between action channels, slowing the rate of evidence accumulation to promote adaptive exploration.

If the environment subsequently remains stable, then this learning process accelerates the rate of evidence accumulation for the optimal decision by increasing the strength of action representations for the new optimal choice. This highlights how these macroscopic systems promote flexible, effective decision-making under dynamic environmental conditions.

We simulated neural dynamics and behavior using a biologically based, spiking CBGT network model (Dunovan and Verstynen, 2019; Vich et al., 2022). The network representing the CBGT circuit is composed of nine neural populations: cortical interneurons (CxI), excitatory cortical neurons (Cx), striatal D1/D2-spiny projection neurons (dSPNs/iSPNs), striatal fast-spiking interneurons (FSI), the internal (GPi) and external globus pallidus (GPe), the subthalamic nucleus (STN), and the thalamus (Th). All the neuronal populations are segregated into two action channels with the exception of cortical (CxI) and striatal interneurons (FSIs). Each neuron in the population was modeled with an integrate-fire-or-burst-model (Wei et al., 2015), and a conductance-based synapse model was used for NMDA, AMPA, and GABA receptors. The neuronal and network parameters (inter-nuclei connectivity and synaptic strengths) were tuned to obtain realistic baseline firing rates for all the nuclei. The details of the model are described in our previous work (Vich et al., 2022) as well as in the ‘Neuron model’ section below for the sake of completeness.

Corticostriatal weights for D1 and D2 neurons in striatum were modulated by phasic dopamine to model the influence of reinforcement learning on network dynamics. The details of STDP learning are described in detail in previous work (Vich et al., 2020), but key details are shown below. As a result of these features of the CBGT network, it was capable of learning under realistic experimental paradigms with probabilistic reinforcement schemes (i.e., under reward probabilities and unstable action–outcome values).

A decision between the two competing actions (‘left’ and ‘right’) was considered to be made when either of the thalamic subpopulations reached a threshold of 30 Hz. This threshold was set based on the network dynamics for the chosen parameters with the aim of obtaining realistic reaction times. The maximum time allowed to reach a decision was 1000 ms. If none of the thalamic subpopulations reached the threshold of 30 Hz, no action was considered to be taken. Such trials were dropped from further analysis. Reaction times were calculated as time from stimulus onset to decision (either subpopulation reaches the threshold). The ‘slow’ and ‘fast’ trials were categorized as reaction times ≥75th percentile (314.5 ms) and reactions time <50 th percentile (196.0 ms), respectively, of the reaction time distributions. The firing rates of the CBGT nuclei during the reaction times were used for prediction analysis as discussed in our description of single-trial response estimation below.

The corticostriatal weights are modified by a dopamine-mediated STDP rule, where the phasic dopamine is modulated by reward prediction error. The internal estimate of the reward is calculated at every trial by a Q-learning algorithm and is subtracted from the reward associated with the experimental paradigm to yield a trial-by-trial estimate of the reward prediction error. The effect of dopaminergic release is receptor dependent; a rise in dopamine promotes potentiation for dSPNs and depression for iSPNs. The degree of change in the weights is dependent on an eligibility trace, which is proportional to the coincidental presynaptic (cortical) and postsynaptic (striatal) firing rates. The STDP rule is described in detail in Vich et al., 2020 as well as in the appendix.

We follow the paradigm of a two-arm bandit task, where the CBGT network learns to consistently choose the rewarded action until the block changes (i.e., the reward contingencies switch), at which point the CBGT network relearns the rewarded action (reversal learning). Each session consists of 40 trials with a block change every 10 trials. The reward probabilities represent a conflict of (75%, 25%); that is, in a left block, 75% of the left actions are rewarded, whereas 25% of the right actions are rewarded. The inter-trial-interval in network time is fixed to 600 ms.

To maximize the similarity between the CBGT network simulations and our human data, we randomly varied the initialization of the network such that neurons with a specific connection probability were randomly chosen for each simulated subject, with the background input to the nuclei for each simulated subject as a mean-reverting random walk (noise was drawn from the normal distribution ${\displaystyle \mathcal{N}(0,1)}$). These means are listed in Supplementary file 6.

Four neurologically healthy adult humans (two female, all right-handed, 29–34 years old) were recruited and paid $30 per session, in addition to a performance bonus and a bonus for completing all nine sessions. These participants were recruited from the local university population.

All procedures were approved by the Carnegie Mellon University Institutional Review Board (Approval Code: $2018\mathrm{\_}00000195$). All research participants provided informed consent to participate in the study and consent to publish any research findings based on their provided data.

The experiment used male and female Greebles (Gauthier and Tarr, 1997) as selection targets. Participants were first trained to discriminate between male and female Greebles to prevent errors in perceptual discrimination from interfering with selection on the basis of value. Using a two-alternative forced choice task, participants were presented with a male and female Greeble and asked to select the female, with the male and female Greeble identities resampled on each trial. Participants received binary feedback regarding their selection (correct or incorrect). This criterion task ended after participants reached 95% accuracy. After reaching perceptual discrimination criterion for each session, each participant was tested under nine reinforcement learning conditions composed of 300 trials each, generating 2700 trials per participant in total. Data were collected from four participants in accordance with a replication-based design, with each participant serving as a replication experiment. Participants completed these sessions in randomized order. Each learning trial presented a male and female Greeble (Gauthier and Tarr, 1997), with the goal of selecting the gender identity of the Greeble that was most rewarding. Because individual Greeble identities were resampled on each trial, the task of the participant was to choose the gender identity rather than the individual identity of the Greeble, which was most rewarding.

Probabilistic reward feedback was given in the form of points drawn from the normal distribution ${\displaystyle \mathcal{N}(\mu =3,\sigma =1)}$ and converted to an integer. These points were displayed at the center of the screen. For each run, participants began with 60 points and lost one point for each incorrect decision. To promote incentive compatibility (Rosenzweig and Evenson, 1977), participants earned a cent for every point earned. Reaction time was constrained such that participants were required to respond within between 0.1 s and 0.75 s from stimulus presentation. If participants responded in $\le 0.1$ s, $\ge 0.75$ s, or failed to respond altogether, the point total turned red and decreased by 5 points. Each trial lasted 1.5 s and reward feedback for a given trial was displayed from the time of the participant’s response to the end of the trial. To manipulate change point probability, the gender identity of the most rewarding Greeble was switched probabilistically, with a change occurring every 10, 20, or 30 trials, on average. To manipulate the belief in the value of the optimal target, the probability of reward for the optimal target was manipulated, with p set to 0.65, 0.75, or 0.85. Each session combined one value of p with one level of volatility, such that all combinations of change point frequency and reward probability were imposed across the nine sessions. Finally, the position of the high-value target was pseudo-randomized on each trial to prevent prepotent response selections on the basis of location.

Statistical analyses and data visualization were conducted using custom scripts written in R (R Foundation for Statistical Computing, version 3.4.3) and Python (Python Software Foundation, version 3.5.5). Scripts are publicly available (Bond et al., 2023).

Binary accuracy data were submitted to a mixed effects logistic regression analysis with either the degree of conflict (the probability of reward for the optimal target) or the degree of volatility (mean change point frequency) as predictors. The resulting log-likelihood estimates were transformed to likelihood for interpretability. RT data were log-transformed and submitted to a mixed effects linear regression analysis with the same predictors as in the previous analysis. To determine if participants used ideal observer estimates to update their behavior, two more mixed effects regression analyses were performed. Estimates of change point probability and the belief in the value of the optimal target served as predictors of reaction time and accuracy across groups. As before, we used a mixed logistic regression for accuracy data and a mixed linear regression for reaction time data.

To assess whether and how much the ideal observer estimates of change point probability ($\mathrm{\Omega }$) and the belief in the value of the optimal target ($\mathrm{\Delta }\text{B}$) (Nassar et al., 2010; Bond et al., 2021) updated the rate of evidence accumulation ($v$), we regressed the change point-evoked ideal observer estimates onto the decision parameters using hierarchical drift diffusion model (HDDM) regression (Wiecki et al., 2013). These ideal observer estimates of environmental uncertainty served as a more direct and continuous measure of the uncertainty we sought to induce with our experimental manipulations. Using this more direct approach, we pooled change point probability and belief across all conditions and used these values as our predictors of drift rate and boundary height. Responses were accuracy-coded, and the belief in the difference between targets values was transformed to the belief in the value of the optimal target ($\mathrm{\Delta }{B}_{\text{optimal(t)}}=B_{\text{optimal(t)}}-B_{\text{suboptimal(t)}}$). This approach allowed us to estimate trial-by-trial covariation between the ideal observer estimates and the decision parameters.

To find the models that best fit the observed data, we conducted a model selection process using deviance information criterion (DIC) scores. A lower DIC score indicates a model that loses less information. Here, a difference of ≤ 2 points from the lowest-scoring model cannot rule out the higher scoring model; a difference of 3–7 points suggests that the higher scoring model has considerably less support; and a difference of 10 points suggests essentially no support for the higher scoring model (Spiegelhalter et al., 2002; Burnham and Anderson, 1998). We evaluated the DIC scores for the set of fitted models relative to an intercept-only regression model (DIC_intercept - DIC_modeli).

Neurologically healthy human participants (N = 4, two females) were recruited. Each participant was tested in nine separate imaging sessions using a 3T Siemens Prisma scanner. Session 1 included a set of anatomical and functional localizer sequences (e.g., visual presentation of Greeble stimuli with no manual responses, and left vs. right button responses to identify motor networks). Sessions 2–10 collected five functional runs of the dynamic two-armed bandit task (60 trials per run). Male and female ‘Greebles’ served as the visual stimuli for the selection targets (Gauthier and Tarr, 1997), with each presented on one side of a central fixation cross. Participants were trained to respond within 1.5 s.

To minimize the convolution of the hemodynamic response from trial to trial, inter-trial intervals were sampled according to a truncated exponential distribution with a minimum of 4 s between trials, a maximum of 16 s, and a rate parameter of 2.8 s. To ensure that head position was stabilized and stable over sessions, a CaseForge head case was customized and printed for each participant. The task-evoked hemodynamic response was measured using a high spatial (2 mm³ voxels) and high temporal (750 ms TR) resolution echo planar imaging approach. This design maximized recovery of single-trial evoked BOLD responses in subcortical areas, as well as cortical areas with higher signal-to-noise ratios. During each functional run, eye-tracking (EyeLink, SR Research Inc), physiological signals (ECG, respiration, and pulse oximetry via the Siemens PMU system) were also collected for tracking attention and for artifact removal.

fMRI data were preprocessed using the default fMRIPrep pipeline (Esteban et al., 2019), a standard toolbox for fMRI data preprocessing that is robust to variations in scan acquisition protocols and minimal user manipulation.

We used a univariate GLM to estimate within-participant trial-wise responses at the voxel level. Specifically, for each fMRI run, preprocessed BOLD time series were regressed onto a design matrix, where each task trial corresponded to a different column, and was modeled using a boxcar function convolved with the default hemodynamic response function given in SPM12. Thus, each column in the design matrix estimated the average BOLD activity within each trial. In order to account for head motion, the six realignment parameters (three rotations, three translations) were included as covariates. In addition, a high-pass filter (128 s) was applied to remove low-frequency artifacts. Parameter and error variance were estimated using the RobustWLS toolbox, which adjusts for further artifacts in the data by inversely weighting each observation according to its spatial noise (Diedrichsen and Shadmehr, 2005).

Finally, estimated trial-wise responses were concatenated across runs and sessions and then stacked across voxels to give a matrix, ${\widehat{\beta }}_{t,v}$, of T (trial estimations) × V (voxels) for each participant.

A machine learning approach was applied to predict left/right Greeble choices from the trial-wise responses. First, using the trial-wise hemodynamic responses, we estimated the contrast in neural activation when the participant made a left versus right selection. A LASSO-PCR classifier (i.e., an L1-constrained principal component logistic regression) was estimated for each participant according to the below procedure. We should note that the choice of LASSO-PCR was based on prior work building reliable classifiers from whole-brain-evoked responses that maximizes inferential utility (see Meissner et al., 2011). This approach is used in case of over-parameterization, as when there are more voxels than observations, and relies on a combination of dimensionality reduction and sparsity constraints to find the effective complexity of a model.

First, a singular value decomposition (SVD) was applied to the input matrix $X$:

where the product matrix $Z=US$ represents the principal component scores, that is, the values of $X$ projected into the principal component space, and $V^T$ an orthogonal matrix whose rows are the principal directions in feature space. Then the binary response variable $y$ (left/right choice) was regressed onto $Z$, where the estimation of the $\beta$ coefficients is subject to an L1 penalty term $C$ in the objective function:

where $\beta$ and Z include the intercept term, $y_i=\{-1,1\}$, and N is the number of observations.

Projection of the estimated $\widehat{\beta }$ coefficients back to the original feature (voxel) space was done to yield a weight map $\widehat{w}=V\widehat{\beta }$, which in turn was used to generate final predictions $\widehat{y}$:

where $x$ denotes the vector of voxel-wise responses for a given trial (i.e., a given row in the $X$ matrix). When visualizing the resulting weight maps, these were further transformed to encoded brain patterns. This step was performed to aid in correct interpretation in terms of the studied brain process because doing this directly from the observed weights in multivariate classification (and regression) models can be problematic (Winkler et al., 2015).

Here, the competition between left–right neural responses decreases classifier decoding accuracy as neural activation associated with these actions becomes less separable. Therefore, classifier prediction serves as a proxy for response competition. To quantify uncertainty from this, we calculated the Euclidean distance of these decoded responses $\widehat{y}$ from the statistically optimal choice on a given trial, $opt\mathrm{\_}choice$. This yielded a trial-wise uncertainty metric derived from the decoded competition between neural responses.

The same analytical pipeline was used to calculate single-trial responses for simulated data with a difference that trial-wise average firing rates of all nuclei from the simulations were used instead of fMRI hemodynamic responses.

To test whether our key effects were robust to variation in parameter schemes, 300 networks were sampled using Latin Hypercube Sampling (LHS), as also described in Vich et al., 2022. From this, we chose two network configurations with biologically plausible parameters, one with slower and faster reaction times (‘fast’ and ‘slow’ networks; upper and lower quartile of the reaction time distributions). We then repeated our key analyses for these two network configurations as shown in Figure 2—figure supplement 1. We show that these adaptive network effects are robust over a range of parameter configurations, so long as the network generates firing rates are biologically plausible.

We used an integrate-and-fire-or-burst model that models the membrane potential $V(t)$ as

where $g_L$ represents the leak conductance, $V_L$ is the leak reversal potential and the first term $g_L(V(t)-V_L)$ is the leak current; the next term is a low threshold $C{a}^{2+}$ current with maximum conductance $g_T$, gating variable $h(t)$, and reversal potential $V_T$, which activates when ${\displaystyle V(t)>V_h}$ due to the Heaviside function ${\displaystyle H;I_{syn}}$ is the synaptic current and $I_{ext}$ is the external current. This neuron model is capable of producing post-inhibitory bursts, regulated by the gating variable that decays with the time constant ${\tau }_h^{-}$, when the membrane potential reaches a certain threshold $V_h$ and rises with time constant ${\tau }_h^{+}$. However, when $g_T$ is set to zero, the model reduces to a leaky integrate and fire neuron. Currently, we model GPe and STN neuronal populations with bursty neurons and the remaining neuronal populations with leaky integrate-and-fire neurons, all with conductance-based synapses.

The synaptic current $I_{syn}(t)$ consists of three components, two excitatory currents corresponding to AMPA and NMDA receptors and one inhibitory current corresponding to GABA receptors, and is calculated as below:

where $g_{\text{i}}$ represents the maximum conductance corresponding to the receptor $i\in$ {AMPA, NMDA, GABA}, $V_I$ and $V_E$ represent the excitatory and inhibitory reversal potentials, and $s_{\text{i}}$ represents the gating variable for each current, with dynamics given by

The gating variables for AMPA and GABA act as leaky integrators that are increased by all incoming spikes, with an additional constraint for NMDA that ensures that the maximum value of $s_{\text{NMDA}}$ remains below 1.

The values of neuronal parameters for all of the nuclei are listed in Supplementary file 3, external inputs to the CBGT nuclei are listed in Supplementary file 4, synaptic parameter values are listed in Supplementary file 5, connectivity types and probabilities are listed in Supplementary file 6, and the number of neurons in each CBGT population is shown in Supplementary file 7.

The plasticity rule we use is a dopamine-modulated STDP rule also described in Vich et al., 2020. All the values of the relevant parameters are listed in Supplementary file 8. The weight update of a corticostriatal synapse is controlled by three factors: (1) an eligibility trace, (2) the type of the striatal neuron (iSPN/dSPN), and (3) the level of dopamine.

To compute the eligibility $(E)$ for a given synapse, an activity trace of each neuron in the presynaptic and postsynaptic populations is tracked via the equations

where $X_{PRE},X_{POST}$ are spike trains, such that $A_{PRE}$ and $A_{POST}$ maintain a filtered record of synaptic spiking of the pre/post neuron, respectively, with spike impact parameters ${\mathrm{\Delta }}_{PRE},{\mathrm{\Delta }}_{POST}$ and time constants ${\tau }_{PRE},{\tau }_{POST}$.

If the postsynaptic spike follows the spiking activity of the presynaptic population closely enough in time, then the eligibility variable $E$ increases and allows for plasticity to occur. On the other hand, if a presynaptic spike follows the spiking activity of the postsynaptic population, then $E$ decreases. In the absence of any activity and spikes, the eligibility trace decays to zero with a time constant ${\tau }_E$. Putting these effects together, we obtain the equation

The synaptic weight update depends on the dopamine receptor type of the striatal neuron; that is, if the neuron is a dSPN or iSPN. We assume that a phasic dopamine release promotes long-term potentiation (LTP) in dSPNs and long-term depression (LTD) in iSPNs. This factor is indicated by the learning rate parameter ${\alpha }_w$, which is set to a positive value for dSPNs and a negative value for iSPNs. The weight update dynamics is given by

where ${\displaystyle X\in }$ {dSPN, iSPN} with ${\displaystyle {\alpha }_{w-dSPN}>0}$ and ${\displaystyle {\alpha }_{w-iSPN}<0}$. Here, the weights of the corticostriatal synapses are bounded between the maximal value $W_{max}^X$, which depends on the SPN type, and a minimal value of $W_{min}=\mathrm{\hspace{0.17em}0.001}$. The precise values used for all relevant parameters are listed in Supplementary file 3.

In the weight update rule (6), $K_{DA}$ represents the dopamine level present. This quantity changes as a result of phasic release of dopamine (increments of size $D{A}_{inc}$), which is correlated to the reward prediction error encountered in the environment. We define a parameter $C_{scale}$ that sets the scaling between the reward prediction error and the amount of dopamine released, and $K_{DA}$ obeys the equation

where

for reward $r(t)$ and expected value $Q_{chosen}(t)$ of the chosen action. Trial-by-trial estimates of the values of the actions (left/right) are maintained by a simple Q-update rule:

where ${\displaystyle a\in }$ {left, right} and where ${\alpha }_q$ represents the learning rate of the Q-values.

Finally, the function $f_X(K_{DA})$ converts the level of dopamine into an impact on plasticity in a way that depends on the identity $X$ of the postsynaptic neuron, as follows:

where $c$ sets the dopamine level at which $f_{iSPN}$ reaches half-maximum. Supplementary file 8 lists the specific parameters used for the STDP rule.

The funders had no role in study design, data collection and interpretation, or the decision to submit the work for publication.

We thank all members of the CoAx Lab and collaborators for their feedback on the development of this work. This work was funded by the Air Force Research Laboratory, Grant Reference Number FA9550-18-1-0251 and NIH award R01DA053014 as part of the CRCNS program. All neuroimaging data were collected at the Carnegie Mellon-University of Pittsburgh Brain Imaging, Data Generation, and Education Center (RRID:SCR 023356).

All procedures were approved by the Carnegie Mellon University Institutional Review Board (Approval Code: 2018_00000195). All research participants provided informed consent to participate in the study and consent to publish any research findings based on their provided data.

© 2023, Bond et al.

This article is distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use and redistribution provided that the original author and source are credited.