Mechanistic models of α-synuclein homeostasis for Parkinson's disease: A blueprint for therapeutic intervention

1. Introduction

Parkinson's disease (PD) is a debilitating neurodegenerative disorder with currently no disease-modifying therapy or cure [1]. The high prevalence of PD on the global scale and its exponential growth over the last few decades suggest the imminent outbreak of a non-infectious pandemic propelled by aging and industrialization as unconventional vectors of this rapidly growing neurological disorder [2]. Given the enormous strain that a pandemic would put on the healthcare system and the resulting human and economic toll, it is no surprise that investigating PD pathogenesis has rapidly become a serious public health concern. Since 1817, when the “shaking palsy” was first described [3], considerable progress has been made in understanding PD, from the clinical symptoms to the molecular mechanisms. Such progress has yet to translate into effective treatments that can hamper or even reverse PD progression. The current gold-standard therapy (Levodopa or L-Dopa treatment) manages motor symptoms but cannot stop neurodegeneration [4]. It counterbalances dopamine (DA) deficiency in the basal ganglia by relying on the pharmacological input of the DA precursor L-Dopa. Indeed, impairments in the dopaminergic pathway and the subsequent neuronal death in the substantia nigra pars compacta stand out as a primary hallmark of the disease [1].

Human and molecular genetic research has soon uncovered PD's multifactorial nature, thus unveiling a complex molecular network [1, 5]. In addition to impairments in DA metabolism, numerous processes such as protein aggregation, defective degradation, neuroinflammation, and oxidative stress contribute to disease onset and progression through various feedback mechanisms (see Fujita et al. [6] for a comprehensive disease map). The presynaptic amyloidogenic protein α-synuclein (αsyn) emerges from the intricate molecular landscape as a key driver of PD pathogenesis [1]. Here, we focus on two mechanisms governing its homeostasis, i.e., the aggregation and degradation pathways. Figures 1A, 2A display the biological mechanisms involved in these pathways. Their link to PD neurodegeneration is supported by multiple experimental and genome-wide association studies [5, 7, 8]. Even though multiple studies have highlighted the significant contribution of these processes to PD neurodegeneration [9–11], the interplay between the uncontrolled accumulation of αsyn aggregates and impairments in the degradation mechanisms is still unclear [12].

Given the intricate biological picture outlined so far, identifying the specific role of αsyn aggregation and degradation pathways would support therapeutic strategies currently under investigation [13], e.g., small molecule inhibition, αsyn antibodies, and autophagy stimulation, or even pave the way for new pharmacological approaches. To this end, the experimental analysis calls for a quantitative systems pharmacology (QSP) approach at the interface between pharmacological research and systems biology [14]. As recently highlighted in Geerts et al. [15], Abrams et al. [16], Bloomingdale et al. [17], and Bloomingdale et al. [18], QSP models indeed hold great potential in neuroscience as they provide a systems-level understanding of complex biological processes and insights into their response to drugs. Specifically, modeling efforts in PD research are moving toward multiscale QSP models that describe molecular, cellular, whole-brain, and organism levels to assess how the effects of a molecular perturbation can scale up to influence the clinical outcome [19].

Starting from the molecular level, various mechanistic models falling within the systems biology framework have tackled PD complexity. As highlighted by previous studies on systems biology models for neurodegenerative disorders [20, 21], they have mainly focused on protein aggregation, synaptic transmission, apoptosis, oxidative stress, and genetic components rather than protein degradation. The relatively small number of models that include degradation processes is likely to rise along with increasing evidence that supports the critical role of the protein clearance machinery in PD [10, 22]. With this in mind, we collect the advances in mathematical modeling of αsyn aggregation and degradation by the ubiquitin-proteasome pathway, macroautophagy, and chaperone-mediated autophagy. As shown in Figures 1, 2, we divide the mechanistic models into two groups: (i) Single-pathway models of αsyn aggregation representing the chemical kinetics of the process (i.e., chemical kinetic models); (ii) Multiple-pathway models of αsyn homeostasis focusing on the interplay between αsyn aggregation and degradation mechanisms (i.e., degradation models). We discuss their potential improvements and possible implications in PD therapy development.

2. Mechanistic models of αsyn aggregation kinetics

Protein aggregation and misfolding are central features in the molecular network of neurodegeneration [23]. It is therefore no surprise that a large part of the modeling effort targeting PD investigates the chain of microscopic events in αsyn aggregation. These studies fall into a broader line of research devoted to uncovering the generality of amyloid formation by focusing on the chemical and physical aspects of this phenomenon [24, 25]. Here, we analyze mechanistic models that extend the formalism of chemical kinetics to the aggregation process; we thus refer to them as chemical kinetic models (see Figure 1B). These models adopt a master equation approach building on the pioneering work by Oosawa and Kasai [26] to provide a mathematical formulation of each molecular step in the process and derive the corresponding rate constants. To gain novel biological insights and formulate new scientific queries, in silico models evolve along with increasingly advanced experimental procedures in a virtuous cycle of model calibration and validation against in vitro quantitative measurements [23].

2.1. αsyn misfolding and aggregation and PD neurodegeneration

