Multidimensional encoding of movement and contextual variables by rat globus pallidus neurons during a novel environment exposure task.

The basal ganglia (BG) play a critical role in a variety of functions that are essential for animal survival. Information from different cortical areas propagates through the BG in anatomically segregated circuits along the parallel direct and indirect pathways. We examined how the globus pallidus (GP), a nucleus within the indirect pathway, encodes input from the motor and cognitive domains. We chronically recorded and analyzed neuronal activity in the GP of male rats engaged in a novel environment exposure task. GP neurons displayed multidimensional responses to movement and contextual information. A model predicting single unit activity required many task-related behavioral variables, thus confirming the multidimensionality of GP neurons. In addition, populations of GP neurons, but not single units, reliably encoded the animals' locomotion speed and the environmental novelty. We posit that the GP independently processes information from different domains, effectively compresses it and collectively conveys it to successive nuclei.

Introduction

The basal ganglia (BG) are a group of subcortical nuclei that are known to be involved in a variety of functions spanning the limbic, associative and sensorimotor domains. Information from different areas in the neocortex reaches the BG’s primary input structure, the striatum, where it is organized topographically and is processed by segregated circuits with minor overlap (Peters et al., 2021; Hintiryan et al., 2016; Alexander et al., 1986). From the striatum, the information is conveyed in the parallel direct and indirect pathways. Both pathways project to the output structures of the BG; however, information in the indirect pathway undergoes additional processing in the globus pallidus (GP) which projects back to the striatum, and to the subthalamic nucleus (STN) (but see (Kupchik et al., 2015)). The STN also provides additional input to the BG via the hyperdirect pathway (Nambu et al., 2000, 2002; Mink, 1996; Kita and Kita, 2012). This complex projection network and specific connectivity enable the GP to influence both the input to, and the output from, the BG (Ketzef and Silberberg, 2021; Mastro et al., 2014).

It has been shown anatomically that the segregated pattern of information flow described in the striatum is maintained in the GP (Foster et al., 2021) even though the number of striatal efferents is an order of magnitude larger than that of the GP (Oorschot, 1996). By contrast, single GP neurons have been shown to respond to a variety of motor and non-motor task attributes possibly because these neurons have a very large dendritic arbor which can branch across sub-regions and therefore integrate information from several functionally distinct areas (Sadek et al., 2007; Kita and Kitai, 1994). For example, these neurons are known to encode movement direction including the signaling of movement onset, the end of a movement in a movement sequence, and are affected by static and/or dynamic load (Dodson et al., 2015; Gu et al., 2020; Brotchie et al., 1991a, 1991b; Mitchell et al., 1987; Arkadir et al., 2004). They also encode cognitive and limbic information such as the contextual setting, the relative difficulty of the task, strategies for behavioral inhibition, reversal learning and reward probability and prediction (Saga et al., 2013; Schechtman et al., 2016; Joshua et al., 2009; Gu et al., 2020; Brotchie et al., 1991b; Arkadir et al., 2004; Bischoff-Grethe et al., 2015; Lilascharoen et al., 2021). However, it remains unclear whether the GP processes all the features encoded in the striatum, or whether it extracts a few relevant features, and whether, like the striatum, these features are encoded by segregated pathways, or alternatively, are integrated by single GP neurons.

Based on their firing rate and pattern, GP neurons are typically divided into two types: high frequency pausers (HFP) and low frequency bursters (LFB) (Hutchison et al., 1994; DeLong, 1971; Elias et al., 2007; Benhamou et al., 2012). Recent studies have proposed a cell type classification of GP neurons based on their molecular signature, cellular properties and innervation targets, into arkypallidal neurons that project to the striatum and prototypic neurons that project to the STN, BG output nuclei and the thalamus (Ketzef and Silberberg, 2021; Mastro et al., 2014; Mallet et al., 2012; Abdi et al., 2015). Although not confirmed, it is thought that the LFBs are the arkypallidal neurons and that the HFPs are the prototypic neurons (Abdi et al., 2015).

When exposed to an unfamiliar environment, animals typically initiate exploratory behavior; i.e., they engage in rapid movement throughout the environment, typically starting with its perimeters and then extending into the center. This behavior allows animals to gradually familiarize themselves with the environment. Two distinct processes are essential for successful task performance: (1) a cognitive distinction between familiar and unfamiliar environments, and (2) the planning and execution of the physical locomotion during the familiarization process. These components are dissociable because the internally generated locomotion undergoes habituation as the animal familiarizes itself with the environment, whereas the identity of the environment remains unchanged as long as the animal stays in a given space. Hence, the novel environment exposure (NEE) task is well suited for studying how the GP integrates and extracts information encompassing different domains.

Previously, we showed that non-overlapping populations of striatal projection neurons – the medium spiny neurons – reliably encode the locomotion and environment identity during performance of the NEE task (Yamin et al., 2013). Here, we inquired whether during a similar task, GP neurons would encode either the animal’s actions or the environmental identity, or both. We chronically recorded and analyzed the activity from rat GP neurons with respect to locomotion, rearing up on the hind limbs, grooming and the different environments while the animals performed the NEE task. We show that populations of GP neurons reliably encode the animals’ average locomotion and the environmental novelty. In contrast, single neurons encoded both attributes poorly even though the majority of the GP neurons (90%) responded significantly to at least one task variable. The likelihood of a GP neuron to respond to a given task-related variable was independent of whether that neuron responded to other variables. Concomitantly, using generalized linear model, we found that many behavioral variables contributed to the prediction of single neurons’ firing rate traces during the task. Interestingly, more behavioral variables contributed to the prediction of LFB firing than to HFP firing, and the LFBs were better tuned to the animals’ speed than the HFPs. These findings suggest that GP neurons integrate different types of task-relevant information, and distribute it independently across populations of single neurons.

Results

Rats explore the novel environment by locomotion and rearing up on their hind limbs

