State and rate-of-change encoding in parallel mesoaccumbal dopamine pathways

Abstract

The nervous system uses fast- and slow-adapting sensory detectors in parallel to enable neuronal representations of external states and their temporal dynamics. It is unknown whether this dichotomy also applies to internal representations that have no direct correlation in the physical world. Here we find that two distinct dopamine (DA) neuron subtypes encode either a state or its rate-of-change. In mice performing a reward-seeking task, we found that the animal's behavioral state and rate-of-change were encoded by the sustained activity of DA neurons in medial ventral tegmental area (VTA) DA neurons and transient activity in lateral VTA DA neurons, respectively. The neural activity patterns of VTA DA cell bodies matched DA release patterns within anatomically defined mesoaccumbal pathways. Based on these results, we propose a model in which the DA system uses two parallel lines for proportional-differential encoding of a state variable and its temporal dynamics.

Extended Data Figure 1.. Proportional-differential encoding in the peripheral and central nervous systems.

(a) Schematic showing proportional and rate-of-change encoding of a stimulus (e.g., pressure on skin) in the peripheral nervous system. (b) Schematic showing state value and rate-of-change encoding in the central nervous system. The state value and its rate-of-change during a Pavlovian conditioning task in which a reward is partially predicted (e.g., as in the behavioral task used for experiments in Figures 1 and 2). Note that in this schematic the ‘state value’ refers to both expected reward as well as the reward itself.

Extended Data Figure 2.. Neuropixels-based opto-tagging of dopamine neurons.

(a) Anatomical locations of the recorded VTA cells. Although Neuropixels probes allowed to record more than n = 186 cells (n = 6 mice), we only included cells that were simultaneously recorded on the same or neighboring electrodes as the opto-tagged cells. Thus, only VTA cells and not cells from adjacent brain regions were included. Different colors indicate locations of recorded cells in mVTA or lVTA and whether the cells were identified as dopaminergic (‘tagged’) or non-dopaminergic (‘non-tagged’). (b) Representative single-channel recording showing 3 action potentials (APs) of the same unit (red arrows). Two of these APs are induced by a 1 ms light pulse (blue). Insert shows how the AP, and the light pulse are separated in time (scale bar 10 ms). (c) Mean delay, mean delay-variation to the first evoked AP, and waveform correlation between evoked and spontaneous APs for n = 67 opto-tagged cells (n = 6 mice). (d) Peri-event histogram of an opto-tagged cell following a 1 ms light pulse. ‘P’ indicates firing fidelity. Total firing fidelity is the sum of the bins between 0 and 5 ms. (e) Peri-event histogram of the same cell as in (d) but following a 2 ms light pulse, which increased firing fidelity, but at the risk of adulterating the waveforms with a light-offset artifact. (f) Evoked (orange) and spontaneous (blue) waveforms for all opto-tagged cells. The waveform on the channel with the highest amplitude was chosen. (g) Peri-event histogram of an opto-tagged cell following a light pulse at 1% of the max laser power (max laser power: 25 mW at tip of optical fiber at a distance of ~1 mm from the recording site). (h) Same unit as in (g), but laser power was set to 10% of max. (i) Same unit as in (g, h), but following a 1 ms light pulse with laser at max intensity. Note, increase in firing fidelity when comparing (g, h, i). Firing fidelity was always < 1 (i.e., one light pulse (1 ms) did not evoke more than a single AP). (j) AP width and baseline firing frequency (FF) for all recorded cells (n = 186 cells). Insert: Sample APs showing calculation of waveform width (i.e., difference between AP peak and hyperpolarization peak). (k) Mean FF (left) and mean waveform widths for dopamine (tagged; n = 67 cells) and non-dopamine cells (non-tagged; n = 119 cells). Individual data points presented in (j). Significance was calculated by means of two-sided Welch t-test; frequency ***P<0.001. Opto-tagged cells had a significantly lower baseline FF and a longer waveform width compared to non-opto-tagged cells. (l) Mean FF (left) and mean waveform widths for dopamine neurons located in mVTA (n = 40 cells; orange) and lVTA (n = 27 cells; blue). Individual data points shown in (j). Significance was calculated by means of two-sided Welch t-test. All data represented as mean ± SEM.

