Impact of Network Topology on Neural Synchrony in a Model of the Subthalamic Nucleus-Globus Pallidus Circuit

Synchronous neural oscillations within the beta frequency range are observed across the parkinsonian basal ganglia network, including within the subthalamic nucleus (STN) - globus pallidus (GPe) subcircuit. The emergence of pathological synchrony in Parkinson’s disease is often attributed to changes in neural properties or connection strength, and less often to the network topology, i.e. the structural arrangement of connections between neurons. This study investigates the relationship between network structure and neural synchrony in a model of the STN-GPe circuit comprised of conductance-based spiking neurons. Changes in net synaptic input were controlled for through a synaptic scaling rule, which facilitated separation of the effects of network structure from net synaptic input. Five topologies were examined as structures for the STN-GPe circuit: Watts-Strogatz, preferential attachment, spatial, stochastic block, k-regular random. Beta band synchrony generally increased as the number of connections increased, however the exact relationship was topology specific. Varying the wiring pattern while maintaining a constant number of connections caused network synchrony to be enhanced or suppressed, demonstrating the ability of purely structural changes to alter synchrony. This relationship was well-captured by the algebraic connectivity of the network, the second smallest eigenvalue of the network’s Laplacian matrix. The structure-synchrony relationship was further investigated in a network model designed to emulate the action selection role of the STN-GPe circuit. It was found that increasing the number of connections and/or the overlap of action selection channels could lead to a rapid transition to synchrony, which was also predicted by the algebraic connectivity.

characteristic feature of Parkinson's disease (PD).Synchrony in this band is associated with the severity of bradykinetic and rigidity motor symptoms [1], [2], [3], and is reduced by both Levodopa therapy and deep brain stimulation [4], [5], [6].The reciprocally connected circuit comprising the subthalamic nucleus (STN) and Globus Pallidus externus (GPe), modulated by cortical and striatal inputs, is believed to play a central role in the generation and maintenance of pathological beta activity in PD [7], [8], [9].Alterations in neuronal properties [10], [11], and in the coupling strengths between nuclei in the basal ganglia [12], [13] have been linked to the emergence of synchrony, however the influence of neural connectivity pattern is less clear.Changes in synaptic connectivity have been observed in animal models of the Parkinsonian STN-GPe circuit.In particular, an increase in the number of GPe-STN synapses [14], [15], and a reduction in the number of cortico-STN [16], [17] and STN to GPe synapses [18].As changes in connectivity can impact the synchronizability of coupled dynamical systems [19], an understanding of how topology influences synchrony in the STN-GPe circuit is critical in understanding the mechanisms underlying pathological neural activity in PD.This is relevant not just from a physiological perspective, but also for computational modelling as the choice of connection arrangement in spiking neuron simulations implicitly affects the synchronizability of the neural system.Computational models can provide insight into generation of pathological oscillations and also the mechanisms by which they are suppressed by therapies such as deep brain stimulation [20].Understanding the contribution of factors such as network topology to neural synchrony is thus critical in understanding the generation and suppression of pathological oscillations in PD.
The effect of connection architecture on firing patterns in the STN-GPe was previously explored in a landmark computational study by Terman and Rubin [21], which illustrated how the circuit can be tuned to exhibit a range of activity patterns including synchronized bursting, by adjusting the pattern and weight of synaptic connections.However, the networks investigated were small in scale, and only three network structures were considered.It is not clear, therefore, if these results generalize to larger networks or different connectivity patterns, such as small world and scale free topologies which have been identified in micro-scale connectomes [22], [23].The synchronizability of neuronal networks has been explored previously, in contexts such as the respiratory brainstem [24] and suprachiasmatic nucleus [25], and for more general areas [26], [27].However, it is not clear how these findings apply to the synchronization of neural activity in the parkinsonian basal ganglia.
In this study, the effect of network topology on neural synchrony is examined in a spiking neuron model of the STN-GPe circuit.The number of neurons and average number of connections was varied first to investigate the effects of network size and density on synchrony, while maintaining mean firing rates across conditions.To isolate the effects of network structure, the network size and number of connections were then held constant, and the arrangement of connections was varied in each topology.Utilizing the framework provided by the master stability function (MSF), [19], [28] which enables network properties to be investigated independent from nodal (neuronal) dynamics, algebraic connectivity was investigated as a predictor of synchrony across various network topologies.Finally, the structure -synchrony relationship was explored in a network model based on the proposed action selection role of the basal ganglia.