Soluble αsyn is an aggregation-prone protein with an intrinsically disordered conformation that prevents its spontaneous self-assembly [27]. Many factors can trigger αsyn misfolding and accumulation, such as high αsyn concentrations, genetic mutations, increased reactive oxygen species (ROS) levels, and interaction with lipid membranes [28]. The process initiates upon structural changes that expose the highly amyloidogenic region of the protein to the cytosol. Misfolded monomers slowly aggregate into oligomeric species, thus creating nucleation seeds. In return, seeding-competent oligomers can rapidly accumulate and elongate by monomer addition forming protofibrils and insoluble fibrils (Figure 1A). Multiple αsyn species eventually gather into amyloid inclusions known as Lewy bodies (LBs) and Lewy neurites (LNs), which represent one of PD's main molecular hallmarks [29].

A considerable body of research points to the harmful effect of αsyn aggregation and its crucial contribution to PD neurodegeneration. For instance, the SNCA gene, encoding for αsyn, was the first gene associated with PD [30]. Still, there is no agreement on the specific mechanisms of toxicity. Seeding and spreading ability and inherent toxicity are the functional and structural factors that determine the overall toxicity of aggregates [27, 31]. According to in vitro and in vivo experimental evidence, the most cytotoxic species are oligomers with a partially folded structure similar to mature fibrils [28, 31, 32]. Due to their conformation, some oligomeric species can create pores within the cell membrane, thus enhancing its permeability and dysregulating calcium homeostasis [33]. They can also impair mitochondrial function and proteolytic machinery. On the other hand, small fibril fragments have been identified as major-spreading agents in PD because of their high seeding efficiency and internalization ability [28]. According to one of the leading hypotheses for αsyn propagation, pathogenic species are assumed to spread in a prion-like manner [34]; they are released in the extracellular space and then endocytosed by neighboring neurons, where they recruit non-native monomers acting as seeds for further aggregation. As a result, protein aggregates propagate across interconnected brain areas over decades, eventually resulting in PD neurodegeneration [1].

Considerable progress has been made in understanding the kinetic of αsyn misfolding and aggregation, however, there is no agreement on many aspects of this mechanism. Chemical kinetic models can provide valuable insights into understanding, for example, the relative contribution of different microscopic events to the process, the structural and functional features of αsyn aggregate species [32], and the relationship between αsyn aggregation and LB formation and maturation [35, 36].

2.2. The master equation approach

The aggregation pathway of a generic amyloidogenic protein can be represented by a minimal set of fundamental reactions, regulating the transition from soluble monomers to fibrils through various intermediate states. As shown in Figure 1A, we can differentiate between primary reactions, i.e., homogeneous primary nucleation and elongation, and secondary reactions, i.e., heterogeneous secondary nucleation and fibril fragmentation [37, 38]. Homogeneous nucleation initiates the aggregation process by combining monomers into newly formed aggregates that, in turn, elongate by monomer addition. By integrating these reactions, fibril fragmentation and heterogeneous nucleation on aggregate surfaces multiply the aggregate number in a positive feedback loop with elongation, thus ensuring self-replication. Their relative contribution to the aggregation process depends on the specific protein considered [25]. Primary and secondary reactions can incorporate more complex events, such as interaction with lipid membranes and multistep nucleation. The master equation approach [26, 39] translates the reaction network into an infinite system of coupled non-linear differential equations describing, in probabilistic terms, the time course of the species involved, from monomers to filaments identified by their monomer count.

The main goal of this approach is to find an analytical solution for the system to describe its reaction time course. As a result, the analytical solution can be fitted to experimental data to identify kinetic rate constants and provide mechanistic insights into the aggregation process. Since the master equation formulation is often analytically and numerically intractable, specific approximation methods can be employed to reduce the system to two differential equations describing the time evolution of aggregate number and mass concentrations, i.e., readily available observables. Multiple efforts have channeled into a formal analysis of this minimal closed system. Knowles and colleagues [25, 37, 40] proposed fixed-point approaches to derive explicit solutions for the entire time course of aggregates. All the identified parameters had a clear physical meaning related to the underlying molecular events. Furthermore, analytical solutions exhibited scaling laws used in model fitting to determine the importance of primary and secondary reactions [41]. For instance, in the case of breakable fibrils [25], the authors found that the reaction time course mainly depended on a parameter combining elongation and fragmentation rate constants. By analyzing the relationship between this parameter and the observables, such as half-time and maximal growth rate, the authors showed that even small perturbations in the fragmentation reaction significantly affected the aggregation process. This procedure was implemented in a global fitting tool for chemical kinetic models of in vitro protein aggregation (AmyloFit) [38].

Model parsimony and data accuracy are essential requirements for the fitting protocol. A parsimonious model can be fitted to a wide range of experimental data acquired for different initial monomer concentrations and in the presence and absence of preformed aggregate seeds. Techniques for measuring aggregate species vary from fluorescent dyes (e.g., thioflavin T) to mass spectrometry for aggregate mass concentration and single-molecule fluorescence for oligomer concentrations [38, 42, 43]. With improvements in experimental approaches come recent advances in the mathematical framework, including the investigation of stochastic aggregation kinetics [44], the spatial dependence of protein aggregation [45], and the nature of oligomeric species [32, 46].

2.3. Chemical kinetic models of αsyn aggregation

A significant number of chemical kinetic models address αsyn aggregation and its implications in PD.

2.3.1. Nucleation-polymerization models

