Dynein Dysfunction Prevents Maintenance of High Concentrations of Slow Axonal Transport Cargos at the Axon Terminal: A Computational Study

Neurons utilize a complicated “railway” system, composed of microtubule (MT) tracks and molecular motors that move on these tracks, to transport large proteins and vesicles across large distances. The motors themselves belong to the kinesin family, which perform anterograde transport, or the dynein family, which perform retrograde transport [1].

Motor protein dysfunction plays a critical role in many neurodegenerative diseases [2–4]. Dynein mutations have been implicated in many cases of malformations of cortical development, spinal muscular atrophy with lower extremity dominance, and congenital muscular dystrophy. Another example is Charcot-Marie-Tooth (CMT) disease, a disorder of peripheral motor and sensory neurons [5]. CMT type 2O disease is caused by a mutation in gene DYNC1H1 that encodes the dynein heavy chain [5–8]. Although the involvement of molecular motors in CMT pathogenesis is not surprising (long sensory and motor neurons require a large anterograde flux of various cargos from the soma to support their remote terminals [9]), an intriguing question is why degeneration of these neurons is associated with dynein dysfunction. Previous hypotheses as to why neurons can be negatively affected by dynein dysfunction include an inability to deliver dysfunctional organelles/proteins to the soma for degradation or a failure of retrograde signaling [10]. Bidirectional transport may be also needed to correct errors if cargo is delivered to a wrong location [11].

Here, we consider another possibility that is associated with our hypothesis that dynein-driven transport is necessary for cargo transport against its concentration gradient. Consider a mathematical argument: when cargo diffusivity is small, a continuum model of anterograde transport is described by a first-order differential equation, which allows for the imposition of a boundary condition only at the axon hillock. The equation does not allow for the imposition of the second boundary condition (prescribing a higher cargo concentration) at the axon tip. To simulate an increased cargo concentration at the axon tip, both anterograde and retrograde components are needed [12]. Kuznetsov and Kuznetsov [12] used this argument to explain the presence of a retrograde component in slow axonal transport. Note that according to our hypothesis dynein-driven transport does not only serve for signaling to the soma to correctly modulate protein synthesis and to regulate forward kinesin transport. Rather, anterograde and retrograde transport are coupled at any location along the axon to form a physical system that allows specifying cargo concentrations at both ends of the axon, at the hillock and at the terminal.

In this paper, we attempt to answer the question of how dynein-driven transport can affect the transport of cargos toward the axon terminal. We will focus on slow axonal transport-b (SCb), which is used to transport ∼200 proteins from the soma to the presynaptic terminal [13]. SCb transports cargo at an average velocity of 2–8 mm/day (0.023–0.093 μm/s), much less than the 1 μm/s observed for kinesin and dynein motors. This slower SCb transport velocity is explained by significant time cargos spend in the pausing state. Another interesting feature of SCb is that during the rapid phase of motion the cargos move bidirectionally with a bias toward anterograde motion [14,15]. An intriguing question raised in Ref. [10] and by many other researchers is then: why is slow axonal transport bidirectional if cargoes just need to be transported to the axonal tip? Since molecular motors use ATP, bidirectional transport is much more energetically expensive than unidirectional transport. One possible explanation is that bidirectional cargo movements are caused by opposite polarity motors that are attached to the same cargo and that are in tug-of-war [16,17]. Here, we suggest that bidirectional transport is required to maintain an increasing cargo concentration along the length of the axon, in the situation when the cargo concentration is small at the soma and high at the presynaptic terminal. This explanation first reported in Ref. [12] holds for the case when cargos are transported as polymers or large multiprotein complexes that have small diffusivity, which is the case for many SCa and SCb proteins.

We focus on one particular slow axonal transport cargo, α-synuclein (α-syn), mostly known for its involvement in Parkinson's disease [18–23]. Since in a healthy neuron α-syn predominantly exists in the monomeric form [24], hereafter we denote α-syn monomer as α-syn.