II. METHODS
The model consisted of equal sized populations of STN and GPe neurons.The STN receives excitatory inputs from the motor cortex and inhibitory inputs from the GPe, and the GPe receives excitatory inputs from the STN and inhibitory inputs from the both the striatum and GPe, similar to previous studies [29], [30], [31] Fig. 1(a).Conductancebased neuron models, originally presented in [29] and [32] were used to simulate the STN and GPe dynamics, whereas cortical and striatal outputs were modelled with Poisson spike generators.Several distinct classes of neurons are present in the GPe [33], but only type A, parvalbumin expressing neurons are considered here, which project predominantly to the STN [18].Connections projecting from the cortex and STN were modelled as excitatory AMPA synapses, while projections from the GPe and striatum were inhibitory GABA A .The weight of each synapse was scaled as a function of the number of connections incident on the postsynaptic neuron to maintain the firing rates observed in PD rats of 34.0 ± 3.4 Hz for the STN [34] and 14.1 ± 0.5 Hz for the GPe [35].Network sizes varied from n = 200, 1000, and 2000 neurons in the combined population.

A. Neuron Models
Neurons were represented using the Hodgkin Huxley formalism where membrane potential (V m ) is related to membrane capacitance (C m ) and ionic currents through the differential equation: I ion is the set of ionic currents listed for each neuron in supplementary materials Section A, I leak is the membrane leakage current, I x syn is the synaptic input current from the x th neuron.I bias is a static bias current and I noise is a time varying Gaussian noise current with zero mean and a standard deviation of 0.05 pA, with a value switch occurring every 2 ms.I async is a small charge balanced oscillating current applied to the STN and present for the first 1000ms of the simulation to initialize out of phase neural firing.Parameter values were identical for each neuron within a population.I x y syn represents the synaptic current due to both the Poisson spike generators and the conductance-based neurons.The dynamics of these AMPA and GABA A currents are described by [36]: where R x y represents a first order process that captures the kinetics of the onset and decay of post-synaptic current following a presynaptic spike at the sending neuron x, E r ev is the reversal potential of the synapse, and w x y is the coupling strength from neuron x to neuron y.Delays of 4 ms were included between STN and GPe [37], and 2 ms for inter-GPe connections.No delays were present between cortex and STN or striatum and GPe as these spikes were generated from a time invariant process.

B. Synchrony
The global synchrony of each neural population was quantified using the χ p measure introduced in [38]: Authorized licensed use limited to the terms of the applicable license agreement with IEEE.Restrictions apply.
where σ 2 V is the time-averaged variance of the ensemble averaged membrane potential of all neurons in a population, and σ 2 V i is variance of membrane potential of the i th neuron.This captures the average fluctuations of the neuronal membrane potentials, normalized to lie within the interval [0, 1], and can be used to infer synchrony from bursting neurons where spikes are not necessarily time locked [38].For activity in the STN-GPe circuit to be considered synchronous, synchrony was required to be present in both subpopulations, and thus synchrony was estimated as the geometric mean of individual STN χ ST N and GPe χ G Pe synchrony: Although neural synchrony and neural oscillations can be distinct phenomena, in this model synchrony is present as correlated bursts at approximately within the beta frequency range.Here, χ thus captures the phenomenon of beta oscillatory activity observed in the STN-GPe circuit in PD.