NP models describe aggregation as fibril formation by monomer addition, which follows a nucleation event and escalates by secondary reactions. This process is represented by a characteristic sigmoidal growth curve, identified by a lag phase, a growth phase, and a plateau [25]. By relying on the kinetic theory of aggregation that revolves around the master equation approach, multiple studies have combined experimental analysis with these mechanistic models to investigate underlying microscopic events and the effect of varying experimental conditions.

Buell et al. [47] focused on the effect of pH on αsyn aggregation in the presence and absence of preformed fibril seeds. Homogeneous and heterogeneous nucleations and fragmentation were negligible in purely monomeric solutions at neutral pH. Adding preexisting aggregates to the solution enhanced the aggregation rate and shortened the lag phase of the sigmoidal growth curve. To obtain the elongation rate constant, the authors used a minimal model involving only an elongation reaction saturating at high monomer concentrations (i.e., two-step elongation). On the other hand, a mildly acidic pH uncovered a much faster aggregation, that is, a longer lag phase and a more rapid growth phase, due to either surface-catalyzed nucleation or fibril fragmentation. Building on Buell et al. [47] and Gaspar et al. [48] determined the pH-dependent role of these two events by relying on the theoretical analysis of a model that includes primary and secondary reactions. According to an experimental analysis supported by a global fitting procedure [38], the autocatalytic nature of aggregation was mainly determined by secondary nucleations of monomers on fibril surfaces. These studies provided insightful information on physiological mechanisms; for example, they showed that even a slight pH change associated with different intracellular environments might enhance the production of toxic aggregates.

Small lipid vesicles can alter αsyn aggregation as much as pH levels [49]. Upon interaction with lipid membranes, the protein adopts a stable α-helical conformation [50]. Galvagnion et al. [49] showed that an imbalance between free and membrane-bound αsyn states can trigger amyloid formation by enhancing heterogeneous primary nucleation on lipid surfaces. To confirm and quantify these results, the authors employed an NP model including this nucleation reaction, two-step elongation, and the conversion between intermediate species [42]. Relying on Buell et al. [47] and Galvagnion et al. [49], Flagmeier et al. [51] investigated how mutations in αsyn influence the fate of amyloid fibrillation. The authors performed three experiments to separately analyze initiation, elongation, and multiplication steps for all αsyn variants. Each experiment was associated with a mechanistic model: the NP model from Galvagnion et al. [49] for lipid-induced aggregation; and linear polymerization models for strongly and weakly seeded aggregation processes at neutral and mildly acidic pH. Experimental and theoretical results agreed that the rate of surface-catalyzed nucleation on lipid membranes and fibril surfaces could vary by several orders of magnitude across different mutations.

2.3.2. Nucleation-conversion-polymerization models

The NP models are coarse-grained to the extent that they do not explicitly account for αsyn intermediate species. The modeling choice mainly relates to the absence of accurate measurements of oligomer concentrations in the aggregation process. Oligomers are indeed transient and highly heterogeneous species. The development of new methods to report their levels has motivated the transition to NCP models [46, 52]. These models include the structural interconversion from disordered clusters to oligomers with a partially formed fibrillar structure, eventually growing into fibrils.

Multiple studies on αsyn aggregation have provided experimental and theoretical support to this updated kinetic theory [53–56]. Single-molecule Förster Resonance Energy Transfer (FRET) experiments can distinguish two oligomeric populations associated with low and high FRET levels. Type-A oligomers are early-formed, disordered, highly degradable, and mildly toxic species, whereas type-B oligomers are later-forming, more compact, degradation-resistant, and ROS-promoting species [42, 55]. Both oligomeric species appear as required subsequent steps in the aggregation process [57]. As displayed in Figure 1A, NCP models tailored for αsyn aggregation include primary nucleation, which generates type-A oligomers from monomers, and the structural conversion from type-A to type-B oligomers and then fibrils [42, 55]. Both fibrillar and oligomeric species increase in size by monomer-dependent elongation. In this reaction network, the conformational change between oligomeric species emerges as a crucial step for initiating αsyn aggregation. Upon fitting the model to experimental data, Cremades et al. [42] showed that the corresponding kinetic rate is slow, thus providing a sufficient lag time for the protein quality control machinery to hamper the pathological process. A potential therapeutic strategy may employ molecular chaperones or small molecules to regulate this reaction. A follow-up study [55] used the NCP model to quantify seeding ability and conditions for type-B oligomers and fibrils, with implications on prion-like spreading. Model predictions showed that fibrils were more prone to seed aggregation, whereas oligomers mainly contributed to the increase of oxidative stress.

The studies reported above focused on the early stages of αsyn aggregation; thus, they did not require full-time analytical solutions. To extend the analysis to the entire time course, Dear et al. [46] provided a general mathematical framework for NCP models. The authors proposed a reaction network with a single oligomeric state. Such a coarse-grained model is suitable to avoid overfitting in the analysis of experimental data on total oligomer concentration that, as such, does not account for differences between oligomeric populations, e.g., see Zurlo et al. [56]. Moreover, the model is instrumental in classifying oligomers according to theoretically defined metrics: persistence or half-life, abundance, and productivity, i.e., the propensity to form fibrils rather than dissociate. The authors proved that oligomers were more prone to dissociation than fibrillation as a universal feature. Results focusing on αsyn oligomers showed that they were very kinetically stable, more persistent than monomers, relatively abundant, and highly productive compared to other proteins.