Our model is a cargo-level model that simulates the behavior of cargos rather than the behavior of motors, as a specific cargo can be driven by several motors. We simulate dynein dysfunction as a decrease of dynein velocity. This is equivalent to specifying the velocity of retrograde cargos during the fast phase of their movement on MTs. Note that cargos can also pause when the motors that drive them temporarily disengage from MTs.

Our goal is to simulate how dynein dysfunction affects slow axonal transport of α-syn. Although we perform computations using α-syn parameters, we expect our results (the need for the retrograde component) to be generalizable to other slow axonal transport cargos.

A schematic representation of an axon is shown in Fig. 1. We start with the simplified model of bidirectional transport suggested in Ref. [12]. Since long-distance transport of α-syn is propelled by kinesin and dynein motors [25], the simplified model includes only two kinetic states, which simulate anterograde (driven by kinesin) and retrograde (driven by dynein) cargos (Fig. 2(a)).

where asterisks denote dimensional quantities, $x∗$ is the Cartesian coordinate that protrudes from the axon hillock (⁠$x∗=0$⁠) to the axon tip (⁠$x∗=L∗$⁠). Sticky speaking, the positions of axon tips are shown in Fig. 1; however, we do not resolve complicated processes occurring in the terminal, our model is only valid between $0≤x∗≤L∗$⁠. $va∗$ is the velocity of rapid motions of α-syn on MTs propelled by kinesin motors in slow axonal transport, and $γar∗$ and $γra∗$ are kinetic constants defined in Fig. 2(a). The first term on the left-hand side of Eq. (1) characterizes a change of the number of α-syn molecules in the control volume (CV) due to anterograde transport of α-syn, the second term characterizes a decrease of the number of α-syn molecules due to their transition to the retrograde state, and the third term characterizes an increase of the number of α-syn molecules due to their transition to the anterograde state (Fig. 2(a)).

where $vr∗$ is the velocity of rapid motions of α-syn on MTs propelled by dynein motors in slow axonal transport. The first term on the left-hand side of Eq. (2) characterizes the effect of retrograde transport of α-syn due to the action of dynein motors while the second and third terms on the left-hand side are the kinetic terms characterizing the effects of α-syn transitions to the anterograde and retrograde kinetic states, respectively.

where $jtot,x=0∗$ is the flux of α-syn at the axon hillock.

where C is the integration constant.

Equation (22) shows that dynein dysfunction prevents the maintenance of a higher concentration of α-syn at the axon tip than at the soma.

We now move on to the full slow axonal transport model (kinetic diagram is shown in Fig. 2(b)). The continuum slow axonal transport model for simulating neurofilament transport was previously developed in Ref. [38]. The model was extended in Refs. [39] and [40] to simulate slow axonal transport of cytosolic proteins by adding a kinetic state for proteins that freely diffuse in the cytosol as well as a degradation term for the destruction of freely moving cytosolic proteins in proteasomes.

Equation (23) is similar to Eq. (1), but kinetic constants are different and are now defined in Fig. 2(b). The first term on the left-hand side of Eq. (23) describes the effect of kinesin-driven transport of α-syn while the second and third terms describe transitions between various kinetic states.

where the first term describes the effect of dynein-driven transport and the second and third terms are kinetic terms describing α-syn transitions between kinetic states. This equation is identical to Eq. (2), but kinetic constants are different.

Because α-syn does not move in the pausing states, Eqs. (25) and (26) include only kinetic terms, which describe transitions between different kinetic states (Fig. 2(b)).