C. Determining Scaling Functions for Connection Weights
To maintain physiological neural firing rates, and facilitate comparison between networks without the confounding effect of firing rate differences, synaptic inputs for each neuron were scaled as a function of in-degree (number of incoming connections).The weighting functions for each connection set, STN to GPe (F SG (k in )), GPe to STN (F G S (k in )), and GPe to GPe (F GG (k in )) were determined by fitting rational polynomials to data generated from a simulated annealing algorithm [39] optimized to produce weights which maintained mean firing rates within both populations as in-degree was varied.
For a given k in , candidate sample weights W * were generated using a log normal distribution with mean equal to the previous accepted weight W (u−1) , with indicating the previous iteration and initial values W 0 SG = 0.1, W 0 G S = 1, W 0 GG = 0.1.For each iteration, three simulations were run with these weights on a k-regular random network with k = k in , and their average cost C (W * ) calculated.The cost function was based on the differences between the desired mean firing rates (34 Hz, 14.1 Hz) and the actual mean firing rates (F S , F G ) for the STN and GPe, along with a term S that penalised low levels of synchrony: where χ * is target level of synchrony that increases to a maximum of 0.45 at k = 10 then remains constant.This target value was chosen as it was typically the highest synchrony value observed during parameter sweeps.Without this moving target, the algorithm prioritizes synchrony at the expense of firing rates for low k.c = 5 is a constant that ensures that S is similar in magnitude to F s and F G .
The probability α of accepting a new weight is: where T is the annealing temperature that decreases to 0 as the iteration number u reaches its maximum u max = 500: The probability of accepting transitions to weights with higher costs decreases with T , and the resulting collection of weights converges to a minimum in the parameter space.
A non-linear least squares algorithm was used to fit a rational polynomial of the form: to the weights

D. Network Topologies
Five network generation rules, used to specify connections between STN and GPe neurons, were chosen to produce qualitatively distinct topologies, and capture features observed in neuronal networks.Insights gained from their analysis were explored in a sixth rule based on a physiologically plausible network structure specific to the STN-GPe which can emulate the action selection role of the basal ganglia [40].The number of nodes n and average degree k were varied for each rule, along with a separate wiring parameter that controlled the pattern of connections within a rule.The lack of inter-connectivity among STN neurons places a topological constraint on applying traditional network models to the full STN-GPe circuit.To account for this, each arrangement of pre-and post-synaptic connections SG, GS, and GG, was treated as a sub-network, and the network generation rule applied to that sub-network, Fig. 1(b).These sub-networks were then combined together to represent the full STN-GPe network, Fig. 1(c), from which the algebraic connectivity was calculated.When generating the networks care was taken to ensure that the networks did not contain self-or multi-edges.
1) Watts-Strogatz Networks: The Watts-Strogatz small world generation rule interpolates between a structured ring network and a random graph by adjusting a rewiring parameter β ∈ [0, 1] [41].Each node is connected to its k 2 neighbors on either side, and every edge is randomly rewired with probability β to have a different target node selected uniformly at random from all other nodes.For low but non-zero values of β the network displays "small world" properties, where local clustering is high, but the average distance between neurons is reduced due to the addition of shortcuts.These small-world properties are commonly observed in neuronal networks across varied species and scales [42].
2) Spatial Nearest Neighbor Networks: Within each spatial nearest neighbor network each node is assigned a point in twodimensional space, the Euclidean distances between all points are calculated, and these distances are indexed in ascending order.A probability p(i) is associated with each index i: This probability distribution is then used to sample kn edges, giving a network with an average degree of k. γ controls the non-linearity of p: at γ = 0 the distribution is uniform with no spatial dependence, and as γ increases, the selection of nearby edges becomes prioritized.Due to the energy cost of maintaining long axons, brain states are often characterized by connections between spatially neighboring neurons [43].
3) Preferential Attachment Networks: Preferential attachment networks were examined due to their distinct qualities: heavy tailed degree distributions and the presence of highly connected hubs.The generation method is based on a modified version of the Barabasi-Albert model [44], [45].An initial cluster of k + 1 all-to-all connected nodes is created and then grown to the desired size n by connecting new nodes to existing ones with each existing node having probability p i of receiving a new connection: k i is the degree of existing node i.The exponent α renders the attachment non-linear and allows control over the tails of the distribution.For 0 < α < 1 the preferential attachment is said to be "sub-linear", and the degree distribution tends towards a stretched exponential.For α = 1 the attachment it is "linear", and the original Barabasi-Albert model is recovered with a scale free degree distribution which follows the power law P(k) ∼ k −3 .For α > 1 the preferential attachment is "super-linear" and the hubs receive the majority of new edges [45].
4) Stochastic Block Networks: The stochastic block topology was investigated to examine the effects of community structure on synchrony [46].The generation rule involved dividing the network into M = 4 subgroups, with N M nodes per group, and an average of k connections per node.The parameter p e controlled the probability of edges occurring between nodes of different groups.As p e is varied the network goes from containing only edges between nodes of the same group at p e = 0, to a random network at p e = 1 − 1 M = 0.75, to a network containing only edges between different groups at p e = 1.
5) K-Regular Random Networks: K-regular random networks, are a subclass of the more general Erdős-Rényi random networks subjected to the condition that each node has the exact same degree [47].Because of this property, they were used to investigate the relationship between synaptic strength and number of connections for constant firing rates.
6) Action Selection Networks: For the action selection network model, an approach based on the network used in [29] was adapted whereby connections could be generated stochastically.The model was parameterized by the average number of connections k, and the channel spread σ .The size of the channels was controlled by a variable d = 40, the distance between the centres of excitation and lateral inhibition, Fig. 2. In [29], the STN projects to GPe neurons within its own channel, while the GPe projects to STN neurons in adjacent channels and to GPe neurons in the current and adjacent channels.This structure was implemented in the stochastic model by specifying that for STN to GPe connections, presynaptic neurons i were chosen at random to connect to postsynaptic neurons j, such that j ∼ N i, σ d 6 , with the process repeated until the average degree k was reached.A similar process was used for GPe to GPe connections but with the channel spread increased j ∼ N i, σ d 3 , to account for overlap with adjacent channels.GPe to STN connections were established with a bimodal distribution for the postsynaptic connection probabilities: j ∼ N i + x, σ d 6 where x = ±d is a variable equal to d or −d with probability of one half.Changing k changes the density of connections within the channels, whereas changing σ changes the amount of channel overlap, Fig. 2(b).The connectivity structure is shown in Fig. 2(c) for connections going to and coming from a single GPe neuron.