Overall, chemical kinetic models contribute to a detailed mechanistic understanding of αsyn aggregation by combining theoretical results and in vitro data acquired through advanced techniques. By exploring a wide range of in silico experimental conditions, these models return coherent predictions that can be useful for studying PD. Further investigation in this direction should focus on whether and how these results apply to in vivo settings to provide insights into the mechanisms of neurodegeneration and, possibly, inform therapeutic strategy design [58].

3. Mechanistic models of αsyn aggregation and degradation pathways

Given the central role of degradation in regulating αsyn homeostasis, multiple models have investigated this process in relation to PD neurodegeneration, including models solely focused on αsyn aggregation and degradation or miscellaneous models that include additional pathogenic processes (e.g., ROS production and defective dopamine metabolism). These degradation models offer a heterogeneous picture in terms of specific mechanisms they consider, granularity, and modeling approaches [21] (see Figure 2B). The biological processes can be represented by varying degrees of detail, and the modeling formalism can range from ordinary differential equations (ODEs) and algebraic equations to stochastic simulations, partial differential equations, or flux balance analysis [20, 21].

3.1. Defective protein degradation machinery and PD neurodegeneration

The proteolytic machinery combined with the molecular chaperone system ensures the renewal of functional components and regulates the fate of altered proteins to avoid the build-up of intracellular damage. In PD and other neurodegenerative diseases, degradation efficiency declines with age and oxidative stress [59], thus increasing misfolded αsyn and aggregate levels. On the other hand, αsyn aggregation can hamper the proteolytic machinery causing the accumulation of other cellular components [60]. The interplay of αsyn aggregation and impairment of degradation in multiple negative feedbacks poses a causality dilemma [61]: which pathological process comes first? There is still no clear answer. Degradation defects and αsyn aggregation can be either primary or secondary events of neurodegeneration [12, 62]. The investigation of their role in PD is even more challenging due to interfering factors such as increased oxidative stress and mitochondrial dysfunctions [1, 59]. Systems biology models can help elucidate these contributions by providing in silico benchmarks to assess these hypotheses.

3.1.1. The ubiquitin-proteasome pathway

UPP is the dominant clearance route for mislocated, misfolded, damaged, and short-lived proteins. Target substrates are selectively tagged by polyubiquitin chains through multiple enzymatic reactions and then transferred to the proteasome. Ubiquitinating and deubiquitinating enzymes such as PARKIN and UCH-L1 play a central role in modulating protein degradation in PD, as indicated by UPP-related familial cases [8, 60]. For instance, UCH-L1 maintains a stable pool of monomeric ubiquitin; an imbalance of this pool can lead to poor substrate labeling with consequent impaired proteasomal degradation. Furthermore, the presence of LBs as a PD molecular hallmark corroborates the correlation between proteasomal dysfunctions and PD pathogenesis. Indeed, LBs mainly consist of ubiquitin-tagged αsyn aggregate, ubiquitinated proteins, and UPP components [60, 63, 64]. According to multiple studies [11, 65], their composition suggests a potential cytoprotective role for LBs, which may act as sinks for poorly degraded proteins to reduce the toxic effect of protein accumulation and slow down neurodegeneration. Still, it is unclear whether these cytoplasmic inclusions represent a harmful or defense molecular mechanism for the neuron. UPP dysfunctions and αsyn aggregation create a negative feedback loop. Genetic and age-related impairments in UPP may hamper the proper degradation of misfolded aggregate-prone αsyn. On the other hand, αsyn species refractory to proteasomal degradation, such as mutant αsyn forms, oligomers, and fibrils, can inhibit proteasomal activity by sequestering the proteasome. In addition, the interaction with oligomers and oxidatively damaged αsyn may directly alter proteasome subunits, as shown in the green box in Figure 2A.

3.1.2. The autophagy-lysosome pathway

ALP removes long-lived proteins, lipids, and dysfunctional organelles, such as aging mitochondria (i.e., mitophagy), thus maintaining neuronal homeostasis. It includes three degradation mechanisms: macroautophagy, chaperone-mediated autophagy (CMA), and microautophagy. Experimental findings indicate that both UPP and ALP degradation of αsyn can occur [66]. The autophagic mechanisms are involved in physiological αsyn turnover and αsyn aggregates degradation (see the purple box in Figure 2A). Specifically, wild-type αsyn monomers are degraded mainly by CMA [67], which relies on molecular chaperones to recognize a specific motif on the target substrate. Upon chaperone binding, the cargo interacts with the membrane receptor LAMP2A, which assembles into multimeric translocation complexes. The substrate then enters the lysosome lumen, where it is rapidly degraded. CMA also clears ALP and UPP-related proteins, such as LRRK2, UCH-L1, DJ-1, and GCase [68]. Furthermore, macroautophagy represents the main clearance route for αsyn aggregates due to the inaccessibility of the CMA degradation motif [62]. This degradation mechanism consists in the formation of autophagosomes, i.e., double-membrane vesicles surrounding the cargo and fusing with lysosomes to allow hydrolase-mediated degradation in the lysosomal lumen. In selective macroautophagy, specific autophagosome membrane receptors target protein aggregates. The link between impaired autophagy and PD neurodegeneration conceals a complex interaction network between αsyn aggregates and the lysosomal proteolytic machinery. Growing evidence supports the negative effect of αsyn aggregation on autophagy acting on multiple levels [62, 69], as displayed in Figure 2A. At the same time, defective ALP cannot handle harmful αsyn species, thus leading to their accumulation. Such defects can originate from various sources, such as increased oxidative stress, aging, and dysfunctional proteins related to genetic PD forms, namely, LRRK2, DJ-1, and GBA [5], and with a physiological role in ALP.