Note that in Eqs. (23)–(27) the degradation term is present only in the free state. Indeed, to enter a proteasome, α-syn must be detached from MTs. Thus, proteins moving on MTs are not subject to degradation in proteasomes, unlike those that are freely diffusing in the cytosol. This must be the case because in long neurons the time that it takes for a protein to travel from the soma to the axon tip exceeds the lifetime of free proteins. Therefore, it is likely that during their transit in axons proteins are protected from degradation [43,44]. If proteins can only be destroyed in the free state, this would explain how this protection is accomplished. Another possibility is that since cytosolic proteins are transported in SCb in the form of large cargo structures (multicargo complexes), they may be too large to enter proteasomes and cannot be degraded at all during transport. In this case, the half-life of α-syn complexes in the free state should be close to infinity.

If α-syn diffusivity $Dnfree∗$ in Eq. (27) is not zero, Eqs. (23)–(27) must be solved subject to four boundary conditions. At the axon hillock, we imposed the total concentration of α-syn at the hillock, $ntot,x=0∗$⁠, and the flux of α-syn entering the axon, which equals to $ntot,x=0∗$ times the average velocity of α-syn (which includes pauses), $vav,estimate∗.$

Unless proteins are MT-bound, such as tau or MAP2, slow axonal transport cargos are usually not transported along the axon as single particles but in protein complexes [10,13,41,46–50]. Because diffusivity of such multiprotein complexes is expected to be very small, we investigated how Eqs. (23)–(27) can be simplified for the limiting case of vanishing diffusivity of free proteins. Our goal is to understand whether after such simplification the model is capable of simulating cargo transport against a cargo concentration gradient.

it can be viewed as the ratio of the convection time to diffusion time.

Note that the solution given by Eq. (66) has the form $jtot,x=0∗va∗exp(−ax∗)$⁠, where the concentration decay constant $a>0$⁠, which means that $na,0∗(0)$ cannot describe an increasing concentration toward the terminal. If a is small, the solution is $na,0∗(0)≈jtot,x=0∗va∗$⁠, which is a constant value.

The total concentration, $ntot,0∗(0)$⁠, can be obtained by finding the sum of concentration components given by Eqs. (59) and (66)–(69).

where $Δvr∗=10−3vr∗$ is the step size. The independence of the sensitivity to the step size was tested by varying the step sizes.

where $k∗$ is the logistic growth rate, a parameter that characterizes the steepness of the S-shaped curve used in approximating the distribution of α-syn along the axon length. We used $k∗=10−3$ 1/μm. The distribution given by Eq. (74) is shown by hollow circles in Fig. 3.

The purpose of the second term on the right-hand side of Eq. (73) is to simulate the difference between the numerically predicted α-syn velocity, $vav,j∗$⁠, and experimentally reported α-syn velocity, $vav,estimate∗$⁠. We used 0.05 μm/s for $vav,estimate∗$⁠, a value which is in the middle of the experimentally reported range of SCb velocity, 2–8 mm/day [13]. In the second term, we also used $ω1$=1 s²/μm². This value was selected based on numerical experimentation, to avoid overfitting either α-syn concentration or its average velocity, see Ref. [65].

The third term on the right-hand side of Eq. (73) is set to a large value, $ω2=108$⁠, if any of α-syn concentration components, $na,j$⁠, $nr,j$⁠, $na0,j$⁠, $nr0,j$⁠, $nfree,j$⁠, or α-syn flux, $jtot$ (j = 1,…,$Nfit$⁠) becomes negative. The latter condition comes from the assumption that α-syn is not synthesized at the terminal, and that the half-life of α-syn is finite; therefore, the terminal must be supplied by α-syn from the soma.

Best-fit values of kinetic constants for the anterograde-retrograde transport model are given in Table S4 available in the Supplemental Materials, and the best-fit values of kinetic constants for the full slow axonal transport model are given in Table S5 available in the Supplemental Materials.