E. Algebraic Connectivity
The MSF provides a framework for relating graph spectral properties to the synchronizability of a system of coupled oscillators.In particular, the eigenvalues of a networks Laplacian can be used to determine the stability of the synchronous state [28].Algebraic connectivity, λ 2 , was calculated as the smallest nonzero eigenvalue of the indegree normalized Laplacian matrix describing the STN-GPe network.Let L denote the Laplacian matrix of the STN-GPe network, with entries L i j = −1 if a connection exists from node i to node j, diagonal entries L ii = k i equal to the out degree so that the row sums are zero, and L i j = 0 otherwise.The in-degree (k i ) normalized Laplacian Gis a matrix such that G i j = L i j k i [48].λ 2 was then calculated by ordering the Authorized licensed use limited to the terms of the applicable license agreement with IEEE.Restrictions apply.
eigenvalues of G by real part and taking the first non-zero term.

F. Implementation
Simulations were performed using the NEURON software [49] backend with the Pynn API package [50] on the University College Dublin Sonic high-performance computing cluster, with a timestep of 0.02 ms.Network edge lists were created using custom scripts for the preferential attachment, spatial topologies, and action selection topologies, and using the network [51] python library for the Watts-Strogatz, k-regular random, and stochastic block topologies.These edge lists were used in conjunction with Pynn to generate connections for the modelled STN-GPe circuit.
To explore the relationship between synchrony and average degree, the wiring parameters for each topology were first set to fixed values.β = 0.01 was chosen as it places the Watts-Strogatz networks in the small world region [41], and P e = 0.01 was selected for stochastic block networks as it has a similar interpretation to β.It was found during initial runs that α = 6 induced a reasonable amount of connection locality in the spatial networks.For preferential attachment networks, α = 1 gives the canonical scale free network.Average degree was varied in logarithmically spaced values from k = 2 to k = 90, 350, or 500 for n = 200, 1000 or 2000 respectively.Simulations were conducted for 10,000 ms with the first 1,000 ms discarded to avoid transient effects.To determine the effect of wiring pattern on synchrony, sample networks were generated with a fixed average degree of k = 10 and size n = 1000 by varying the wiring parameters [β, γ , α, P e ] for each topology over a range of values.The frequency of oscillations was estimated as the frequency at which the coherence between the STN and GPe signals was maximum from data pooled across all simulations with χ ≥ 0.45.