3.2. Degradation models

3.2.1. UPP models

Relying on experimental evidence pointing to the essential contribution of UPP to neurodegeneration [70, 71], Sneppen et al. [72] provided a minimal theoretical model of the interplay between impairments in proteasomal degradation and αsyn fibrillation. The proteasome ensured protein homeostasis by keeping the amount of αsyn aggregates low. On the other hand, insoluble αsyn fibrils sequestered this proteolytic machinery and formed a slowly-degraded complex, thus reducing proteasome function. Computational simulations of this ODE system supported by stability analysis showed a bifurcation behavior. When the fibril level exceeded a particular threshold related to proteasome production capacity, the double-negative feedback broke αsyn homeostasis in favor of a pulsatile oscillatory regime characterized by spikes and long periods of low proteasome concentrations. As a result, the system alternated between the hypothetical recovery of proteolytic activity and the uncontrolled accumulation of harmful aggregates and fibrils. Several oscillations eventually built up toxicity over decades, thus explaining the slow progression of PD neurodegeneration. The authors therefore associated proteasome overload and the subsequent oscillatory regime with the onset of sporadic PD phenotype. In a follow-up study [73], the stochastic analysis of the core negative feedback confirmed and generalized the above phenomenon for the whole protein clearance machinery. The bifurcation behavior characterized by a healthy homeostatic state and an oscillatory disease state indeed emerged for both UPP and ALP degradation. In addition, stochastic simulations captured dynamical differences between the two degradation systems, modeled as enzymatic reactions slowly reducing aggregate concentrations. Lysosomes were assumed to be more efficient but less abundant than proteasomes. Fibril-mediated sequestration of the proteasome therefore led to regular and almost deterministic oscillations, whereas lysosome overload resulted in irregular fluctuations.

In contrast to the high-level representation of the UPP in Sneppen et al. [72], Pigolotti et al. [73], and Proctor et al. [74, 75] proposed a detailed description of this pathway by including mono- and polyubiquitination, proteasome binding, ATP-dependent degradation, and proteasome inhibition by protein aggregates. These studies were devoted to the mathematical modeling of the aging process, which contributes to the decline of proteolytic activity and the build-up of oxidative stress due to flawed mitochondria [76, 77]. Both models indeed investigated the role of the proteasomal degradation pathway in age-related neurodegeneration. In Proctor et al. [74], Proctor and colleagues first presented a stochastic model for aggregation and proteasomal degradation of a generic protein under normal homeostasis and aging. The authors then tailored the system to assess the UPP role in αsyn aggregation during PD neurodegeneration [75]. The final model included additional modules representing αsyn metabolism and UCH-L1 deubiquitination activity in proteasomal degradation. This stepwise approach relied on the idea that the molecular landscape of neurodegeneration can be arranged in building blocks, as most biological processes are assumed to be modular [20, 78]. Both general and specific models proved to be functional for investigating the negative feedback of protein aggregate formation and UPP impairment. Employing the more general model [74], the authors explored the impact of various experimental procedures, e.g., shutting down a ubiquitin-protein ligase or the proteasome activity. The results suggested that depleted ubiquitin pools might contribute to proteostasis impairment. In addition, relying on the model tailored for PD [75], the authors analyzed the effects of both mutated and oxidatively damaged UCH-L1. Stochastic simulations indicated that proteasome inhibition might not be the triggering event for neurodegeneration but rather an exacerbating factor accelerating αsyn aggregation and inclusion formation. Both models were validated against data gathered in ad-hoc designed experiments, in a virtuous cycle of model refinement and in vitro experiments.

3.2.2. UPP and ROS models

Since oxidative stress and impairments in UPP degradation are interconnected, multiple studies have analyzed their synergistic effect on αsyn homeostasis. Cloutier et al. [79] provided a mathematical model accounting for the interplay between αsyn aggregation, ROS production, and proteasome-mediated degradation. Given the variety of molecular processes involved, the authors were able to evaluate the impact of multiple risk factors on αsyn metabolism. For instance, perturbations of system parameters related to oxidative stress and protein clearance capacity reproduced the high variability of LB concentration in PD forms [79, 80]. Moreover, supported by experimental verification [81], deterministic simulations of the model uncovered a bistable process whereby the system irreversibly switched from a healthy to a disease state, identified by high misfolded αsyn and ROS levels, in response to toxin exposure, genetic mutations, and aging. On the other hand, the model could not recapitulate the oscillatory regime due to decreased levels of available clearance machinery obtained in Sneppen et al. [72] and Pigolotti et al. [73]. To identify the primary mechanisms determining the bifurcation behavior, the authors reduced the system to focus only on the double-positive feedback between αsyn misfolding and ROS production [82]. This minimal version retained the bistable behavior, thus associating PD onset with a switch-like transition in response to enhanced ROS concentrations. In line with [72, 73], the central idea remained that neurodegeneration results from the transition between two dynamical behaviors.