Seven male Long-Evans rats engaged in the NEE task (see STAR Methods) once a week for up to four consecutive weeks (n = 27 sessions). Each experimental session comprised 15 min in the home cage, defined as baseline activity (BL), after which the animal was transferred to the novel environment (NE), a cage identical to the home cage but lacking familiar odors or objects. After a period of 15 min in the novel environment, the animal was transferred back to its home cage (HC) for an additional 15 min period (Figure 1A). Upon transition to the NE, the rats initiated exploratory movement that consisted primarily of locomotion and rearing up on their hind limbs and leaning against the cage wall (see example session in Figure 1B). The exploration of the NE gradually ceased and the animals tended to select a preferred corner and remained motionless while sporadically initiating short epochs of exploratory behavior. When transferred back to their home cages the animals reinitiated exploratory behavior; however, the exploration time of the HC was significantly shorter than for the NE (326 ± 53 s and 475 ± 66 s for HC and NE, respectively, paired t-test, t = 2.37, df = 26, p = 0.026, see STAR Methods). Along with the increased exploration time, the fraction of time the animals spent walking and rearing up was significantly higher in the NE as compared to the BL and the HC (Figure 1C, locomotion: 5.5 ± 0.8%, 9.9 ± 1.0% and 7.1 ± 0.8% for BL, NE and HC, respectively; ANOVA, F = 7.06, df = 2, p = 0.0015; Figure 1D, rearing up: 3.9 ± 0.8%, 16.2 ± 2.4% and 7.1 ± 1.6% for BL, NE and HC, respectively; ANOVA, F = 13.08, df = 2, p = 1.8 × 10−5, post-hoc testing with multiple comparison adjustments). These findings suggest that rats use a strategy of interleaved locomotion and rearing up on their hind limbs to explore the NE.

Figure 1

Rats explore the novel environment by locomotion and rearing up on their hind limbs

(A) Experimental timeline: rats spent 15 min in their home cage, defined as the baseline activity, then were transferred to a novel environment for 15 min, and finally transferred back to their home cage for another 15 min.

(B) A single animal’s trajectory (gray lines) and rearing up locations (circles) during BL (left), NE (center) and HC (right) of an example session. Rectangles represent the cage (35 × 46 cm2).

(C) Boxplots of the fraction of time (in %) the animals spent in locomotion in each environment. In each box, the central mark indicates the median, and the bottom and top edges of the box indicate the 25th and 75th percentiles, respectively. The whiskers extend ±2.7 SD. Colored dots denote individual sessions, and each color represents values obtained from a single animal.

(D) Boxplots of the fraction of time (in %) the animals spent rearing up in each environment. Boxplot conventions as in (C).

(E) Animals’ speed in 5 s bins averaged across sessions from the same exposure week (1–4). Gray scale denotes the exposure number. Black vertical lines represent transition times between experimental conditions.

(F) Animals’ speed in 5 s bins averaged across all sessions (black) and SEM (gray shading). Black vertical lines represent transition times between experimental conditions.

(G) Boxplots of the animals’ speed for each experimental condition. Boxplot conventions as in (C).

In all subplots, the star denotes an ANOVA, p < 0.01 with post-hoc testing adjusted for multiple comparisons.

As reported in mice (Yamin et al., 2013), the rats did not habituate to the task and displayed a similar pattern of locomotion speed in all sessions irrespective of the number of exposures (Figure 1E). Hence, locomotion speed was averaged across sessions (Figure 1F). Movement speed was significantly higher in the NE and the HC as compared to BL (Figure 1G; 10.4 ± 0.3 cm/s, 12.3 ± 0.2 cm/s and 11.6 ± 0.3 cm/s for BL, NE and HC, respectively; ANOVA, F = 10.56, df = 2, p = 9 × 10−5, post-hoc testing corrected for multiple comparisons). The pattern of the average locomotion speed was similar to that previously reported in mice, suggesting that epochs of rearing up occurred sporadically, and thus, were canceled out.

GP neurons can be classified into high frequency pausers and low frequency bursters

We chronically recorded the activity of 78 single neurons in the GP of 7 rats during performance of the NEE task (n = 21 sessions). The electrodes’ position in the GP was verified histologically for all recorded neurons (Figures 2A and 2B; see STAR Methods). As in previous work, we tested whether the recorded neurons could be classified into two distinct subpopulations based on their electrophysiological properties during BL (Benhamou et al., 2012). First, we calculated three firing features that have been shown to differentiate HFPs from LFBs: the Fano Factor (FF), the post-spike suppression (PSP) and the mode of the inter-spike interval (mode ISI), and plotted them on a 3D graph (Figure 2C; see STAR Methods). This presentation revealed two distinct clusters, each corresponding to a different cell type as identified by its properties: the HFP cluster contained 50 neurons (marked in blue in Figure 2C) and the LFB cluster contained 28 neurons (marked in red in Figure 2C). Figures 2D–2G depict raw traces and their corresponding autocorrelation of representative HFP (Figures 2D and 2E) and LFB (Figures 2F and 2G). Comparisons of other firing properties across these two clusters revealed that all of them except for firing rate were significantly different, thus strengthening the classification of GP neurons into HFPs and LFBs (see Table 1).

Figure 2

GP neurons can be classified into high frequency pausers and low frequency bursters

(A) A 40 μm coronal slice stained by AChE (see STAR Methods). The arrows denote the electrolytic lesion indicating that the recording electrodes were placed in the GP.

(B) Electrode placement of all neurons included in the dataset overlaid on a planar section. Electrodes at different depths are projected on the displayed section. Black and green dots indicate electrodes placed in the right and left hemispheres, respectively.

(C) A 3D scatterplot of GP neurons’ firing properties generated two distinct clusters each corresponding to a different cell type: the HFPs (n = 50; marked in blue) and the LFBs (n = 28; marked in red). Arrows denote the example neurons shown in panels (D–G).

(D) A representative raw trace of a neuron identified as an HFP. Tick marks above the trace represent sorted action potentials.

(E) An autocorrelation of the spike train of the neuron shown in (D).

(F) Same as D for a neuron identified as an LFB.

(G) Same as (E) for the neuron shown in (F).

Table 1

Firing properties of neurons classified as LFBs and HFPs