Extended Data Figure 3.. Four distinct activity patterns describe neural activity of VTA cells during behavior.

(a) Sparce PCA embedding of all recorded VTA cells (n = 186 cells; n = 6 mice); color coding indicates that four different activity patterns were obtained using k-means clustering. (b) Graph showing action potential (AP) width and baseline firing frequency of all recorded VTA cells (n = 186 cells); color coding indicates cluster identity (cluster #1: grey, cluster #2: green, cluster #3: blue, cluster #4: red). (c) Bar graphs showing mean baseline firing frequency (left) and action potential width (right) for all recorded VTA cells (n = 186 cells) based on cluster identity. Significance was calculated by means of 1way ANOVA. The annotations ‡ and # indicate significant difference in Tukey’spost hoc testing; ***P<0.001. Data represented as mean ± SEM. (d) Percentage of non-opto-tagged (i.e., non-dopaminergic neurons) neurons in lVTA (left) and mVTA (right) based on cluster identity. (e) Box plot (median, quantiles and outliers) depicting the anatomical location of all opto-tagged (i.e., DA) neurons along a 45° (ventro-medial to dorso-lateral) axis trough the VTA. Cluster #2 DA neurons were significantly more likely to be located in dorso-lateral part of the VTA compared to DA neurons in the other clusters. Cluster 1: n = 100 cells, Cluster 2: n = 33 cells, Cluster 3: n = 35 cells, Cluster 4: n=18 cells. Significance was calculated by means of 1way ANOVA; ***P<0.001. Data represented as median, interquartile range (shading) and total range (error bars).

Extended Data Figure 4.. VTA DA cell body activity is sufficient to explain DA release patterns in distinct NAc subregions.

(a) AAV-DIO-ChR2-mCherry was injected bilaterally into the VTA of DAT-Cre mice. The same mice were injected with dLight into NAcLat and NAcMed. Optical fibers were implanted dorsal to the VTA, in the NALat and NAcMed (optical fibers in NAcMed and NAcLat were counterbalanced across all animals, n = 6 mice). (b) 100 trials of 1 ms 470 nm light stimulation of VTA DA neurons were used to calculate DA release (auROC) in NAcLat (left) and NAcMed (right). Blue line indicates 1 ms light pulse in VTA. (c) Response kernels in (b) were used to convolute single-unit recordings. Top: Sample unit, each vertical line denotes a single action potential. Bottom: Corresponding inferred dLight trace when the sequence of action potentials is convoluted with the NAcLat kernel in (b). (d) Mean normalized (auROC) single-unit recordings of all cluster #2 neurons in the lVTA; color code indicates random interval length. CS+ onset was at t = 0 sec. (e) Activity traces of the neurons in (d) were convoluted as in (c). Trials were normalized (auROC) and averaged over all cells as in (d). (f) Recorded NAcLat DA release patterns normalized (auROC) and averaged across all mice (n = 6 mice). The 2 sec interval (orange) is the same as in Fig. 2j. Note the qualitative and quantitative similarities with (e), including the below-baseline reduction in the interval between CS+ onset and reward delivery. (g-i) Same as in (d-f), but for NAcMed. Note that this experiment does not unequivocally prove that mVTA cluster #3 single-unit activity is responsible for NAcMed DA release patterns. However, it does suggest that mVTA cluster #3 single-unit activity is sufficient to explain NAcMed DA release patterns. All data represented as mean ± SEM (shading).

Extended Data Figure 5.. A linear model to predict dopamine transients across all NAc subregions.

(a) Each NAc recording was normalized, and a linear regression model was fitted to the discrete task events to predict NAc DA release. Left: Animal behavior was binned to 100 ms and shifted in time (1 sec in both directions, 20 coefficients were fitted for each task event). Task events were encoded as ‘1’ or ‘0’ depending on occurrence in any given time bin except for the trial duration, which was coded from 0 to 1 over the 0–2500 ms interval. Right: Data were split into a training (80%) and test (20%) set. The task events in the training set were fitted to the dLight data using linear regression. The sample trace shows a high overlap between recorded response (blue) and model prediction (orange). Bar graph shows model performance on training and test set (mean ± SEM, n = 29 recordings). Note that the model performed notably better in NAcLat; probably because it performed well at fitting the large transients at CS+ onset and reward delivery. Lower right: Recorded data (left) and model prediction (right); color coded based on trial duration (averaged over all trials). (b) Task-event kernels used to fit the behavior to the fiber photometry signal. (c) AUCs of the kernels in (b), hexagons denote individual recording locations. The color map was selected to highlight the direction (up or down) of the kernels (aca: anterior commissure). (d) Because the anatomical location of individual recordings was histologically verified, we used Gaussian interpolation to calculate the task-event kernels for each pixel in a coronal section of the NAc. Note the gradient from NAcMed to NAcLat with NAcMed showing positive AUC for kernels related to reward seeking, reward consumption and trial duration. Conversely, NAcLat kernels were specifically tuned to trial start (CS+ onset) and reward delivery. The NAc core represents an intermediate structure with signals tuned either towards NAcMed or NAcLat based on the position more medial or more lateral to the aca, respectively. All data represented as mean ± SEM (error bars or shading).

Extended Data Figure 6.. Topography of the mesoaccumbal dopamine system in mice and rats.

(a) Schematics showing representative fluorescent retrobeads locations (green) in different NAc subregions (NAcMed, NAc core, NAcLat) along the ventro-medial to dorso-lateral axis for different mice (n = 3 mice for each injection location; aca: anterior commissure). (b) Corresponding sample fluorescent images of coronal brain sections from the NAc (tyrosine hydroxylase (TH): red). (c) Sample fluorescent images of coronal midbrain sections from the same mice shown in (a, b). (d) Corresponding schematics highlighting locations of retrogradely labeled (i.e., beads-positive, green) neurons in the medial VTA (mVTA) and lateral VTA (lVTA) (IF: interfascicular nucleus, PN: paranigral nucleus, PBP: parabrachial pigmented nucleus, IPN: interpeduncular nucleus; scale bars (b, c) 500 μm). (e) Schematic showing location of injection sites for fluorescent retrobeads (red) in the NAcLat and FluoroGold (grey) in the NAc core (n = 3 rats; CPu: Caudate Putamen). (f) Sample fluorescent images showing coronal brain sections of the NAcLat and NAc core at two different anterior-posterior coordinates (TH: blue, beads: red, FluoroGold: white; Scale bar 1 mm). (g) Fluorescent images showing coronal brain sections of the VTA from a rat that was injected with fluorescent retrobeads (red) into the NAcLat and FluoroGold (white) into the NAc core (scale bar 500 μm). Note, that mice and rats show a similar anatomical topography of NAcLat- and NAc core-projecting DA neurons (compare panels (c) and (g)). In both cases, NAcLat-projecting DA neurons are located in the lVTA, while NAc core-projecting DA neurons are located in the mVTA).