Raichur et al. [83] also presented a model targeting the interplay between oxidative stress, αsyn aggregation, and UPP dysfunctions. Compared to Cloutier et al. [79], this system displayed a more granular description of the UPP pathway, from ubiquitination and proteasome recognition to degradation. The close molecular detail of ROS production and UPP degradation provided in the model inevitably hampered the theoretical analysis. On the other hand, the explicit representation of proteasome-mediated degradation enabled the authors to explore different scenarios of UPP-related genetic predisposition: mutations in UCH-L1, PARKIN, and DJ-1 genes, combined with age-related oxidative stress. These familial PD forms determined the accumulation of ubiquitinated proteins, leading to a more rapid αsyn aggregation than sporadic forms. Furthermore, the synergistic combination of ROS production and UPP-related genetic mutations fueled αsyn aggregation. These results pointed to therapeutic strategies involving restored UPP functionality or enhanced ALP activity in combination with antioxidant treatments.

3.2.3. UPP and ALP models

The modeling picture outlined so far focuses on the ubiquitin-proteasome pathway, while the autophagy-lysosome pathway is under-represented despite its critical role in αsyn turnover; when considered, autophagic clearance has often been incorporated as a pseudo-first-order reaction. This is not the case of Ouzounoglou et al. [84], who provided a detailed representation of selective macroautophagy, CMA, and proteasome-mediated degradation, focusing on the inhibitory feedback between αsyn oligomers and the CMA receptor LAMP2A. Upon model calibration and validation against αsyn overexpression data from human neuroblastoma cells, the authors performed stochastic simulations to explore different interventions with neuroprotective potentials, such as shutting off DA production, increasing LAMP2A levels, and reducing αsyn synthesis. The resulting predictions of increased cell viability were in agreement with experimental evidence and literature [85, 86].

To provide valuable insights into the role of αsyn degradation in neurodegeneration, a mechanistic model should integrate ALP and UPP degradation mechanisms, ideally accounting for the crosstalk between these pathways. In this direction, Sass et al. [87] proposed a model targeting asyn aggregation, UPP and ALP degradation, and DA metabolism and based on the biochemical system theory. A granular representation of these mechanisms and their interactions facilitated the analysis of specific system perturbations in both disease and treated states (e.g., enhanced αsyn fibrillation and aggregation rate). Including DA metabolism in the model also enabled the authors to simulate the impairment of DA vesicle packaging and available treatments related to its shortage, such as L-Dopa therapy combined with monoamine oxidase inhibitors. Increasing model complexity even further, Büchel et al. [88] presented an extended model of a dopaminergic neuron that includes the highest number of molecular components among the available PD models of αsyn homeostasis [20]. The model described 11 biological processes, including DA metabolism and transport, αsyn aggregation, oxidative stress, UPP and ALP clearance mechanisms, and specific details on PARKIN and DJ-1 proteins. Here, the modeling approach adopted was crucial to tackle model complexity. The authors employed flux balance analysis to investigate the steady-state behavior of the model rather than deterministic or stochastic analysis based on mass action kinetics. The qualitative prediction of multiple in silico experiments pointed to increased ROS, neurotoxin MPTP, and αsyn levels as crucial factors in PD neurodegeneration.

3.2.4. αsyn degradation and axonal transport models

The degradation models analyzed so far do not account for the spatial dependence of αsyn aggregation and degradation mechanisms. However, αsyn homeostasis and disease progression may strongly depend on intraneuronal protein propagation. Along this line, Kuznetsov and Kuznetsov [89, 90] investigated αsyn axonal transport and how its perturbation between the somatic and synaptic compartments affected the aggregation and degradation processes. By relying on a two-compartment ODE model [90], the authors predicted axonal death due to αsyn accumulation at the synaptic terminal, pointing to defective protein clearance machinery as the primary cause for αsyn aggregation onset compared to reduced αsyn production and transport [90]. Differences in the degradation of αsyn monomers and aggregates were determined by assuming distinct half-lives. Furthermore, building on recent findings on the absence of αsyn fibrils in LB inclusions [36], the authors proposed a mathematical model of αsyn aggregation, degradation, and transport to test the hypothetical autocatalytic nature of LB formation and its impact on protein accumulation [91]. Here, the authors assumed spatial-dependent proteasomal and autophagic degradation mechanisms. The study pointed to the central role of the autophagic clearance pathway in LB formation and PD progression; defective autophagy may favor the aggregation of membrane-bound organelles catalyzed by αsyn aggregates and, thus, LB formation. Most degradation models do not differentiate between αsyn fibrillation and LB formation, whereas chemical kinetic models do not even include these cytoplasmic inclusions. The detailed implementation of LB formation and maturation would allow for the investigation of LB neuroprotective role [65, 91].

The wide variety of degradation models reflects the multifactorial nature of PD. Such a modeling effort has provided multiple insights into the interplay between αsyn aggregation and protein degradation by spanning various levels of molecular detail, from coarse-grained to highly granular representations (see Figure 2B). Minimal models have been used to uncover dynamical features of disease progression (e.g., Sneppen et al. [72], Pigolotti et al. [73], and Cloutier et al. [82]). More detailed models (e.g., Proctor et al. [75] and Ouzounoglou et al. [84]) have proven to be suitable for testing the effects of potential treatments. Overall, degradation models move toward an integrative representation of PD pathogenesis. Future work may focus on the crosstalk between UPP and ALP pathways, the interplay between autophagy and apoptosis in PD neurodegeneration [92], and the integration of other control systems for protein homeostasis such as the molecular chaperone system [93].