Property	HFP (n = 50)Median (min – max)	LFB (n = 28)Median (min – max)
Baseline firing rate (spikes/s)	4.47 (2.8–30.90)	3.76 (1.1–20.78)
Fano Factor (FF)	0.71 (0.31–1.01)	1.31 (0.98–2.25) ∗
Coefficient of variation (CV)	0.67 (0.36–1.15)	1.38 (0.87–2.62) ∗
Post spike suppression (PSP)	57.5 (14–128)	8.5 (0–24) ∗
Mode ISI distribution	0.08 (0.02–0.29)	0.02 (0.002–0.05) ∗
Percent spikes in bursts (%)	5.14 (0–42.81)	36.63 (2.52–64.03) ∗

∗ denotes significant difference between groups (Kruskall-Wallis test, p < 0.01, Bonferroni corrected).

Single GP neurons display multidimensional responses during the NEE task

Throughout the task, the neurons displayed tonic firing with complex modulations (see two example neurons in Figure 3A (for LFB) and 3B (for HFP); see Video S1). Examination of these firing rate modulations with respect to different animal behaviors (locomotion, rearing up and grooming, denoted in the figure by different shades of green) indicated that the neuron shown in Figure 3A increased its firing rate for locomotion and rearing up and decreased its firing rate during grooming as compared to its firing rate during quiet awake (QA) periods (Figure 3C; see STAR Methods). The neuron shown in Figure 3B increased its firing rate during grooming and rearing up but not during locomotion as compared to the QA periods (Figure 3D). Examination of firing rates measured in these example neurons during the QA periods in the different environments showed that the neuron in Figure 3A displayed a higher firing rate in the NE and the HC as compared to BL, whereas the neuron shown in Figure 3B did not distinguish between environments (Figures 3E and 3F, respectively). Overall, the fraction of HFPs that responded to at least one action (39/50 (78%)) was similar to that of the LFBs (25/28 (89%), χ2 test, χ2 = 1.6, df = 1, p = 0.21). Similarly, the fraction of HFPs that encoded the environment (30/50 (60%)) was similar to that of the LFBs (13/28 (46%), χ2 test, χ2 = 1.3, df = 1, p = 0.25). Hence, we avoided decreasing statistical power and performed the analysis described below on the whole population.

Figure 3

Single GP neurons display multidimensional responses during the NEE task

(A) An example firing rate trace of a neuron identified as an LFB (black) throughout the NEE task, calculated in 0.5 s bins and smoothed using a Gaussian window of 40 bins. Background colors represent the activity of the animal and are denoted in the figure legend. Gray vertical lines represent transition times between experimental conditions. See also Video S1.

(B) Same as in (A) for a neuron identified as an HFP.

(C) The average firing rate of the neuron plotted in (A) while the animal was engaged in the different actions. Error bars are SEMs. Color scheme as in (A).

(D) Same as (C) for the neuron plotted in (B).

(E) The average firing rate during quiet awake periods for the different experimental conditions. Same neuron as in (A). Error bars are SEMs.

(F) Same as (E) for the neuron plotted in (B).

(G) Pie chart of the proportion of neurons encoding different combinations of tested actions.

(H) The expected (dotted black line) and the observed (purple line) percentages of units encoding different possible combinations of actions. The expected value for these combinations was calculated from the marginal probabilities of the tested actions assuming independent processes. Gray rectangles represent two-way 99% confidence intervals calculated from 10,000 repetitions of random draws of responses from the marginal probabilities of each action, and were plotted for demonstration sake only.

(I) The expected (dotted black line) and the observed (purple line) percentages of units encoding both the environmental novelty and any of the analyzed actions. The expected values were calculated from the marginal probabilities to encode an action and the NE assuming independent processes. Gray rectangles as in (H) but for the NE and the different actions.

In all subplots the star denotes ANOVA, p < 0.01 with post-hoc testing adjusted for multiple comparisons.

Video S1. Activity of an example neuron during different rat behaviors

Same neuron as in Figure 3A, related to Figure 3.

Quantification of firing rate alterations in the whole population during the different motor behaviors showed that the majority of GP neurons displayed significant firing rate alterations in response either to one (21.8%), two (37.2%) or three behaviors (23.1%) whereas only a minority (17.9%; Figure 3G) of the neurons did not respond to any of these motor behaviors. To examine whether the likelihood of a GP neuron to encode one motor behavior co-varied with its likelihood to encode another behavior, we calculated the probability of GP neurons to encode different combinations of behaviors based on their marginal probabilities to encode single behaviors. The use of the marginal probabilities sufficed to predict the probabilities of GP neurons to encode a single behavior, two behaviors or three behaviors (Figure 3H; 2 test for independence, 2 = 3.4, df = 4, p = 0.49; see Table S1 for data per cell type). This indicates that the likelihood of a GP neuron to encode a given motor behavior was independent of whether or not it encoded other behaviors. That is, even though locomotion and rearing up are part of the exploration process and rearing up and grooming have a common body posture, their joint probability was similar to that expected from independent processes.

Quantification of firing rate alterations of GP neurons during the QA periods in the different environments (experimental conditions; see STAR Methods) showed that the majority (55%) of the neurons encoded information about the environment, out of which 46.5% (20/43) distinguished the NE from the other two environments. Only 8/78 (10.3%) of the neurons did not respond to any task attribute. We then examined whether the likelihood of a GP neuron to encode any of the tested motor behaviors co-varied with its likelihood to encode the NE. To that end, we calculated the likelihood of GP neurons to encode the NE and/or any motor behavior by using the marginal probabilities of the different behaviors and the NE. The marginal probabilities sufficed to predict the probability of the GP neurons to encode any behavior and/or the NE (Figure 3I; 2 test for independence, 2 = 0.15, df = 1, p = 0.70, 2 = 0.23, df = 1, p = 0.63, 2 = 0.87, df = 1, p = 0.35, for grooming, locomotion and rearing up, respectively). Thus, single GP neurons appear to be able to encode multiple task-related variables, and the propensity to encode one variable is independent of the neuron’s propensity to encode other variables.

Many behavioral predictors contribute to the activity of single GP neurons