Extended Data Figure 7.. Bayesian model selection showing that the behavior of the animals can be best described by a Q learning model with a single learning rate.

(a) Results of Bayesian model selection. Model 2, a Q learning model with a single learning rate, has the highest protected exceedance probability, and is therefore the best descriptor of the data of the three given models. Model 1 is a random choice model (in which every choice has a probability of 0.5); model 3 is a model with separate learning rates for learning from reward versus reward omission. (n = 10 mice) (b) Best-fit model parameters for the experimental mice shown in Fig. 3; each dot indicates an individual mouse/session. The two panels show the relation between learning rate α (left) and explore/exploit parameter β (right) versus the number of achieved reversal per 100 trials; lines indicate least-squares linear fit and its 95% confidence intervals; significance was tested by means of two-sided t-test.

Extended Data Figure 8.. Dopamine release in NAcLat, but not NAcMed, encodes negative RPE.

(a) Graphs showing licking behavior of a representative sample animal during the task. Note that licking generally started ~500 ms after CS+ onset with only minimal variability. Inset shows schematic of experimental design: Head fixed mice were trained on a behavioral assay, in which the presence of a CS+ predicted the availability of a reward (3 μl of a 1% sucrose solution). The reward was delivered to the animal after one single lick (i.e., it assumes minimal effort to the animal to obtain a reward). After training, mice were subjected to a session in which 10% of the trials were not rewarded (omission trials). DA release was measured simultaneously in the NAcLat and NAcMed using dLight1.2. (b) Top: Heatmaps showing auROC normalized dLight responses in the NAcLat and NAcMed for all animals (n = 9 mice) in standard trials in which reward delivery occurred after the first lick during CS+ presentation. Bottom: mean DA response in NAcLat (blue) and NAcMed (orange) as well as the mean lick rate (dashed line) for all mice. Note that NAcLat DA release tracks CS+ onset and reward delivery, while NAcMed DA tracks licking behavior, including reward consumption. (c) Same as in panel (b) but for omission trials in which no reward was delivered. Note that NAcMed DA release returns to baseline concurrent with the termination of licking behavior. Conversely, NAcLat DA release showed a negative deflection below baseline level during reward omission, which reflects a negative RPE. (d) Quantification (box plot, median and quantiles) of panels (b) and (c). Significance was calculated by means of a two-sided Students t-test; ***P<0.001, n = 9 mice. All data, except panel (d), is represented as mean ± SEM (shading).