4. Discussion

Mechanistic modeling grounded in molecular-level knowledge is a powerful tool for tackling PD biological complexity. By exploring the intricate interaction network of neurodegeneration through in silico experiments, mathematical models can investigate the key drivers of PD pathogenesis, discover promising targets, and even simulate the effect of potential single or combination therapies. Given the pivotal role of αsyn homeostasis in PD, this review has provided a collection of available mathematical models of αsyn aggregation and degradation pathways.

Chemical kinetic models focus on a single molecular pathway and rely on a unifying kinetic theory that translates into the consistency of modeling formalism (i.e., master equations approximated to ODE systems). In contrast, degradation models describe multiple biological processes and adopt a wide range of approaches to analyze the dynamical features of the system. The latter are modeled as deterministic and stochastic simulations based on mass action kinetics; yet, there are some exceptions, such as flux balance analysis, mainly related to the complexity of the reaction network considered. Stochastic approaches [94, 95] have been employed to analyze biological processes with low molecular counts and capture their inherent fluctuations. Future investigation may explore the role of intrinsic noise in PD onset and progression since it can unpredictably affect the system behavior, for example, by triggering the transition from a healthy to a disease state leading to PD onset [21, 79].

4.1. Integrating chemical kinetic and degradation models

Both chemical kinetic and degradation models include a representation of αsyn aggregation. As shown in Figures 1A, 2A (light-blue box), chemical kinetic models account for each molecular step separately, whereas degradation models employ macro-variables for different αsyn aggregate species according to their size and functional features. The master equation approach returns parameters that physically relate microscopic events to macroscopic observables. On the other hand, a higher-level representation of the process enables the targeting of additional pathways and the investigation of their interplay in the molecular network. A modeling opportunity lies at the intersection of the two groups: degradation models can be integrated with a chemical kinetic representation of the aggregation process. In this way, they can benefit from the kinetic rate constants derived from fitting to the wide variety of in vitro data available. At the same time, chemical kinetic models can be extended to include UPP and ALP pathways, converging toward an in vivo representation. These integrated models can be employed to determine the relative importance of each microscopic event in relation to degradation mechanisms (e.g., the structural interconversion between type-A and type-B oligomers), thus suggesting specific target species and reactions. Recent studies have moved a few steps toward model integration [73, 96]. Here, clearance mechanisms are approximated by pseudo-first-order reactions to facilitate theoretical analysis. Further investigation should focus on an explicit representation of the degradation processes, especially the ALP pathway, considering the growing interest in the role of autophagy in neurodegeneration [7, 22, 97]. However, this choice translates into an increased complexity of the system itself. Future work should favor an integrative representation of pathological mechanisms that balances high granularity and analytical or computational tractability [98].

4.2. In vitro/in vivo translation

Integrating chemical kinetic and degradation models is part of the modeling efforts currently chasing the mechanistic representation of αsyn aggregation in vivo. To this end, both groups of models would benefit from calibration and validation against measurements of aggregate levels in living systems. However, the availability of this quantitative information is limited by technical difficulties in monitoring the time evolution of aggregates and addressing the heterogeneous and elusive nature of oligomers [58]. As a result, chemical kinetic models have been trained mainly on in vitro data. Among the limitations of using this type of data is the unclear relevance of the results on a disease time scale, moving from the few hours or days of test-tube experiments to the decades of PD progression [21]. In addition, in vitro aggregation assays typically rely on initial protein concentrations that far exceed αsyn levels detected intracellularly and in post-mortem brain extracts from PD patients and controls [99–101]. On the other hand, concentrations used in vitro are comparable to those estimated in vivo at presynaptic terminals [102–105]. Such a potential discrepancy in αsyn concentrations may hamper translatability due to the strong dependence of aggregation rates on initial protein levels [58, 106]. Nevertheless, studying aggregation mechanisms in silico can bridge in vitro and in vivo settings. For example, in vitro studies paired with chemical kinetic models showed that adding lipid vesicles to the solution accelerates in vitro aggregation, which occurs at lower bulk concentrations [47, 48]. This result may suggest a similar mechanism in living systems, given αsyn intracellular localization at presynaptic terminals and its physiological role in vesicular trafficking. Moreover, according to a recent perspective paper [58], the late onset of amyloid formation in sporadic PD may emerge from the intrinsic properties of aggregation kinetics, specifically from a combination of dominant secondary processes and prion-like spreading. Age-related defects in protein degradation may also work in concert with secondary events in the aggregation process to drive amyloid formation. Qualitative similarities between in vitro and in vivo aggregation mechanisms support this scenario. Indeed, biophysical features of amyloid formation detected in test-tube experiments seem to be preserved in living systems [58, 107]. On the other hand, a quantitative analysis of the aggregation kinetics in vivo is lacking. Current studies focus on the applicability of in vitro molecular-level mechanistic studies to living systems and the related modeling limitations and opportunities, as extensively discussed in Meisl et al. [108]. Molecular models combining the kinetic of αsyn aggregation and multiple biological pathways underlying PD neurodegeneration can assist in vitro to in vivo translation [109].

4.3. Modulating the protein aggregation pathway