To directly assess the contribution of the different behavioral variables to the activity of single GP neurons, we examined whether the firing rate profile of each neuron could be predicted by a generalized linear model (GLM) which took into account a variety of task-related variables (Engelhard et al., 2019). The model included continuous and categorical variables (see STAR Methods). The continuous variables were the speed and the acceleration of the animal’s center of mass. The categorical variables consisted of the type of behavior and the experimental conditions. The behavioral variables included bins tagged according to whether the animal was engaged in grooming, rearing up, locomotion or QA. The remaining bins were marked as ‘others’. The experimental condition variables were composed of bins tagged according to the environment in which the animal was located (i.e. BL, NE or HC). The variables were fed into the GLM to predict the firing of single GP neurons recorded during the NEE task (Figure 4A). We used R2 as a metric to assess the model’s accuracy. The significance of each model was evaluated by testing whether its accuracy exceeded the confidence interval of an accuracy distribution calculated from control models in which the neuron’s firing was predicted from the same behavioral data that were shifted forward in time (see STAR Methods). The GLM of 40/78 (51.3%) of GP neurons were significant, of which 17 (61%) were LFBs and 23 (46%) were HFPs (χ2 test, χ2 = 1.55, df = 1, p = 0.21). Two examples of the estimated firing traces (blue and red for HFP and LFB, respectively) overlaid on the actual firing traces of their corresponding GP neurons (black) are shown in Figures 4B and 4C. Figure 4D depicts the distributions of the R2 values of the significant (blue and red for HFP and LFB, respectively; range: 0.03–0.5) and non-significant (gray; range: 0.001–0.34) models. The R2 values were similar across the two subpopulations (t-test performed on the Fisher transform of the correlation values, t = 1.62, df = 76, p = 0.11). Adding interaction-between-variables terms to the GLM did not substantially increase the number of significant models (40 when using linear terms and 42 when including interaction terms) and did not improve the quality of approximation as measured by the R2 values of the significant models (a range of 0.03–0.5 in both cases). Thus, inclusion of interaction terms in the GLMs increased the complexity of the model (degrees of freedom) without improving the outcome.

Figure 4

Many behavioral predictors contribute to the activity of single GP neurons

(A) Schematic illustration of the behavioral and environmental variables used by the GLM to predict single neuron’s firing rate.

(B) An example of a GLM fit (blue) to the actual firing rate trace of a neuron identified as an HFP (black).

(C) An example of a GLM fit (red) to the actual firing rate trace of a neuron identified as an LFB (black).

(D) The distribution of R2 values for all GLMs. Blue and red denote the significant models of the HFPs and the LFBs, respectively. Gray denotes the non-significant models.

(E) The number of predictors in the true GLM (x axis) versus the mean ± SEM number of predictors in a set of shifted GLMs (y axis; see STAR Methods). Only significant models are plotted. HFPs and LFBs are marked in blue and red, respectively.

The number of contributing predictors to the significant GLMs was high (6.32 ± 0.28) considering that each model could potentially have a maximum of 9 coefficients. In the vast majority (33/40) of the significant models, the number of contributing predictors of the neurons’ firing rates was higher than that of the shifted models (Figure 4E; 4.81 ± 0.07 predictors in the shifted models, KW test, H = 24.62, df = 1, p = 6.9 × 10−7). The number of contributing predictors of the LFBs’ firing rates was higher than that of the HFPs’ firing rates (Figure 4E, blue and red markers for HFP and LFB, respectively, 7.47 ± 0.28 and 4.46 ± 0.22, t-test, t = 4.2, df = 38, p = 0.0001). Observation of the significant models’ coefficients revealed high variability in the relative contribution of the different predictors (Figure S1). Interestingly, the speed variable contributed to a significantly larger number of LFB models than to the HFP models (12/17 and 3/23 of the significant models for LFB and HFP, respectively, χ2 test, χ2 = 13.81, df = 1, p = 0.0002). Overall, many task-relevant variables appeared to contribute to the firing of each GP neuron. This outcome is consistent with the result described in the previous section that GP neurons modulate their firing rate significantly in response to behavioral and environmental conditions.

The animals’ locomotion and the environmental novelty are encoded by populations of GP neurons but not by single neurons

The fact that many task-relevant variables were required for successful modeling of single neuron activity implies that input about these variables may be distributed across the population rather than be encoded by single neurons. To investigate this possibility, we used a linear regression model to approximate two task attributes that were common to all animals: the average locomotion speed and the environmental novelty (see STAR Methods). The population of GP neurons yielded high R2 values for the animals’ average speed (R2 = 0.77, Figure 5A) and the environmental novelty (R2 = 0.70, Figure 5B). To control for the degrees of freedom, we randomly shuffled the firing rate of the neurons 100 times and recalculated the models’ accuracy (see STAR Methods). The probability to achieve the population’s R2 values or higher using the shuffled data was 0 for both the locomotion speed and the environmental novelty (Figures 5C and 5D, respectively) suggesting that the high performance of the regression model was not due to the number of predictors (n = 78 neurons).

Figure 5

The animals’ locomotion and the environmental novelty are encoded by populations of GP neurons but not by single GP neurons

(A) The animals’ speed in 5 s bins averaged across all sessions (black), and the estimated speed (purple) calculated by a linear regression model using the firing rate of all recorded neurons (n = 78).

(B) Same as (A) for NE identity.

(C) The distribution of R2 values of linear regression models applied to shuffled firing rates to predict the animals’ average speed (see STAR Methods). Purple line is the R2 value of the true model.

(D) Same as (C) for estimating the NE identity.

(E) The best example of a linear regression model applied to a single GP neuron to estimate the animal’s speed. The estimated speed (orange) overlaid on the animal’s speed in 5 s bins (black).

(F) Same as (E) for estimating the NE identity.

(G) The distribution of R2 values of linear regression models applied to single units to predict the animal’s speed in a session. Purple line is the R2 value of the model using all neurons (n = 78). The example shown in E is marked by an orange arrowhead. Inset shows the R2 value distributions of HFPs (blue bars) and LFBs (red bars).

(H) Same as (G) for estimating the NE identity. The example shown in (F) is marked by an orange arrowhead.

(I) Boxplots of R2 values calculated from 100 repetitions of the linear regression model for animals’ speed using a growing number of neurons drawn from either the HFPs (blue) or the LFBs (red) (see STAR Methods). In each box, the central mark indicates the median, and the bottom and top edges of the box indicate the 25th and 75th percentiles, respectively. The whiskers extend ±2.7 SD and outliers are marked by a + sign. (J) Same as (I) for the novel environment.