Equations (1) and (2) were solved using MATLAB's BVP5C solver (MATLAB R2020b, MathWorks, Natick, MA). For the full slow axonal transport model, we first eliminated $na0∗(x∗)$ and $nr0∗(x∗)$ from Eqs. (23)–(27) using Eqs. (25) and (26). The resulting system of ordinary differential equations of the fourth order was again solved using MATLAB's BVP5C solver. To determine values of kinetic constants (Tables S4 and S5 available in the Supplemental Materials), we used Multi-Start with the local solver fmincon. These routines are part of MATLAB's Optimization Toolbox. We used 1000 random points plus the initial point given in Table S3 available in the Supplemental Materials as starting points in searching for the global minimum. We followed the numerical procedure described in Ref. [66].

We used the anterograde-retrograde model (Fig. 2(a)) to simulate transport of an SCb protein that is predominantly localized in the presynaptic terminal. α-syn data were used for our simulation. We found that the concentrations of proteins pulled by anterograde and retrograde motors remain constant along most of the axon length and then sharply increase near the axon tip (Figs. 4(a) and 4(b)). A decrease of the dynein velocity results in a decrease of the boundary layer thickness (i.e., the region where cargo concentrations change from the constant low concentration which is maintained along most of the axon length to the high concentration seen at the axon tip; see Figs. 4(a) and 4(b)). If the dynein velocity approaches zero, the anterograde and retrograde concentrations are represented by horizontal lines (no concentration increase toward the axonal tip, see the curves corresponding to the analytical solution in Figs. 4(a) and 4(b)).

A similar decrease of the boundary layer thickness with a decrease of the dynein velocity is observed in Fig. 3 which displays the total cargo concentration. In the model given by Eqs. (1) and (2) the boundary layer thickness is very sensitive to the dynein velocity. The dimensionless sensitivity, $Svr∗δ∗$⁠, calculated using Eq. (72) for this model is 9.21. Therefore, in addition to the base case, corresponding to $vr∗=0.7$μm/s, we show two more concentration distributions where the dynein velocity is decreased by factors of 0.8 and 0.5 (see Figs. 3 and 4). Notably, further decreasing dynein velocity causes the boundary layer thickness to become small, leading to an unphysiological near-vertical transition in concentration distribution.

The values of kinetic constants used for computing all curves in Figs. 3, 4, and S1 (the latter is available in the Supplemental Materials) were determined using $vr∗=0.7$μm/s. Otherwise, if kinetic constants were optimized for each value of $vr∗$⁠, all three numerical curves would be the same and would give a perfect fit with the synthetic data. Numerically computed curves in Figs. 3 and 4 exhibit a constant concentration along most of the axonal length which is followed by a sharp concentration increase near the axon tip. This feature is due to the concentration distribution of synthetic data (see Eq. (74)). The sharp increase of the total concentration near the terminal exhibited by synthetic data (see hollow circles in Fig. 3) simulates the presynaptic localization of α-syn. Since values of kinetic constants in the model are determined such that numerically predicted total cargo concentration would fit the synthetic data distribution, the shapes of the numerical curves in Fig. 3 correspond to the synthetic data distribution. The important feature in Fig. 3 is the inability of the model to simulate the increase of cargo concentration toward the axonal tip if dynein velocity tends to zero (see the curve that represents the analytical solution).

The total flux of cargos (defined in Eq. (5)) is constant independent of the position in the axon (Fig. S1 available in the Supplemental Materials on the ASME Digital Collection) because the cargo flux at the hillock is imposed by Eq. (3) and there is no cargo degradation in the model given by Eqs. (1) and (2).