A. Scaling of Synaptic Connection Weights
The optimization algorithm identified connection weights which maintained STN and GPe mean firing rates within the target ranges as the degree was varied Fig. 3. Mean firing rates lay within 1 Hz of the target for 4 < k < 240 in the STN, Fig. 2(a), and k > 3 in the GPe, Fig. 3(b).The level of synchrony within the network remained close to target for k > 10, Fig. 3(c).The strength of connections from STN to GPe and GPe to STN decreased monotonically with increasing degree, Fig. 3 (d) and 3 (e).In contrast, the strengths of GPe to GPe connections increased with degree initially, reaching a maximum at k = 14, after which they decreased, Fig. 3(f).The parameters governing each weighting function are presented in Table I.
Using the connection weights identified in section A, the average degree and the connectivity patterns were then varied, while keeping the mean firing rates approximately constant.

B. Influence of Average Degree
Global synchrony generally increased with average degree for the topologies considered, however the exact nature of

TABLE I COEFFICIENTS OF RATIONAL POLYNOMIALS GOVERNING CONNECTION WEIGHTS
these relationships varied greatly depending on the network size and topology Fig. 4(a).The effect of network size on synchrony was evident for Watts-Strogatz and stochastic block topologies, with the transition to the synchronous state occurring at higher degrees for larger networks, Fig. 4(a).In contrast, preferential attachment and spatial topologies showed little dependence on network size, reaching their maximum values near k = 30 for all networks.
The asynchronous state was characterized by uncorrelated activity in the raster plots which averaged to a steady aperiodic mean membrane potential Fig. 4(b).The highly synchronous state consisted of correlated firing in the STN and GPe, shown in Fig. 4(c) for spatial network with n=200 k=30.The frequency of oscillations within highly synchronized networks, χ > 0.45, was found to be 20.7 ± 0.2 Hz, independent of network topology or degree.The synchronous oscillations thus lay within the beta frequency range, consistent with pathological beta oscillations observed within the parkinsonian basal ganglia network.
In the larger networks examined, low synchrony states were observed for Watts-Strogatz and stochastic block networks up to relatively large degrees, Fig. 4(d) and 4(e).These exhibited distinct activity patterns which were locally, but not globally, correlated.Due to the dependence of connectivity on local neuron index, waves of activity emerge in the raster plots of Watts-Strogatz networks which average out in the mean field resulting in minimal global synchrony, shown in Fig. 4(d) with n=200 k=5.Similarly, stochastic block networks conferred a local 'block' structure to the spiking activity.Firing was correlated within a block, but the frequency and phase of different blocks were uncorrelated due to the limited connections between them.For synchronous and asynchronous block networks the dominant frequency was 20.7 Hz, but in asynchronous block networks this frequency varied with a mean range of 0.32 Hz which led to the blocks moving in and out of phase with one another, limiting global synchrony.The lack of between block interaction was also evident in the mean pairwise coherence of the blocks mean Vm, which was 0.101 at 20.7 Hz in asynchronous networks, compared to 0.993 in synchronous networks.

C. Connectivity Pattern and Algebraic Connectivity
It was possible to tune each topology to display a range of synchrony values by adjusting their respective wiring parameters, though the extent of this range varied over each topology.The exception was the k regular random network which does not have an additional rewiring parameter.Watts-Strogatz networks became more synchronous as the network structure evolved from ring-like at low β to uniform random at β ≈ 1, seen in Fig. 5(a).Spatial networks became less synchronous as the distance exponent increased, and local connections were prioritized over global ones Fig. 5(b).Preferential attachment networks always displayed a high level of synchronous beta activity even without preferential Authorized licensed use limited to the terms of the applicable license agreement with IEEE.Restrictions apply.attachment in the growth process, α = 0, Fig. 5c.These α = 0 networks still contained hubs, but these hubs were less influential and received a lower proportion of the total connections.As inequality in the degree distribution grew with increasing α, synchrony increased to levels higher than those observed in the other network models.Finally, higher proportions of inter-block connections p e resulted in greater levels of synchrony for stochastic block networks, shown in Fig. 5(d).
When algebraic connectivity was compared against the measured synchrony over all networks, synchrony was observed to increase monotonically and displayed a consistent trend across topologies.A sigmoid function sig (x) = a 1+exp(−b+cx) was fit to the data and gave a coefficient of determination across all network topologies and wiring parameters of R 2 = 0.89, with a = 0.4, b = −3.77,and c = 2.08, Fig. 5.