The high performance of the regression model on the population could have emerged either from a few single neurons that effectively predicted the two task attributes, or alternatively, that the minor contribution of many single neurons was summed. To distinguish between these possibilities, we applied linear regression models to single GP neurons to assess how well they approximated the animal’s speed and presence in the novel environment. Unlike the population, single neurons approximated both task attributes poorly (the best example models are shown in Figures 5E and 5F, for speed and environmental novelty, respectively). Figures 5G and 5H depict the R2 distribution of single units (white bars) as compared to the population (purple bar). Interestingly, the R2 values of single LFBs were higher than those of single HFPs when approximating the animal’s speed (see inset in Figure 5G; t-test performed on the Fisher transform of the correlation values, t = 2.77, df = 76, p = 0.007) but not when approximating the environmental novelty (inset in Figure 5H; t-test performed on the Fisher transform of the correlation values, t = 1.59, df = 76, p = 0.1).

These results indicate that each neuron by itself made a minor contribution to the population model and support the claim that individual neurons may convey non-overlapping information that tends to accumulate when integrated across neurons. We tested this possibility by recalculating the regression model on an increasing number of randomly selected groups of LFBs and HFPs (see STAR Methods). The averaged R2 gradually increased with the number of units, suggesting that the information was distributed across the whole population (Figures 5I and 5J). Overall, populations of LFBs outperformed the populations of HFPs not only when approximating the animals’ speed but also when approximating the environmental novelty, as quantified by measuring the effect size using Cohen’s d (see STAR Methods). For the NE estimation, LFB had a big effect size as compared to the HFP starting from 10 units, and for all unit numbers for the estimation of speed (see Table 2 for Cohen’s d values). These results indicate that information about many task variables is distributed independently across GP neurons, and may be restored and conveyed to successive nuclei by integrating information from the population activity.

Table 2

Effect size of R2 between LFB and HFP as measured by Cohen’s d

Number of units	5	10	15	20	25
Speed LFB - Speed HFP	1.38	2.29	3.45	4.50	6.47
NE LFB - NE HFP	0.53	1.03	0.90	0.82	0.97

Cohen’s d values between the R2 distributions calculated for models with different number of units. Values higher than 0.8 are considered big effect (see STAR Methods).

The gain in regression performance when taking into account the population as compared to single units could have emerged from averaging across correlated signals with random noise, thus enhancing the signal to noise ratio. To control for this possibility, we calculated the cross-correlograms (see STAR Methods) of all simultaneously recorded pairs of GP neurons. We found that a minority of the tested pairs (12/147 (8%)) were significantly correlated thus ruling out this possibility.

Discussion

This study explored the encoding of diverse attributes of the NEE paradigm by GP neurons. We showed that rats explored the novel environment by locomotion and rearing up on their hind limbs. We classified single GP neurons into two cell types and showed that irrespective of their type, GP neurons responded during behaviors performed by the rats such as grooming, rearing up and locomotion and distinguished between environmental conditions including identifying the novel environment. The probability of a GP neuron to respond to one or more of these task attributes was independent of whether the neuron responded to other attributes, and the activity of simultaneously recorded pairs of neurons was uncorrelated. Prediction of the single neuron firing rate by GLM demonstrated that many behavioral variables contributed to each model. Nevertheless, the information carried by single GP neurons was insufficient to approximate either the animal’s locomotion or the animal’s presence in the novel environment. By contrast, both task attributes were reliably approximated by a population of GP neurons. Consistent with the fact that the speed variable contributed to the GLM of more LFBs than HFPs, populations of LFBs were better than populations of HFPs in approximating the animals’ average speed. Populations of LFBs outperformed the HFPs as well when approximating the environmental identity. These data indicate that task-relevant information is distributed across GP neurons in a non-overlapping manner, and can therefore be extracted from the population, rather than from single GP neurons.

It should be noted that although GP neurons may encode variables such as anxiety and sensory novelty (Lever et al., 2006), we analyzed only overt, quantifiable behavioral variables and did not address variables that may have triggered the different behaviors. The fact that GP neurons responded to many behavioral variables raises the possibility of a common latent variable. While we cannot reject this possibility, its likelihood declines by the finding that the probability of GP neurons to encode different behavioral variables could be calculated from their marginal probabilities assuming independence. A latent variable would have increased the probability of each neuron to encode specific behaviors beyond chance and by that reject the null hypothesis of independence. In addition, a latent variable would have decreased the number of contributing predictors to the GLM due to increased dependence.

It has been shown that single neurons in the globus pallidus external segment (GPe) of primates, which is homologous to the rodent GP, can encode more than one task-relevant variable. For example, GPe neurons were reported to encode the direction of reaching movement and reward during a probabilistic visuomotor task (Arkadir et al., 2004), and could encode effort and reward in a task where a different force was required to obtain different amounts of reward (Nougaret and Ravel, 2018). These task-relevant parameters were integrated linearly, suggesting their independent encoding by GP neurons (Arkadir et al., 2004). Moreover, GP neurons were reported to generate dynamic ensembles which at any given time encoded the most relevant information for best task performance (Nougaret and Ravel, 2018; Adler et al., 2012; Saga et al., 2013). Our results in the rodent GP are in agreement with those described in the primate GPe; and extend them by showing that two distinct cell types in the GP, the HFPs and the LFBs, encode multiple actions as well as the nature of the environment simultaneously and independently, thus encompassing both the motor and cognitive domains.

During the NEE task, striatal MSNs encoded either the animal’s locomotion or the environmental novelty but not both (Yamin et al., 2013). The identity of MSN encoding each task attribute was not addressed. Hence, it remains unclear whether these task attributes were encoded selectively by MSNs from either the direct or the indirect pathway, or both. The fact that GP neurons responded to both task attributes indicates that MSNs from the indirect pathway contributed to task performance. Furthermore, it rules out the possibility that the GP processed these task attributes in segregated parallel pathways. It remains unclear whether MSNs from the direct pathway also process the two task attributes, and whether the GP selects specific task-relevant features to maximize performance or optimally compresses all the converging data from the striatum using dimensionality reduction (Bar-Gad et al., 2000, 2003b). The NEE paradigm may be too simplistic to differentiate between these two alternatives because both the motor and the cognitive features are essential for task performance. Hence, deciphering whether GP neurons are tuned to extract functionally relevant features or to compress all converging data from the striatum may require tasks utilizing more complex feature representations.

