{'fraction_non_alphanumeric': 0.03359139220574728, 'fraction_numerical': 0.03853387569435332, 'mean_word_length': 4.653808110781404, 'pattern_counts': {'":': 0, '<': 0, '<?xml version=': 0, '>': 0, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Dopamine depletion leads to pathological synchronization of distinct basal ganglia loops in the beta band

Motor symptoms of Parkinson’s Disease (PD) are associated with dopamine deficits and pathological oscillation of basal ganglia (BG) neurons in the β range ([12-30] Hz). However, how the dopamine depletion affects the oscillation dynamics of BG nuclei is still unclear. With a spiking neurons model, we here captured the features of BG nuclei interactions leading to oscillations in dopamine-depleted condition. We found that both the loop between subthalamic nucleus and Globus Pallidus pars externa (GPe) and the loop between striatal fast spiking and medium spiny neurons and GPe displayed resonances in the β range, and synchronized to a common β frequency through interaction. Crucially, the synchronization depends on dopamine depletion: the two loops were largely independent for high levels of dopamine, but progressively synchronized as dopamine was depleted due to the increased strength of the striatal loop. Our results highlight the role of the interplay between the GPe-STN and the GPe-striatum loop in generating sustained β oscillations in PD subjects, and explain how this interplay depends on the level of dopamine. This paves the way to the design of therapies specifically addressing the onset of pathological β oscillations. Author summary Parkinson’s Disease is associated to the death of neurons generating a particular neurotransmitter: the dopamine. Motor symptoms of PD, on the other hand, are known to be due to dysfunctions in a particular subcortical area of the brain, the BG network. In particular, the BG network develops pathological oscillations in a specific frequency range (β: [12-30] Hz). What is unclear is how dopamine depletion leads to these oscillations. In this work we developed a BG network model and we found the actual reason for these abnormal oscillations is the synchronization of two loops within the network that are individually oscillating in the β range. For healthy level of dopamine the two loops are decoupled and the oscillation power is low. When dopamine is depleted (as in PD) the two loops synchronize and originate the pathological oscillations associated with motor symptoms.


Introduction
The third hypothesis is the one of a prominent role of the interaction between GPe 48 and STR: experimental recordings in rats [38][39][40][41][42] and mice [43][44][45][46] suggest that the GPe 49 nucleus and the pallido-striatal pathway could play a major role. In vivo observations of 50 pauses in fast spiking neurons FSN support this hypothesis as transient drop in the 51 activity might arise from interactions with GPe [47] (see also [48,49] for an extended 52 review). From a computational point of view, a model of the GPe-striatal circuit has 53 been developed by Corbit et al [50], who investigated the onset of pathological β 54 oscillations in a subnetwork consisting by the inhibitory synapses FSN → D2-MSN → 55 GPe. 56 In this work we propose a comprehensive model including all the loops with the aim 57 of addressing the aforementioned hypothesis on the mechanisms underlying abnormal β 58 oscillations in PD. In line with this purpose, the starting point of our work is the 59 network developed by Lindahl and Kotalewski [51] which includes all major nuclei and 60 connections of the BG. In this context, we identify the major sources of β activity, 61 investigate the role of their interaction and highlight the effects of dopamine depletion. 62 Materials and methods 63 1.1 Basal ganglia network model 64 Our model of the basal ganglia (Fig 1) includes: 65 • three striatal populations: D1-type dopamine receptor Medium Spiny Neurons 66 (D1); D2-type dopamine receptor Medium Spiny Neurons (D2) and a population 67 of Fast Spiking (inter-)Neurons (FSN); 68 • the external Globus Pallidus divided into two populations that will be labelled as 69 GPe-TA (characterized by a lower discharge rate and by a negligible input from 70 striatal populations) and GPe-TI (with a higher activity and receiving input from 71 D2); 72 • the SubThalamic Nucleus (STN). 73 Each population presents specific size (see Table 2) and neuron model (see section 74 1.2). Single neurons within connected populations are randomly associated according to 75 specific connection probabilities (see Table 3). Each neuron receives inputs from within 76 the network and from other brain regions (mainly the cortex) which have not been 77 explicitly included in our model. These external inputs have been modeled as 78 independent poissonian trains of pulses. The mean rates ν ext of these signals have been 79 adjusted for each population in order to ensure realistic population firing rates: FSN 80 [10][11][12][13][14][15][16][17][18][19][20] Hz ( [52], [48]), D1 and D2 (MSN) [0.05-2.5] Hz ( [53]), GPe-TI  Hz 81 ( [38]), GPe-TA [5][6][7][8][9][10][11][12][13][14][15] Hz ( [38]) and STN [12][13][14][15][16][17][18][19][20] Hz ( [24]).

82
Note that the GPi and the SNr nuclei were not included in the simulated model Architecture of the adopted Complete Model of the basal ganglia. FSN, D1 and D2: striatal Fast Spike Neurons, and medium spiny neurons with D1 and D2 dopamine receptors; GPe-TA and GPe-TI: globus pallidus externa type A and type I; STN: subthalamic nucleus; GPi: globus pallidus interna; SNr: substantia nigra pars reticulata; ext: external poissonian input. By convention red/green arrows are inhibitory/excitatory projections. Dashed elements have not been included in our model.

86
All neurons in the adopted network are adaptive, conductance based, point neurons ( [54][55][56][57]). The STN, GPe-TI and GPe-TA populations are modelled as adaptive October 7, 2022 3/21 exponential neurons (aeif cond exp model in the code implementation) [58] and the dynamics of their membrane potential is governed by: The striatal populations are modelled as adaptive quadratic neurons (aqif cond exp or aqif2 cond exp model) [59] and their membrane potential evolves according to: In both cases excitatory and inhibitory conductances present exponential decays: The evolution of the adaptation variable w is governed by ( [60,61]): in the aeif cond exp and aqif cond exp model, while by: in the case of the aqif2 cond exp model.

87
To introduce variations in each neural population, the synaptic weight of the 88 external poissonian input is not equal in all neurons. Rather, it is assumed to be 89 uniformly distributed around its central value with an amplitude equal to the 90 dev-ext-weight parameter. The adopted values of the parameters of each neuronal 91 population are summarized in Table 2. The connectivity properties (delays, connection 92 probabilities and synaptic weights) are reported in Table 3. These values have been 93 adapted from the work of Lindhal and Kotaleski [51], here simplified neglecting synaptic 94 plasticity, spatial restrictions on the connected populations and dopamine effects 95 different from the ones described in section 1.5.

96
Note that the results of the model are not critically dependent on the exact values of 97 the network parameters.

98
In the following, we will refer to the described network of the basal ganglia as the 99 Complete Model. In order to identify the more relevant structures generating β oscillations, the following 102 strategy has been pursued. Starting from the Complete Model (Fig 1), the relevance of 103 each connection in generating oscillations in the β regime has been quantified referring 104 to the decrease in the mean intensity of β oscillations following its elimination. In   Table 2. Adopted values of the parameters of each neuronal population.
• an auxiliary subnetwork S* has been introduced with the same neuron model and 108 size of S but with no input sources other than the external one;

109
• the parameter I e of the S* subnetwork has been adjusted so that the mean  As from our analysis emerged that two oscillators are mainly responsible for the onset of 123 the pathological activity, we developed an ad hoc network to analyze the properties and 124 the consequences of their interaction (see Fig 4A).

125
The GPe-TI nucleus is present in both oscillators, hence its neurons have been   In contrast, in order to modulate the coupling between the two oscillators, inter -loop 132 connections (thin continuous line in Fig 4A) present connection probabilities which are 133 modulated by a coupling parameter ε. Particularly, for each pair (S,T ), the connection 134 probability of each S neuron with each T neuron is set equal to: where p 1,S→T is the probability characterizing the corresponding connection in the 136 Complete Model (see Table 3).

137
For every connection, synaptic weights are the same as in the Complete Model (and are not affected by ε). As a result, for ε = 0 the two loops are completely independent, while for ε = 1 the network presents realistic connectivity properties. Further, in order to ensure that the overall inputs to each nucleus remain the same when ε is changed, the auxiliary populations GPTI*, D2* e STN* have been introduced. Neurons in these nuclei receive constant and independent external input such that their mean spiking rates are equal to the mean spiking rates of the corresponding non-auxiliary populations. The auxiliary connections are shown as dashed lines in Fig 4A. While their synaptic weights are the same as in the original network, the connection probabilities S * → T are set equal to: Hence, for ε=1 the auxiliary nuclei do not play any role; otherwise, they become study the process of synchronization in these states, it has been necessary to modulate 149 the intensities of the connections constituting the two loops.

150
The complete list of the introduced adaptations are listed in Table S1. In the 151 following, and in contrast to the original Complete Model (see Section 1.1), the network 152 obtained after the listed assumptions and simplifications will be labelled Simplified In order to investigate the role of dopamine and accounting for its effects on the mean 156 spiking rate of neurons into striatal populations ( [10]), the severity of dopamine 157 depletion has been modelled in the Simplified Network by the modulation of the 158 intensity of the external input rate towards the D2 population. More precisely, the 159 external input rate for this population (and the corresponding auxiliary population D2*) 160 has been modulated according to: where ν ext,1 is the reference value of the input rate (see Table 2) and D d is the 162 parameter regulating the severity of the condition of dopamine depletion (the higher is 163 D d , the more severe is the condition). Throughout the analysis, D d varies in an interval 164 such that the D2 activity ranges from ∼ 0.05 to ∼ 4 Hz. The effects of dopamine 165 depletion will be also studied in the Complete Model applying the same modulation. • its mean frequency of oscillation: • its mean value in the β regime: where P i (f ) is the PSD of the nucleus activity. The values of m and M have been fixed 179 to 8 and 24 Hz respectively, within the β regime and coherently with the natural 180 frequencies of the two oscillators.

181
The Mean PSD i is a measure of the intensity of β oscillations within each considered nuclei. We highlight that this measure is biased : in the presence of constant activity (no β oscillations), fluctuations around the mean are present and the related spectrum is affected by noise. In order to eliminate this bias we considered the corrected quantity: where Q(ν 0,i , N i ) is the Mean PSD of the activity of a population of N i neurons with constant mean activity ν(t) = ν 0,i . The value of Q(ν 0,i , N i ) has been estimated computing the PSD of a fictitious activity in which the number {n i j } j of spiking neurons in each bin is given by independent extractions from a Bernoulli distribution: In order to verify whether the properties that we obtained persist in the limit of large populations (the realistic size of each population is approximately 10 3 times the adopted sizes [63,64]) we repeated all the analysis with each population size multiplied by a factor n ∈ [1,2,4]:  Since both loops contain the GPe-TI population, in the Complete Model of the basal 219 ganglia the two loops interact with one another. The consequences of this interplay will 220 be the objective of the next sections.

Loops interplay 222
The results highlighted in the previous section suggest that the main responsible for the 223 onset of β activity in the adopted model of the basal ganglia are the two oscillators set 224 up by the STR and STN loops. In order to deepen the properties of the interaction 225 between these oscillators, we now consider the Simplified Model (Fig 4A). In this 226 section, we particularly focus on the consequences of increasing the value of the coupling 227 ε which regulate the intensity of the connections between the two oscillators (see 228 Methods section 1.4).

229
Interplay variability resulted in fundamental changes in the spectral dynamics of the 230 network (Fig 3): for very low values of ε the Str and STN loops, measured through the 231 D2 and STN nuclei respectively, exhibit a power spectrum which is not dissimilar from 232 the case of the independent oscillators (ε = 0) (Figs 2B and 2C). As the value of the ε 233 parameter is increased, the PSDs of both populations become affected by the presence 234 of the other resonance and present a broader peak which increasingly include the  In order to highlight this behaviour, we computed the average frequency of the STR 240 and STN loops (with specific measurements respectively at the D2 and STN nuclei) as a 241 function of ε (see Methods section 1.6). As expected, the two frequencies start from the 242 natural frequencies of the two oscillators (∼ 13 Hz and ∼ 19 Hz in the region of ε = 0) 243 and, as the interplay between the two loops increases, converges to an intermediate 244 common frequency f ∼ 16 Hz (ε = 1, see Fig 4B). Even more interestingly, the process 245 of synchronization is associated to a strong increase in the intensity of the activity in 246 the β regime ( Fig 4C).

247
In order to verify whether the properties highlighted above persist in the limit of 248 large populations, we repeated the analysis with the size of each population multiplied 249 by a factor n ∈ [2, 4] (see Methods sec 1.7). The results of this analysis highlight two 250 different regimes (Fig 4D):

251
• in case of low synchronization (ε ∼ 0) the intensity of β oscillations goes to zero in 252 the large n limit;

253
• in case of high synchronization (ε ∼ 1) the intensity of β activity is preserved.

254
On the one hand, these results confirm that the establishment of a state of complete 255 synchronization and the consequent onset of prominent β activity in the model are 256 stable for variations in populations size; on the other hand, they suggest that single loop 257 models are not capable of explaining strong β activity in the limit of large populations. 258 A further insight on the process of synchronization is finally captured by the value of 259 PSD STN (f STR ) (i.e., the power spectral density of the STN population computed at the 260 natural frequency of the STR loop; Fig 4E). Particularly, for low coupling, this value is 261  low as STN uniquely oscillates at the natural frequency of the STN loop (Fig 3A). In 262 the intermediate regime (ε ∼ 0.5) higher values are registered since the activity of the 263 STN nucleus becomes more affected by the STR oscillator ( Fig 3B). Finally, for ε ∼ 1, 264 the value of PSD STN (f STR ) decreases again as the two oscillators synchronize to a novel 265 frequency which is different from both the natural ones (Fig 3C).  nucleus moves from ∼ 16 Hz, that is mean value between 8 and 24 (see Methods and Eq. 288 (4)), to the natural frequency of the D2 oscillator at ∼ 13 Hz (see Fig 5C).

289
The role of dopamine in the Simplified Model We continue considering the 290 Simplified Model but we now introduce the coupling between the two oscillators by 291 setting the value of ε = 0.75.

292
In contrast to the previous case, the increase in the intensity of the STR Loop leads 293 the two oscillators to undergo a process of synchronization: for low values of D d the two 294 oscillators interact poorly (the STR oscillator is almost absent and only the STN 295 resonance at ∼ 18 Hz characterizes the status of the system; see Fig 6A); as the value of 296 D d is increased, the enhancement in the intensity of the STR oscillator implies higher  Similarly to the increase of ε (section 2.2), the increase of D d determines the 300 convergence of the mean frequencies of oscillation of the STN and D2 populations (Fig 301  7B) and a remarkable increase in the intensity of the pathological oscillations (Fig 7C). 302 Further, the analysis of the limit of large populations shows that:

309
In accordance with previous results, this proves that the intensity of β activity in the

312
Beside that, the comparison between the uncoupled (ε = 0.00) and coupled 313 (ε = 0.75) cases highlights that the increase in the intensity of β activity is stronger in 314 the latter case (Fig 7E).   Model of the basal ganglia (see Fig 8A). Particularly, even in this model, as the value of 330 D d is increased, the mean frequencies of the D2 and STN populations converge (Fig 8B) 331 and the intensity of β oscillations increases and presents the characteristic speed up due 332 to completion of the synchronization process ( Fig 8C).

333
Finally, and in accordance with previous results, the size analysis confirms that our 334 model is capable of capturing prominent β activity also in the limit of large populations 335 (Fig 8D). 336 Please note that in all the studies regarding the effects of dopamine depletion (Figs 5, 337 6, 7B, 7C, 8B and 8C) we showed the results related to the case n = 8. This has been 338 necessary to capture the relevant properties of the system in a clearer way (for example, 339 the reader can compare the speed up in the increase of β activity in case of n = 1 and 340 n = 8 shown in Fig 7E). This choice does not affect the results of the large n limit 341 analysis.

346
One of the most important feature of the pathological β activity in the basal ganglia 347 is that it presents a strongly phasic characterization ( [19,20]). We characterized the interplay between two loops generating β oscillations within the 361 basal ganglia: one between GPe and striatum and one between GPe and STN. While for 362 high values of dopamine the two loops are both weak and act largely as independent, 363 when dopamine is depleted the striatopallidal loop increases its strength, and this leads 364 to synchronization and mutual reinforcement between the two loops and eventually to 365 strong oscillations in the whole basal ganglia network. We suggest that this mechanism 366 plays a key role in the generation of abnormal β oscillations in Parkinson's Disease.

367
Our findings are coherent with state-of-the-art experimental results, and shed light 368 on the underlying mechanisms. D2 neurons are known to increase their firing activity 369 and to strongly synchronize in the β range in dopamine-depleted conditions [13]. Our 370 modeling work suggests a mechanism in which the former effect leads to the latter due 371 to the internal resonances of the BG network. The relevance of GPe in the generation of 372 β oscillations was recently proved by a work of Mallet group [42]. First, the cortical 373 origin of β oscillations was ruled out by observing that optogenic cortical silencing did 374 not affect β power. Interestingly, also STN activity suppression had a small effect on β 375 power, suggesting that STN might only have, in the words of the authors "a supportive 376 role" in the generation of abnormal β oscillations, which where instead suppressed by The key difference is that in our 385 view this loop is not simply modulated by the intensity of the striatal inhibition to GPe, 386 but rather by the intensity of the β oscillations of the striatopallidal loop.

387
Most previous computational models on the onset of β oscillations (e.g., [26,34,67]) 388 described how the STN-GPe loop was originating the pathological oscillations but 389 focused only on perturbations induced by modulations of the intensity from striatal 390 population to pallidum [29], without taking into account the pallidostriatal feedback 391 that leads to the striatopallidal β loop. In a relevant work, Corbit and colleagues 392 introduced a detailed small scale model of the striatopallidal loop [50] and showed how 393 it generated β oscillations. Later models, including [51] which was the model laying the 394 ground for this work, implemented then this loop in their network descriptions 395 ( [68,69]), but did not explicitly address how the two loops interacted in originating β 396 oscillations, rather focusing on action selection. Interestingly, a recent study by the 397 same group focusing on transient responses of the basal ganglia network [70] also 398 captured β activity, but did not address the question of its origins.

399
In a key paper, Tachibana et Al [71] studied the effects of (locally) blocking different 400 connections of the STN-GPe circuit on the intensity of β activity. As accurately 401 described in [72], the strong constraints highlighted in [71] rejected most of the previous 402 hypotheses on the origin of β oscillations (see Introduction). Differently from those 403 cases, thanks to the presence of two loops and to the effects of their synchronization 404 (which are not affected by local blocks of the connections), the results of our model are 405 coherent with the evidences in [71]. Further, and differently from the mechanisms 406 proposed by Pavlides et al [72], our model does not need to include cortical populations. 407 As a final remark, the presented size analyses highlight that the intensity of β 408 activity is preserved when the two oscillators are in the synchronous regime. On the one 409 hand this confirms that our model is cable of reproducing β activity also in the limit of 410 large populations; on the other hand it suggests that single loop models ( [29,34,50]) 411 may not be capable of explaining prominent β activity.

412
In future works we will investigate how our understanding of the mechanisms agonists, but to proper model this latter aspect we have to take into account several 417 other pathophysiological aspects [73]. Moreover, it would be interesting to simulate the 418 action of Deep Brain stimulation as in [67,[74][75][76] also in our network to evaluate the 419 input of different targets on suppressing the β oscillations and compare it to 420 experimental results [77][78][79].