D. Action Selection Networks
Finally, the analysis was extended to a more specific and biologically plausible model of STN-GPe connectivity which can support action selection.Synchrony varied with channel overlap and average number of connections, forming a circular quadrant of low synchrony in the k = σ parameter space.The results further confirm the strong effect of the number and arrangement of connections in a topology on the level of synchrony within the network.
Algebraic connectivity performed well in capturing the network changes which contributed to variations in synchrony within the action selection network, Fig. 6(b).A step like transition to the synchronous state was observed, Fig. 5(e), Fig. 6(b), which was sharper than seen in the other networks examined.In the action selection network, the asynchronous state was composed of alternating bands of active and inactive channels, Fig. 6(c).In Figure 6(d) presents an action selection network structure just outside the bounds of synchrony where segregated action channels could not be maintained, and beta oscillations emerged.

IV. DISCUSSION
In this work, we have shown that the number and arrangement of connections in the STN-GPe circuit plays a crucial role in determining its level of intrinsically generated synchrony.To isolate the effects of network topology on synchrony, changes in net synaptic drive were controlled for by establishing a functional relationship between connection strength and post-synaptic degree.Synchrony generally increased as the number of connections increased.However, the nature of this relationship was highly dependent on the network generation rule.The pattern of connections was shown to be a key determinant of synchrony, as networks with the same number of connections could display asynchronous or synchronous activity depending on their topology, Fig. 5.This relationship was well approximated within and across different network topologies using the algebraic connectivity λ 2 .Synchrony could also be induced in action selection topologies proposed for the basal ganglia, through network changes corresponding to increases in λ 2 .

A. Scaling of Synaptic Connection Weights
The weighting functions identified for the STN-GPe and GPe-STN were found to be monotonically decreasing functions of the number of connections, while the inter-GPe scaling displayed a more complex relationship, Fig. 3. Scaling the weights enabled changes in synchrony to be attributed to structural changes rather than variations in mean firing rate or other internal neuron properties such as ion channel dynamics.These results show that synchrony can be achieved at realistic firing rates across a large range of degrees.While constraining mean firing rates and scaling inputs may seem artificial, similar synaptic scaling mechanisms can be observed in various brain regions [52], [53], including the STN [15], [16].It has been suggested that such mechanisms may exacerbate the motor symptoms of PD [16].

B. Influence of Network Size and Number of Connections
Network size had an impact on the level of synchrony for only two of the topologies investigated, small world and stochastic blocks, Fig. 4. The small world models were assigned a low rewiring parameter β = 0.01 and behaved similar to ring networks, with distances growing in proportion to network size [19].This may be observed in Fig. 4a, where the onset of synchrony occurs at a higher degree value for large networks.Stochastic block networks with small inter-group connection probability P e experienced a qualitative change in their degree-synchrony relationship at values of degree that depended on network size k = n m=4 , Fig. 4(a).Below this threshold, most new connections were added between neurons in the same block which led to local synchrony, Fig. 4(e), Authorized licensed use limited to the terms of the applicable license agreement with IEEE.Restrictions apply.
while above the threshold new connections occurred between different blocks, facilitating the global synchrony.
Motter et al. demonstrated that synchronizability depends only on the mean degree when using input normalization with large sufficiently random networks [48].We similarly observe a monotonic relationship between average degree and synchrony for the k-regular random, spatial, and preferential attachment topologies, Fig. 4. For the remaining topologies the assumption of "sufficiently random" breaks down and the pattern of neuron-to-neuron connections becomes important, even for relatively large networks comprising up to 2,000 neurons.This assumption of large, sufficiently random networks is core also to mean field models of the STN-GPe circuit which are widely used to investigate the emergence of pathological oscillations in Parkinson's disease [54], [55].These models absorb connection strength, number, and pattern into a single variable, and in doing so make implicit assumptions about network topology.As the arrangement of connections in the STN and GPe nuclei is unknown, uniform randomness may be a conservative assumption that introduces minimum bias.However such uniform randomness is not observed in the available connectomes [56] and the resulting networks display a high level of intrinsic synchronizability [57].This assumption may therefore limit the generalizability of results from such models.