We classified the GP neurons into HFPs and LFBs based on electrophysiological properties measured during the baseline period. We found task-related differences between these cell types. First, more behavioral predictors were required to approximate the firing rate of single LFBs than single HFPs. This outcome may indicate a difference in the role played by the two subpopulations, and may reflect the fact that the two cell types utilize different encoding schemes to convey information to their target nuclei; i.e., bursts rather than pauses. Second, the LFBs contained more information about the animals’ speed than the HFPs. This outcome is consistent with a previous report showing that the motor cortex directly drives the arkypallidal neurons (Karube et al., 2019), thus supporting the claim that the LFBs are the arkypallidal neurons and that the HFPs are the prototypic neurons, as previously suggested by Abdi et al. (2015).

The firing rate of GP neurons recorded in this study was lower than that described in other studies conducted in vivo, including a study from our own laboratory (Abdi et al., 2015; Benhamou et al., 2012; Gu et al., 2020). Although the source of this difference remains unclear, one potential explanation has to do with the use of recording electrodes with different diameters and coating thicknesses, which can influence the electrodes’ impedance. The electrodes’ impedance influences the noise level and signal attenuation and thus may have created a sampling bias (Viswam et al., 2019). Nonetheless, the parameters used in this study for classification of the GP neurons into HFPs and LFBs were similar to those described elsewhere (Benhamou et al., 2012). In both studies, the neurons’ firing rates did not differ significantly across cell types, and thus did not contribute to their classification. This outcome and the fact that firing patterns rather than firing rates are preserved across species (Benhamou et al., 2012) highlight the importance of the unique encoding schemes utilized by each subpopulation.

A recent study of BG networks showed that the mouse GP faithfully follows the striatal topography and is characterized by a small percentage of converging inputs arising from different regions of the striatum (Foster et al., 2021). This regional specificity, which is more segregated than the coarse division into limbic, associative and sensorimotor domains, highlights the expectation of a limited overlap in converging afferents onto the GP. Hence, it remains unclear how GP neurons display multidimensional responses to task attributes in both primates and rodents. GP neurons have a large dendritic arbor extending across large portions of the GP volume that is aligned perpendicularly to MSN axons (Kita and Kitai, 1994; Park et al., 1982; Jaeger and Kita, 2011; Percheron et al., 1984). This anatomical organization may allow GP neurons to integrate information from relatively distant striatal domains but at the same time yield a substantial overlap of inputs in neighboring GP neurons. The high overlap in converging inputs is likely to generate interactions between the probabilities of GP neurons to encode different task attributes and enhance neuronal cross-correlations. However, as shown here and elsewhere (Bar-Gad et al., 2003a), GP neurons display weak to nonexistent correlations and their probabilities to respond to a given action is independent of whether they responded to another action. Two distinct scenarios may account for these phenomena: the first utilizes the high convergence ratio of 40:1 between MSNs to GP neurons (Oorschot, 1996). In this case, MSN projections onto GP neurons are highly selective so that despite the high spatial overlap, each GP neuron receives afferents from a distinct group of MSNs. The second scenario suggests that although overlap in GP afferents exists, their common input is canceled out by other mechanisms. For example, by the recurrent collateral axons of GP neurons (Kita et al., 2004) which have a complex spatial spread in the GP (Matsumura et al., 1995) and/or by the strong reciprocal connections between the GP and the subthalamic nucleus (STN) (Atherton et al., 2013; Baufreton et al., 2009; Magill et al., 2000).

The encoding scheme utilized by GP neurons is substantially different from the one utilized by the striatal projection neurons, the MSNs: GP neurons display multidimensional responses whereas the MSNs respond typically to a single event by a short bursting activity (Fobbs et al., 2020; Berke et al., 2009; Shidara et al., 1998; Atallah et al., 2014; Jog et al., 1999; Jin and Costa, 2010). Further research is needed to determine why the GP converts the temporally precise, highly specific information present in single MSNs into spatially distributed information where each component carries little information about the task. Additional experiments are required to determine how information from the striatum and the GP is integrated in the basal ganglia output structures, the entopeduncular nucleus and the substantia nigra pars reticulata, and which task attributes are represented in these structures.

Limitations of the study

Two potential caveats should be taken into account when interpreting the results reported in this study: the relatively small sample size of GP neurons and their lower-than-expected firing rate. We recorded the activity of 78 single GP units that were subdivided into HFPs (n = 50) and LFBs (n = 28). Although we used appropriate statistical methods to evaluate the results’ significance, we cannot rule out the possibility that results reported here as insignificant may become significant when tested on a larger sample size.

The second caveat is that in the current study the recorded GP neurons displayed firing rates lower than those reported in previous studies. We have verified carefully the quality of recorded neurons by repeating offline sorting and electrode placement analysis by two independent individuals. However, it remains plausible that our results characterize a subpopulation of GP neurons having a low firing rate, and thus, caution is in place when generalizing these results to the high firing neurons of the GP.

STAR★Methods

Key resources table

Resource availability

Lead contact

Further information and requests for resources and reagents should be directed to and will be fulfilled by the lead contact, Dana Cohen (dana.cohen@biu.ac.il).

Materials availability

This study did not generate new unique reagents.

Experimental model and subject details

Animals

All procedures were approved by the Bar Ilan University Institutional Animal Care and Use Committee and were performed in accordance with the National Institutes of Health guidelines. Data were collected from 7 Long Evans male rats (3–7 months old). All animals were housed two in a cage, separated by a divider after surgery, on a 12/12 h light/dark cycle and had ad-libitum access to food and water. Experiments were performed during the light phase.

Method details

Surgery