Extended Data Figure 9.. Dopamine release patterns for all mice used in experiment shown in Figures 2g–o.

NAcLat (left), NAcMed (middle) DA release was Z-score normalized and trials were sorted by interval length. The derivative of NAcMed DA release (dNAcMed/dt, right) was calculated by determining the slope of a line that was fitted (linear regression) over a 500 ms interval at each time point. Mouse #4 is the sample animal shown in Fig. 2i; mouse #1 is the sample animal shown in Fig. 2o. Although there is substantive similarity in the pattern of DA release in the NAcLat with the temporal derivative of DA release in the NAcMed, there are two important limitations with such a comparison. First, the release of DA in the NAc reflects a complex response to activity in DA cell bodies in the VTA. As shown in Extended Data Fig. 4, NAc DA release is a convolution of VTA cell body activity. Second, there is significant between animal variation. For example, mouse #6 lacks a transient increase in the dNAcMed/dt trace in response to reward delivery. Additionally, while most mice show a dip in NAcLat DA release below the baseline between the onset of CS+ and reward delivery, only mouse #2 and mouse #5 exhibit a downward slope in NAcMed DA release during the same period (see also Fig. 2j, n).

Extended Data Figure 10.. Dopamine release patterns in NAcMed and NAcLat in response to aversive stimuli can be conceptualized in the context of state and rate-of-change encoding.

(a) Hypothesis: the decrease in NAcLat DA release in response to a mild electric tail shock reflects the derivative of the downward slope in NAcMed DA release. If this is true, a sustained decrease in NAcLat DA release should be preceded by a brief increase that reflects the derivative of the onset of DA release in the NAcMed. (b) Experimental data: Head-fixed mice (n = 12 mice) were subjected to 10 unpredictable tail shocks while DA release in NAcMed and NAcLat was recorded using dLight 1.3. Left: Sustained increase in DA release in the NAcMed. Right: Sharp increase at shock onset that is followed by sustained below-baseline reduction of NAcLat DA release. While the gradual downward slope in the NAcMed and the transient increase in the NAcLat at shock onset suggest that our hypothesis is correct, there are two limitations that may be due to variability between experimental animals. First, the NAcMed DA response appears to be delayed compared to the NAcLat DA increase. Second, there is a sharp dip in NAcLat DA release that did not directly correspond to a sharp decrease in NAcMed DA release. All data represented as mean ± SEM (shading).

Figure 1.. Four distinct activity patterns describe neural activity of VTA DA and non-DA neurons during behavior.