C. Influence of Connection Pattern
The results presented demonstrate the critical role that connection pattern plays in the synchronizability of the STN-GPe circuit.Previous studies have relied on changes in external synaptic drive [29], [58], changes in connection strengths between populations [59], [60], or both [9], [61], [62] to explain the emergence of synchronous oscillations within the circuit.Connectivity is often implemented in these studies using one of two common patterns.The first entails connecting each post-synaptic neuron with a fixed number of presynaptic efferents chosen at random for each of the valid SG, GS, and GG synapse groups [31], [63], [64].This specification is useful for maintaining a fixed level of synaptic current in each neuron but comes with the same problems as uniform random networks discussed in the previous section.The second pattern covers connectivity structures that emulate the proposed action selection role of the basal ganglia, and stem from the Gurney, Prescott, and Redgrave model [29], [40], [59].However, it should be noted that the precise arrangement of neuron to neuron connections in the basal ganglia remains to be determined.It is possible to induce synchronous parkinsonian activity in both connectivity models by increasing coupling strengths or changing the activity of external inputs [64].However, the relationship between these changes and the resulting level of synchrony will depend on the underlying network structure and may not generalize.A main aim of this study was to uncover an underlying relationship between the arrangement of connections in arbitrary networks and the level of synchrony present.It was important to ensure that this relationship could be maintained across networks with distinct features.To examine this, five different network topologies were selected for investigation.The selected networks incorporated distinct features including the heterogenous degree distributions of preferential attachment networks [44], the block structure of stochastic block networks [65], and the small world structure of the Watts-Strogatz networks [41], the two dimensional locality of the spatial networks, and the on-center off-surround structure of the action selection networks [40].Algebraic connectivity, derived from the MSF, was found to capture this structure synchrony relationship for the networks considered, Fig. 5.For the MSF to apply to a networked system, the dynamics of the nodes must be identical and the coupling must be diffusive [28].Under these conditions, the first non-zero eigenvalue will determine the stability of the full synchronized state in which all nodes have the same dynamical trajectory.Issues arise in applying this framework to neural systems as neurons are rarely identical, and chemical synapses are not a two-way diffusive process.Coupling with chemical synapses exacerbates the problem of non-identical nodes as the summation of synaptic currents from different numbers of sources alter neural dynamics by inducing bursts in addition to changing firing rates, which further intertwines structural and dynamical properties.However, despite the limitations in applying the MSF to neuronal systems, the level of synchrony in the STN-GPe circuit correlated with algebraic connectivity for all four topologies, Fig. 5 This highlights the utility of algebraic connectivity in capturing the structure-synchrony relationship.A better understanding of this relationship could improve on interventions which aim to achieve desynchronization via targeted rewiring of the underlying network structure [66].
The model as examined in this paper did not include excitatory NMDA receptors.The activation of NMDA receptors can play a role in enhancing neural synchrony [67].To examine whether the observed relationships still hold in the presence of NMDA receptors, simulations were repeated with the addition of NMDA receptors incorporated using an NMDA synapse model [68] (see Supplementary Material S4).From the perspective of the MSF, increased excitability from NMDA receptor activation would alter the intrinsic synchronizability of the neurons, but not affect the synchronizability brought about by network properties.Consistent with this expectation, it was observed that the presence of cortico-STN NMDA receptors did not fundamentally alter the structure synchrony relationships observed in this study Fig S3 , S4.

D. Action Selection Networks
The literature on the action selection role of the basal ganglia, suggests that neurons are connected in discrete processing units corresponding to particular actions [40].Connectivity changes have been observed in PD for STN to GPe [18], GPe to STN [14], [15] and CTX to STN [16] connections, as well as within the striatum [69].It is possible that these changes could induce synchrony through a desegregation of the action selection circuits in a manner that increases the algebraic connectivity.The connection between circuit rewiring and aberrant motor symptoms has be suggested previously [70], [71], [72], and the link between beta oscillations and motor symptoms is well-established.Our results show how desegregation of action selection circuits as a result of changes in topology can lead to synchrony, demonstrating a possible path between observed changes in physiological connectivity and pathological state.