For the full slow axonal transport model (Eqs. (23)–(27), Fig. 2(b)) we observed that predicted cargo concentrations are less sensitive to the dynein velocity than in the model given by Eqs. (1) and (2). The dimensionless sensitivity, $Svr∗δ∗$⁠, calculated using Eq. (72) for the full slow axonal transport model is 1.28. Therefore, in Figs. 5, 6, and S2–S4 available in the Supplemental Materials, we have shown results for the base dynein velocity (⁠$vr∗=0.7$μm/s) and two much smaller values of dynein velocity, $vr∗=0.7×10−1$ and $0.7×10−2$μm/s. Again, the values of kinetic constants used for computing all curves in Figs. 5 and 6, and S2–S4 (the latter are available in the Supplemental Materials) were determined using $vr∗=0.7$μm/s. The analytical solution obtained for the case when the cargo diffusivity and dynein velocity both approach zero predicts almost uniform cargo concentrations (Figs. 5, 6, S2, and S3a, the latter two are available in the Supplemental Materials). Indeed, for the results displayed in Figs. 5 and 6 the concentration decay constant a defined in the paragraph after Eq. (66) equals to $3.59×10−13$ 1/μm. This means that the exponential function in Eq. (66) is almost identical to unity and $na,0∗(0)$ is constant. This explains why the model predicts that without a retrograde motor neurons are incapable of supporting high cargo concentration at the axon terminal. It should be noted that even if the soma can sense the low cargo concentration at the axon terminal (if some retrograde signaling is present) and starts producing more α-syn, it will only shift the horizontal “analytical” α-syn distribution displayed in Fig. 6 upward, but the distribution will remain horizontal.

For the full slow axonal transport model, the cargo flux displayed in Fig. S3b available in the Supplemental Materials on the ASME Digital Collection remains constant because of the large value of cargo half-life used in computations (Fig. S3b available in the Supplemental Materials is computed for $T1/2,free∗=5.76×1010$ s, assuming that α-syn is protected from degradation during its transport in the axon). Computations using a smaller cargo half-life (⁠$T1/2,free∗=5.76×104$ s) predict a decay of the total flux of cargos near the axon tip (Fig. S4b available in the Supplemental Materials). This is because the model assumes that α-syn can be destroyed in proteasomes only in the free (cytosolic) state. Since the concentration of free α-syn increases toward the tip (Fig. S4a is available in the Supplemental Materials), there is more α-syn degradation near the tip. Hence, the total flux of α-syn decreases near the tip (Fig. S4b is available in the Supplemental Materials).

The results obtained in this paper suggest that dynein dysfunction results in the inability of neurons to transport cargo against a concentration gradient. This result holds when cargo diffusivity is small, which applies to cargos transported as polymers or as a part of multiprotein complexes. A mathematical argument to explain why dynein motors are needed to support a cargo concentration distribution that increases toward the axon tip is as follows. With the anterograde-only motion of cargo there is no way to impose a boundary condition at the axon tip that is to require an elevated cargo concentration at the tip. Our results, counterintuitively, suggest that the progressive distal to proximal (dying-back) axonal degeneration may result from dynein dysfunction because of the inability of the neuron to support a high cargo concentration at the terminal.

Limitations of our model are related to neglecting local protein synthesis and local controlled protein degradation. These effects should be included in future models and considered in future research.

Our results are testable and falsifiable. Future work could test whether neurons with dysfunctional dynein are indeed unable to support a cargo distribution in an axon that increases with distance from the soma. Ref. [11] reported that dynein knockdown causes arrest of bidirectional peroxisome transport in Drosophila. However, cargo transport can be rescued by artificially attaching a construct containing another retrograde motor (kinesin-14, Ncd) to peroxisomes. It would be interesting to investigate whether dysfunctional dynein can be replaced by another retrograde motor to support slow axonal transport.

IAK acknowledges the fellowship support of the Paul and Daisy Soros Fellowship for New Americans and the NIH/National Institute of Mental Health (NIMH) Ruth L. Kirchstein NRSA (F30 MH122076-01). AVK acknowledges the support of the National Science Foundation (award CBET-2042834) and the Alexander von Humboldt Foundation through the Humboldt Research Award.

• Division of Chemical, Bioengineering, Environmental, and Transport Systems (Award No: CBET-2042834; Funder ID: 10.13039/100000146).

• National Institute of Mental Health (Award No. F30 MH122076-01; Funder ID: 10.13039/100000025).

The datasets generated and supporting the findings of this article are obtainable from the corresponding author upon reasonable request.