(a) Experimental design. (b) Left: Sample recording site (scale bar 0.5 mm). Right: Recording location of two DA cells in different VTA subregions (Unit ‘A’ in lVTA (blue), Unit ‘B’ in mVTA (orange)). (c) Spontaneous and evoked waveforms of the two simultaneously recorded cells shown in (b). (d) Reward-seeking task. Top: all CS+ trials (sorted by interval length) of a sample session (green: licks during reward seeking; red: licks during reward delivery; blue: licks during reward consumption). Bottom: task structure. (e) Peri-event histograms of the neural activity pattern of unit ‘A’ and ‘B’ shown in (b) during behavior. Interval lengths are color coded and grouped together in bins of 300 ms. (f) Normalized (auROC) activity patterns of each VTA cell (n=186 cells, n=6 mice) sorted based on cluster identity and task component. Right: Mean normalized activity patterns of the four clusters during behavior color coded by interval length. Vertical blue bar denotes DA neurons. Bottom: Mean normalized activity responses during different task components for the four clusters. (g) Anatomical location of opto-tagged VTA DA neurons (color coded based on cluster identity). (h, i) Percentage of opto-tagged DA neurons in the four clusters (C1-C4) based on location in lVTA (h) or mVTA (i). (j) Antidromic opto-tagging. Insert: mean delay of light response (n=19 cells, n=7 mice). (k) Sample of collision test. Top: collision trials in which a spontaneous action potential (AP) occurred in the 15 ms leading up to the timepoint of the mean light response. Note how presence of spontaneous (i.e., presumably orthodromic) APs block the light-evoked AP. Bottom: trials in which no spontaneous AP occurred in the 15 ms leading up to the acquisition of the light response (i.e., ‘free’ trials). (l) Normalized (auROC) responses of all NAcLat-projecting DA neurons (n=5 cells) during cue onset (left) and reward delivery (right). Note the homogeneous response patterns across different cells. (m) Same as in (l) but for NAcMed-projecting DA neurons (n=14 cells). All data represented as mean ± SEM (shading).

Figure 2.. VTA DA cell body activity resembles DA release within defined mesoaccumbal subsystems.

(a) Experimental design; same task as in Figure 1. (b) Top: GCaMP6m (green) expression and recording site in mVTA (scale bar 1 mm); bottom: Recording locations in different VTA subregions for all recorded animals (n=5 mice). (c) Normalized photometry signals obtained during the 2.2–2.4 sec interval between CS+ onset and reward delivery at five different VTA recording sites (numbers and colors correspond to (b)). (d) lVTA DA cells: Transient ‘peaks’ at CS+ onset and reward delivery; mVTA DA cells: sustained (‘ramping-like’) activity patterns. (e) Left: Behavioral responses. Right: Calcium transients of lVTA DA neurons. Top: Individual trials (sorted by interval length); bottom: Mean plots grouped by interval length; different colors represent different interval lengths (bin width: 0.3 s). (f) Same as in (e), but for mVTA DA cells. (g) Top: Simultaneous recordings in NAcMed and NAcLat; bottom: dLight expression and recording sites in NAcLat and NAcMed (scale bar 200 μm, representative example from n=6 mice). (h, i) As in (c, f) but for DA release measurements in NAcLat and NAcMed. (j) Group level analysis of the auROC-normalized 2.2–2.4 sec interval for NAcLat and NAcMed (n=6 mice, auROC calculated over 1 session (200 trials), 1 auROC trace for each mouse was included). (k) Box plot (median and quantiles) showing significant differences in mean DA response between NAcMed and NAcLat during 1 sec interval before reward delivery. Significance was calculated by means of two-sided paired t-test; ** P=0.004, n=6 mice. (l) dLight response in NAcLat and NAcMed to CS onset in ‘Go’ trials (left, n=6 mice), ‘No-Go’ trials (middle, n=5 mice) and following CS− onset (n=6 mice). (m) Top: Peri-event response to lick-bout onset (~500 ms after CS+ onset). Dotted line indicates lick rate. Bottom: Rate-of-change of the NAcMed response during lick-bout onset. (n) Same as in (j), but for the rate-of-change response in NAcMed (n=6 mice). (o) Sample simultaneous recordings in NAcLat and NAcMed. Insert: Single-trial response highlighting how the slope of the NAcMed trace is reflected in the NAcLat trace (green: reward seeking licks, blue: reward consumption licks). All data represented as mean ± SEM (error bars or shading).

Figure 3.. Separate mesoaccumbal DA pathways encode value and RPEs.