PD research efforts seek to provide a solid foundation for therapeutic strategy design in PD. Not only can mechanistic models of αsyn homeostasis identify candidate target mechanisms, but they can also support the development of various anti-aggregation treatments currently under study. In this context, chemical kinetic models provide a theoretical framework for analyzing the inhibitory effect of multiple compounds on the aggregation process. Employing a nucleation-polymerization model that accounts for interactions with a new general component, Michaels et al. [110] identified different regimes of inhibition related to specific combinations of kinetic rate constants and binding affinities, thus optimizing the inhibitor efficacy of the compound. Furthermore, a follow-up study simulated the impact of therapeutic interventions on key oligomeric features by tuning specific rate constants [111]. This analysis suggested that inhibiting the primary nucleation event that dominates lipid-induced αsyn aggregation strongly reduces the peak and overall concentration of the oligomeric population. As shown in [112], this mechanism of action characterizes the natural aminosterol compound known as trodusquemine, which interferes with αsyn aggregation by blocking its interaction with lipid membranes and thus inhibits lipid-induced primary nucleation and fibril-mediated secondary nucleation. Also, it can suppress aggregation-related cytotoxicity by hampering oligomer interaction with cell membranes [113]. In addition to small molecules, potential therapeutic strategies targeting different αsyn species and reactions include peptides and peptidomimetics, antibodies, and molecular chaperones [13]. In this respect, chemical kinetic studies have recently investigated the mechanism of action of specific molecular chaperones [114, 115] and antibodies [116] on amyloid beta aggregation in Alzheimer's disease (e.g., Brichos and Aducanumab, respectively); the results were in good agreement with clinical outcomes. These studies should be extended to Parkinson's disease, e.g., supporting clinical trials that involve anti-αsyn antibodies such as prasinezumab and cinpanemab.

4.4. Targeting protein degradation pathways

Inhibiting specific molecular steps of the aggregation pathway is not the only possible intervention targeting impaired αsyn homeostasis. Genetic-based strategies reducing protein synthesis have been proposed to avoid aggregation in the first place, such as employing antisense oligonucleotides, RNA interference, and β2-adrenoreceptor agonists [13]. However, side effects may arise from protein downregulation due to αsyn physiological role in presynaptic vesicle trafficking. Mechanistic models should be employed to determine a critical range for αsyn expression levels, facilitating the transition from preclinical to clinical trials. However, currently available models do not include a genetic level representation; thus, transcription and translation inhibitions can be simulated only by varying the rate constant associated with αsyn production. Another strategy consists in modulating UPP and ALP degradation mechanisms to hamper αsyn aggregation and PD neurodegeneration. Autophagy can be enhanced selectively or as a whole. An overall ALP enhancement driven by chemical compounds such as rapamycin and trehalose has limited applicability in PD since it may interfere with essential physiological processes and the homeostasis of proteins other than αsyn [62]. The same limitations hold for therapeutic strategies focusing on overall UPP degradation [117]. In contrast, selective targeting of specific molecular steps in proteasomal degradation, CMA, and macroautophagy may be a promising strategy [62]. UPP-related candidate treatments include stimulating ubiquitination and regulating the activity of chaperone systems. In addition, selective ALP activators may act by modulating LAMP2A assembly and translocation, enhancing autolysosome fusion, and targeting lysosomal function. Ideally, mechanistic models can support the development of current and future interventions of this type. However, further research is needed to reach this goal. Most degradation models analyzed here lack a granular representation of autophagic pathways that would allow the simulation of these selective strategies. Still, they can provide insightful information on the role of UPP and ALP degradation pathways in neurodegeneration as the first step toward in silico-aided therapeutic strategies against PD.

4.5. Toward multiscale QSP models for PD

The points discussed above highlight the need for multiscale QSP models that can link molecular mechanisms to cellular and tissue-level events and, finally, clinical manifestations. These models would leverage the large amount of available data associated with neuroimaging techniques and disease-rating scales for motor symptoms (MDS-UPDRS) and cognitive impairment (MoCa and MMSE) [1, 118]. This type of data has already been employed by quantitative models that address the heterogeneity of the PD clinical picture by accounting only for physical symptoms [119]. In addition to quantitative tools supporting PD diagnosis, a mechanistic understanding of the molecular processes underlying PD pathogenesis and how they connect to the clinical endpoints is essential to inform drug discovery and development. The models in this review represent a starting point in this sense. Being confined at the molecular level, they have to settle for surrogate endpoints such as ROS, αsyn aggregate, and LB levels to identify disease phenotypes. Moving to the cellular, tissue, and organ levels through network-based, epidemiological, and compartmental PK/PD approaches [120–123], these models can account for αsyn prion-like propagation between neurons and across multiple brain areas, using potential PD biomarkers such as αsyn concentrations in plasma and CSF [124] and magnetic resonance imagining measurements.

Overall, this review makes a case for using mechanistic models of αsyn homeostasis as valuable tools to uncover PD pathogenesis, suggesting their integration into a QSP framework to provide blueprints for therapeutic interventions.

Author contributions

ER and FR designed, wrote, and edited the manuscript. AA provided a preliminary version of some sections. LM and ED discussed the results and edited the manuscript. All the authors reviewed the manuscript.

Acknowledgments

The authors would like to thank Dr. Gianluca Selvaggio for his helpful suggestions, comments, and support in preparing this manuscript.

Conflict of interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Publisher's note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

References