The surgical procedure has been described in detail elsewhere (Benhamou et al., 2012). In brief, 7 adult male Long Evans rats weighing 460–545 gr (bred in-house) were initially sedated with 5% isoflurane and then injected intramuscularly with ketamine HCl and xylazine HCl (100 and 10 mg/kg, respectively). Supplementary injections of ketamine and xylazine were given as necessary. The skull surface was exposed and two 0.5 × 2 mm2 craniotomies were made above the GP. Each craniotomy was shaped as a rectangle centered on (anterior-posterior: 1.65 mm; mediolateral: ± 4.1 mm; specific coordinates were: [-0.7, ±3.7], [-1.0, ±3.3], [-2.6, ±4.5], [-2.3, ±4.9] mm relative to bregma; (Paxinos and Charles, 2014)). Custom made microwire electrodes (microwires of S-isonel-coated tungsten, California Fine Wire Company, 50 μm in diameter arranged in a 2 × 8 array) were lowered from the surface of the brain while recording neural activity. The final placement of the electrodes was determined based on the measured coordinates from the surface of the brain (5.5–6.6 mm) and the quality of neural activity. Electrode placement in the GP was confirmed histologically after electrolytic lesioning, perfusion with 10% formalin, brain fixation with 20% sucrose, sectioning and acetylcholinesterase (AChE) staining. A protocol for AChE staining was obtained from neurosciencecourses.com and adjusted for perfused tissue. The rats were given at least 1 week to recover prior to recording.

Novel environment exposure (NEE) task

24 h prior to the task, each rat was individually housed in a clean cage (35 × 46 × 20 cm3) with a small amount of bedding to familiarize the animal with its cage. Then the cage lid was opened, and the rat was connected to the head stage for at least 30 min before the beginning of the experiment to allow for habituation to the recording wires and the experimenter’s presence. Each recording session comprised 15 min of baseline activity in the home cage (BL), 15 min in the novel environment (NE) which was a clean cage identical to the home cage in size and amount of bedding but lacking familiar odors, and an additional 15 min back in the home cage (HC). Sessions were repeated once a week for up to 4 consecutive weeks. The animals’ behavior was monitored by a digital camera synchronized with the recording system (30 frames/s, Cineplex, Plexon Inc.).

Data acquisition

Neuronal activity was recorded throughout the performance of the NEE task. The activity was amplified, band-pass filtered at 0.5–8000 Hz and continuously sampled at 40 kHz using the OmniPlex data acquisition system (Plexon Inc). The broadband activity was band-pass filtered for spiking activity (300-8000Hz). Offline sorting (offline sorter, Plexon Inc) was performed on all recorded units having a signal-to-noise ratio of more than 3:1 independently by two individuals and only units identified as single cells by both were included in the final dataset. Single units were identified as distinct clusters in a 3D presentation of the first three principal components of the waveforms that remained stable throughout the recording session. Additional analysis was performed in MATLAB (R2013b, MathWorks Inc., Natick, MA).

Quantification and statistical analysis

The short transition periods between environments were removed from all subsequent analyses, thus the durations of all trials were equal across sessions.

Behavior analysis

The animals’ behavior was video recorded during the NEE task. The grooming and rearing up times were manually noted from the video. Rearing up time was measured from the moment the animal rose up on its hind limbs and raised its fore limbs into the air until the front paws touched the floor again. Grooming time was measure from the beginning of the stereotypical grooming movements until these movements ceased. Then, for all other frames, the coordinates of the animal’s center of mass were extracted using CinePlex Studio recording & tracking software (Plexon inc), and convoluted with a 9 frame Gaussian window. The speed of the animal in each frame was calculated by measuring the distance traveled in two consecutive frames and dividing it by the inverse of the frame rate. The animal’s speed was averaged in 0.5 s bins. Bins in which the animal’s speed was above 5 cm/s were considered locomotion bins. Bins in which the animal’s speed was lower than 0.5 cm/s were considered quiet awake (QA) periods. Bins with intermediate speeds were marked as ‘others’.

Exploratory behavior duration

The time the animals were engaged in exploratory behavior was measured from the time they were transferred to the NE/HC until identifying a time-period of 15 s without locomotion or rearing up on hind limbs.

Baseline firing rate

The number of spikes during baseline divided by its length in seconds.

Fano Factor (FF)

The variance of the spike count in non-overlapping 100 ms bins divided by its mean.

Mode of ISI distribution

The mode of the inter spike interval (ISI) distribution.

Post spike suppression (PSP)

The spike train autocorrelation of each unit was calculated during baseline using a ±1000 ms window at a resolution of 1 ms. PSP was defined as the earliest latency at which the autocorrelation rate equaled the baseline firing rate of the unit.

Percent spikes fired in bursts

Calculated in a similar manner to that described in detail in (Benhamou et al., 2012). Briefly, in order to identify bursts in the spike trains, we used the Poisson surprise method (Legendy and Salcman, 1985), which calculates how surprising it is to obtain n spikes in a time interval T given the firing rate r, where the time interval T is set to the mean firing rate divided by 2. Poisson surprise (S) was calculated by where . The burst percentage counts the number of spikes in the bursts compared to the total spikes emitted by the neuron.

Coefficient of variation (CV)

The SD of the inter-spike interval distribution divided by its mean.

Action encoding neurons

For each neuron, the firing rate distribution in 0.5 s non-overlapping bins for each behavior (QA, grooming, rearing up and locomotion) was calculated for the whole session. To determine whether a neuron significantly modulated its firing rate in response to any of these behaviors, the firing rate distributions of all behaviors were compared to its firing rate distribution calculated during QA using a two tailed, 1-way unbalanced ANOVA (p < 0.01).

Experimental condition encoding neurons

We calculated the firing rate distributions for each neuron in 0.5 s non-overlapping bins during the QA periods separately for each environmental condition (BL, NE and HC). To determine whether a neuron significantly changed its firing rate during any of the environments we compared its firing rate distributions during QA periods across environments using a two tailed, 1- way unbalanced ANOVA (p < 0.01).

Neurons showing a significant ANOVA (p < 0.01) after Bonferroni correction for multiple comparisons (n = 156 comparisons for behaviors and environments) were tested post-hoc to determine which action and/or environment they responded to. A neuron that responded differently during NE as compared to BL and HC was considered a novel environment encoding neuron.

Generalized linear model (GLM)