(a) Two-armed bandit task. (b) Performance of individual mice (after training). (c) Sample dLight (green) recording in NAcLat. Reward expectation (blue) and RPE (red) were extracted from trial-by-trial data using a Q learning model fitting procedure (see Extended Data Fig. 7 for model selection). (d) Top: Peri-stimulus time histograms of individual Z-score normalized DA release dynamics during the three main phases of the task. Bottom: Mean response of all mice showing differential DA dynamics in NAcLat during rewarded and unrewarded trials. (e) Same as (d), but for DA release in NAcMed. (f) Mean dLight Z-score 0.5 s after choice (i.e., after nosepoke). Circles indicate individual animals (n=5 mice per group). Significance was calculated by means of two-tailed Student’s t-test different from 0; P values for the four comparisons in order of appearance: 0.0011, 0.0582, 0.0014, 0.0276. (g) Median latency until peak dLight Z-score and median latency until first lick in rewarded trials. NAcLat: DA peak is observed before the animals made first lick; NAcMed: DA peak coincides with first lick. Circles indicate individual animals (n=5 mice per group). Significance was calculated by means of Student’s paired t-test; P=0.0296 for NAcLat, P=0.5527 for NAcMed. (h) Regression plots around the time of choice showing the fraction of trial-to-trial variance of dLight signal that could be explained by trial-to-trial modeled reward expectation. (i) Mean NAcMed DA release during trials with high or low reward expectation highlighting (arrow) the correlation depicted in (h). In trials with high reward expectation, NAcMed DA release is increased in the interval preceding the side-in response. (j) Quantification of the data shown in (i, h). Trial-by-trial variation in NAcMed DA release tracks reward expectation in the interval before choice (n=5 mice per group). Asterisks indicate significant Sidak’s post-hoc test after a significant value component × brain area interaction effect in a two-way RM ANOVA; P<0.0001 for value, P=0.9931 for RPE. (k) Regression plots around the time of choice showing the fraction of trial-to-trial variance of the dLight signal that could be explained by the trial-to-trial modeled RPE. (l) Mean NAcLat DA release highlighting the correlation depicted in panel (k). NAcLat DA release tracks both positive and negative RPE immediately after the choice outcome is made known to the animal. (m) Quantification of data shown in (k) (n=5 mice per group). Asterisks indicate significant Sidak’s post-hoc test after a significant value component × brain area interaction effect in a two-way RM ANOVA; P=0.7318 for value, P<0.0001 for RPE. All data represented as mean ± SEM (error bars or shading).

Figure 4.. State and rate-of-change are encoded by different VTA DA subpopulations.

(a) Schematic of model. The task events of the last 500 ms were fitted to a one-dimensional state vector using linear regression. This state vector and its temporal derivative were fitted to a set of simultaneously recorded neurons using a formula with two free parameters: ‘A’ controlled the amplitude of the predicted activity pattern, ‘B’ controlled the bias of the best fit towards the state vector or its temporal derivative (the encoding bias). (b) Differential distribution of Z-score normalized neural activity above and below the pre-trial baseline suggests implementation of Leaky ReLU in our model. (c) Sample showing state and Δstate/Δt traces averaged over all trials with an interval of 2±0.15 s. (d) Recorded activity (left) and model predictions (right) of 3 representative neurons with varying B coefficients and similar values for the A coefficient as well as the coefficient of determination (R²). (e) Individual cells (grey dots) and means for each cluster. Left: Bias (B) towards encoding either state or Δstate/Δt value; middle: Model amplitude (A); right: R² values as a measure of the quality of fit. Cluster 1: n=100, Cluster 2: n=33, Cluster 3: n=35, Cluster 4: n=18 cells. (f) Anatomical location of DA neurons with different color shades denoting ‘B’ values; marker size refers to ‘A’ values. (g) The model fit significantly lower amplitude parameters to mVTA neurons. Significance was calculated by means of two-sided t-test **P=0.0095. (h) mVTA DA neurons were significantly more likely to encode the state compared to lVTA DA neurons. mVTA DA neurons include cells that are i) biased in either direction, or ii) show mixed responses. Significance was calculated by means of two-sided Welch’s t-test; ***P=0.0008. (i) Same as in (h), but individual data points were linearly weighted for amplitude (i.e., their contribution to the mean, SEM and statistical test scaled with the A parameter, which prevents cells with a small amplitude from confounding the analysis). Significance was calculated by means of (weighted) two-sided Welch’s t-test **P=0.0098. All data represented as mean ± SEM (shading or error bars).