Coupled beta and high-frequency oscillations emerge from synchronized bursting in a minimal model of the parkinsonian subthalamic nucleus

Local field potentials (LFPs) recorded from the subthalamic nucleus (STN) in Parkinson's disease (PD) exhibit a distinctive multiscale spectral signature: exaggerated beta-band oscillations (13-30 Hz) coupled to high-frequency oscillations (HFOs, 200-400 Hz), with HFO amplitude being phase-locked to the beta cycle. This phase-amplitude coupling (PAC) has been identified as a promising biomarker of the parkinsonian state, yet no biophysical model has explained how it emerges, what determines the HFO frequency, or how HFOs can exist without beta modulation in the medicated STN. Here we show that a heterogeneous population of excitatory Izhikevich neurons with recurrent coupling produces three dynamical regimes: (i) asynchronous tonic firing, (ii) asynchronous bursting, in which neurons burst individually producing broadband HFO power but without coherent population-level PAC, and (iii) synchronous bursting, which gives rise to beta-HFO PAC. The regimes are governed by two biophysically interpretable parameters that capture complementary effects of dopamine depletion: one reflecting changes in intrinsic neuronal excitability, the other reflecting changes in synaptic coupling strength. The transition from asynchronous to synchronous bursting in this model emulates the emergence of pathological STN neuronal activity in the parkinsonian state. HFO peak frequency varies continuously across the two-parameter landscape, suggesting a possible mechanism of the clinically observed shift from slow (200-300 Hz) to fast (300-400 Hz) HFOs between medication states. The character of the synchronization transition depends on baseline excitability, ranging from a sharp co-emergence of bursting and synchrony at low excitability to a decoupled transition at intermediate excitability, where the bursting fraction saturates while the population synchronization continues to increase with coupling. The asynchronous-to-synchronous bursting transition and the beta-HFO PAC in the STN population also arise in an extended version of the model with reciprocal coupling to an inhibitory population of globus pallidus externa (GPe)-like neurons, with the beta frequency varying with the delay in the synaptic coupling between the two populations. The model generates testable predictions for future clinical and experimental studies, provides a numerical dissection of how mesoscopic LFP features could map onto microscopic neuronal dynamics, and serves as a computational building block for future circuit-level models that may inform brain stimulation strategies tailored to the patient-specific dynamical state of the STN.