To quantify the contribution of the different variables to neural activity, we used a GLM. We used the MATLAB function stepwiseglm to fit the per trial neuronal firing rate in 5 s bins during the NEE paradigm. stepwiseglm begins with an initial constant model and takes forward or backward steps to add or remove variables, until a stopping criterion is satisfied. We used linear and interaction models and the Akaike information criterion to add or remove predictors and assumed a Poisson distribution of the firing rate. The predictors included continuous and categorical variables. The continuous variables quantified kinematics including speed and acceleration (calculated as the derivative of the speed), that were z-scored prior to modeling. The categorical variables were an action vector corresponding to the per bin main activity of the animal (QA, grooming, rearing up, locomotion, or others as defined in the behavioral analysis), and an experimental condition vector corresponding to the per bin experimental condition (BL, NE or HC). The use of the described input variables can yield a maximum of 9 and 29 coefficients, including the intercept for the linear and interaction models, respectively. The coefficient of determination, was used as a measure of the models’ accuracy. To evaluate the significance of each model, we calculated the of 25 shifted models, in which the predictors’ matrix was shifted forward in time by n bins of 5 s (n = 1, 2, …, 25 bins) relative to the firing rate which remained unchanged. The preservation of the temporal structure reduces the likelihood of identifying significant models erroneously (Harris, 2021). The model was considered significant if its exceeded the 95% confidence intervals of the distribution of the shifted model.

Linear regression analysis

Linear regression (task attribute = A + B⋅FR) was calculated between the firing rates of one individual neuron at a time or for a group of neurons (in 5 s bins) and a task attribute vector – either the speed of a single trial (for single neurons) or the average speed over all contributing trials (for the group of neurons) or the NE condition. The vector for the NE condition contained ones at the times the animal was at the NE and zeroes at all other times (identically for single neurons and groups of neurons). To run the model on groups of neurons of different sizes, we randomly selected 100 different samples of neurons out of the entire population, increasing in sample group size by increments of 5 neurons (i.e., sample unit number: 5, 10, …, modulo(max(number of units), 5) = 0 (nearest multiplication of 5 of the unit count)). The regression performance was estimated using the coefficient of determination . To control for the large number of degrees of freedom when using the whole population of GP neurons as predictors, we compared the model’s calculated value to that of a bootstrap consisting of 100 different shuffles of each neuron’s firing rate.

Cross correlation

Cross-correlograms (CC) were calculated for ±1 s in 1 ms bins for all neuronal pairs that were not recorded from the same electrode (Bar-Gad et al., 2001). A total of 147 CCs from 15 sessions were calculated and analyzed for significance. A CC was considered significant if 3 consecutive bins in the inner part of the CC (±100 ms) crossed the 95% confidence intervals calculated from the outer parts of the CC (>|±500ms|) (Nevet et al., 2007; Oliveira-Maia et al., 2012; Kopelowitz et al., 2014).

Statistical analysis

All the data are presented as mean ± SEM unless otherwise noted. All statistical tests were ANOVAs unless otherwise noted. All post-hoc corrections for multiple comparisons were Bonferroni corrections. To statistically compare the correlation values, correlations were first Fisher transformed, and the statistical comparisons were conducted on the transformed data. To compare R2 values of the linear regression models for different sample sizes of neurons between LFB and HFP, a direct statistical measure was inappropriate due to the dependencies within groups. Thus we used Cohen’s d as a measure of the effect size as follows:

where are the samples mean and where are the samples’ variances. Cohen’s d values exceeding 0.8 were considered a big effect (Cohen, 1988).

REAGENT or RESOURCE	SOURCE	IDENTIFIER
Chemicals, peptides, and recombinant proteins
Ketamine	Zoetis	N/A
Xylazine	Eurovet animal health	N/A
Atropin	SALF laboratory	N/A
Isoflurane	Piramal critical care	N/A
Carperive	Norbrook	N/A
Ringer lactate	Teva medical	N/A
Morphin	Kalcex	N/A
Dental cement	Lang dental www.langdental.com	Cat#1303R5
NaCl	Sigma-Aldrich www.sigmaaldrich.com	Cat#S9888
Formaline	Bio-lab www.biolab-chemicals.com	Cat#06750523F1
Sucrose	Sigma-Aldrich www.sigmaaldrich.com	Cat#S-0389-1KG
S-Acetylthiocholine iodide	Sigma-Aldrich www.sigmaaldrich.com	Cat#A5751-5G
Cupric sulfate CuS045H20	Sigma-Aldrich www.sigmaaldrich.com	Cat#209198-250G
Ethopropazine HCL	Sigma-Aldrich www.sigmaaldrich.com	Cat#E−5406-50MG
Reagent grade alcohol 300 mL 100% ETOH	Bio-lab www.biolab-chemicals.com	Cat#0525050200
Glycine (aminoacetic acid) NH2CH2COOH	Bio-lab www.biolab-chemicals.com	Cat#07132391
Sodium acetate trihydrate CH3COOHNa3H20	Sigma-Aldrich www.sigmaaldrich.com	Cat#S8625
Na2S9H2O	Sigma-Aldrich www.sigmaaldrich.com	Cat#13468-1KG
Glacial acetic acid	Gadot https://he.gadot.com	Cat#64-19-7
HCL	Sigma-Aldrich www.sigmaaldrich.com	Cat#H1758
Xylene	Gadot https://he.gadot.com	Cat#1330-20-71
Experimental models: Organisms/strains
Rats: Long Evans	BIU animal facility	N/A
Software and algorithms
MATLAB 2017b	MathWorks	https://www.mathworks.com/
Offline sorter version 3.3.5	Plexon Inc.	https://plexon.com/products/offline-sorter/
NeuroExplorer version 4.133	Plexon Inc.	https://plexon.com/products/neuroexplorer/
PlexUtil version 4.0.2	Plexon Inc.	https://plexon.com/products/plexutil/
PlexControl version	Plexon Inc.	https://plexon.com/products/cineplex-studio/
CinePlex Studio version 3.5.0	Plexon Inc.	https://plexon.com/products/cineplex-studio/
CinePlex Editor version 3.6.0	Plexon Inc.	https://plexon.com/products/cineplex-editor/
Other
Isonel coated tungsten microwires measuring 35 and 50 μm in diameter	California fine wire https://calfinewire.com	Cat#Cfw-2015376, 2,015,378
Connectors	Omnetics www.omnetics.com	Cat#A79016-001
PCB	A.D.I Electronics Ltd. http://www.adielectronics.com/company/index.htm	PN; ad5254p
Silver wire 0.008″, bare	A-M system	Cat#AM-782000