E. Limitations
While the results demonstrate the importance of the structure-synchrony relationship in the STN-GPe circuit, there are several limitations to be considered.The quantitative relationship between synchrony and algebraic connectivity λ 2 will depend on the neuron model, so the numerical values presented here only apply to this model.Despite the topologies being "stochastic", their randomness is welldefined, their degree distributions are unimodal, and a single parameter governs the network diversity.In addition, the network generation algorithms used here may not exhaustively or uniformly sample the space of possible graphs for a given class of topology.For all but the action selection topologies, the same pattern was used to generate all connection types (STN to GPe, GPe to STN, GPe to GPe), confining the networks to narrower regions in the sample space.The STN and GPe comprised of equally sized neural populations in order to isolate the fundamental dynamics of the system in a simplified network.However, in rodents, GPe neurons are three times more numerous than STN neurons [73].We confirmed that these differences in size did not fundamentally alter the results for Watts-Strogatz and spatial networks (supplementary materials S3).Finally, the synchrony measure captures globally correlated variations in both spiking and subthreshold oscillations but cannot capture local correlations, most evident in stochastic block networks.

V. CONCLUSION
The relationship between network structure and synchronous beta frequency oscillatory activity has been explored in a model of the STN-GPe circuit.Connections between neurons were weighted as a function of the in-degree of the receiving neuron, to isolate the effects of network topology.Subject to this scaling, we have shown that in general, synchrony increases as the number of connections increases, but the nature of this relationship depends on the topology governing the arrangement of the connections.Within a given topology, synchrony can also be tuned by varying the arrangement of connections while holding the average number constant.Algebraic connectivity correlated with synchrony across the topologies considered, providing a means of assessing the synchronizability of arbitrary connection arrangements.A critical role for network topology in determining the synchronizability of the network was also demonstrated in biologically plausible action selection network models.The results have important implications in the design of computational models of neuronal networks in the parkinsonian basal ganglia, as they provide insight into how the choice of network structure can impact the model's propensity for synchronization.They also provide a means of understanding how structural changes may contribute to pathological synchrony observed in the parkinsonian STN-GPe through an understanding of the network connectivity.

Fig. 1 .
Fig. 1.(a) Schematic of the STN-GPe circuit.(B) Example of a connectivity rule for specifying connections between pre and post synaptic populations.(C) The networks for each connection subset are combined, specifying the full STN-GPe circuit.

Fig. 2 .
Fig. 2. Probability distributions for each connection type (STN to GPe, GPe to STN, GPe to GPe) in the action selection topology.The distance between the center of the channel at index 0, and the center of lateral inhibition due to GPe to STN connections is d = 40 neurons.(a,b) Connectivity distributions for channel spread σ = 1, σ = 2 respectively.(c) Schematic of an action selection network with connections going to and from a specific GPe neuron.

Fig. 3 .
Fig. 3. Synaptic weights identified to maintain constant mean firing rates as the number of connections per neuron varied.Each point corresponds to the results from the simulation that minimized the cost function for the given degree.(a) STN firing rates, (b) GPe firing rates, (c) synchrony, (d) STN to GPe connections, (e) GPe to STN connections, and (f) GPe to GPe connections as functions of in degree.Simulations were not run for networks with unitary average degree so F(k=1)=F(2) for all connection types.

Fig. 4 .
Fig. 4. (a) Degree -synchrony relationship for each topology, with networks comprised of 200, 1000 and 2000 neurons.The median value of synchrony across 10 simulations is presented, error bars denote the interquartile range.Raster plots and mean membrane potentials for representative networks are shown in (b),(c),(d) and (e).(b) A typical low synchrony state for a spatial network with k =5.(c) A synchronous state with beta oscillations exemplified in a spatial network with k =30.(d) The asynchronous state for Watts-Strogatz networks with low rewiring is unique in having spatially varying waves of activity which average over time.(e) Stochastic block networks show correlated activity within their groups which averages due to the phase and frequency differences between blocks.

Fig. 5 .
Fig. 5.The effect of network wiring pattern on synchrony for different topologies with fixed network size n=1,000, and average degree k =10.Relationship between synchrony and (a) rewiring probability β for Watts-Strogatz networks, (b) distance dependence γ for spatial networks, (c) attachment inequality σ for preferential attachment networks, and (d) inter-block connectivity p e for stochastic block network ε (e) Synchrony χ vs algebraic connectivity.

Fig. 6 .
Fig. 6.Action selection networks.(a) Synchrony versus mean degree k and channel spread σ.(b) Synchrony vs algebraic connectivity.(c) Channel activity emerges in the network with k =10 and σ = 1.(d) The channel activity cannot be sustained when channel spread is increased, σ = 2, and oscillations emerge as channels desegregate.