## Introduction

Dense core vesicles (DCVs) transport neuropeptides, neurotrophins, and enzymes toward the axon terminal (Wong et al. [1]). DCVs are synthesized in the soma (Li et al. [2]). Malfunctions in DCV transport can lead to devastating neurodegenerative diseases. Bulgari et al. [3] reported that the loss of Huntingtin protein results in more capture of retrogradely transported DCVs and thus more DCV accumulation in the nerve terminal. They also suggested a possible link between DCV transport and Huntington's disease, caused by accumulation of mutant Huntingtin protein and loss of wild type Huntingtin, through modification of synaptic capture of DCVs. Kwinter et al. [4] pointed out that understanding DCV transport may be also important for understanding the molecular basis of Alzheimer's disease and amyotrophic lateral sclerosis.

By performing a careful experimental investigation of Drosophila nerve terminals, Wong et al. [5] found that in type Ib terminals, DCVs circulate through the terminal, moving first anterogradely and then retrogradely, and that only a portion of DCVs are captured when passing each of the en passant boutons (hereafter boutons). This discovery was important to understand how DCVs transported to the nerve terminal end up being equally distributed across the boutons, without being preferentially supplied to the most proximal bouton (Wong et al. [5]). The principle of equal distribution of resources is referred to as synaptic democracy (Bressloff and Levien [6], Bressloff and Karamched [7]).

It is well known that, by stating the conservations of mass, energy, and momentum, useful information can be obtained about a thermodynamic system. Here, we treat the number of DCVs as a conserved property and state the conservation of this property in five separate control volumes (compartments). These compartments include the axon and four boutons composing the nerve terminal (Fig. 1). We previously applied this approach to investigating DCV transport in type Ib nerve terminals in Kuznetsov and Kuznetsov [8]. We extended the model to simulate type III nerve terminals in Ref. [9], and further extended the model to simulate various rates of DCV synthesis in the soma in Ref. [10]. In the present paper, we make two further extensions to our model. The first extension deals with simulating DCV capture in type III boutons. Kuznetsov and Kuznetsov [10] assumed that type III boutons exhibit an increased capture efficiency (compared to type Ib boutons) when DCVs move through the boutons in both anterograde and retrograde directions. However, recent experimental results obtained by Levitan's lab [11] suggest that type III boutons demonstrate the increased capture efficiency only for DCVs that move anterogradely. The capture efficiency for retrogradely moving DCVs is comparable to type Ib boutons. In this paper, we modified the model to account for this recent experimental finding. This required the introduction of three new mass transfer coefficients in our model, which characterize the rates of capture of retrogradely moving DCVs in boutons, and the development of a method for estimating their values.

The second extension deals with considering different possible fates for captured DCVs. In our previous papers [810], we assumed that all captured DCVs are eventually destroyed in boutons. However, this assumption may be difficult to reconcile with the fact that organelles like DCVs are usually destroyed in lysosomes, which are abundant in the soma but not in the terminals. In this paper, we re-formulated our equations and extended the model to allow for the possibility that DCV capture in boutons may be reversible and that at least some DCVs may be re-released from boutons to reenter the circulation. Our new modeling results show that the assumption about the fate of DCVs in boutons has profound effect on the presence or absence of DCV circulation in the nerve terminals, as well as on its strength. This extension allows us to pose questions such as how the fate of DCVs captured in boutons affects the presence or absence of DCV circulation in nerve terminals. With currently employed experimental techniques, measurements alone make it difficult to distinguish between possible fates of DCVs in boutons (destruction in the boutons or leaving the boutons, e.g., to be destroyed in somatic lysosomes, Levitan [11]). This means that the development of our model so that it can make predictions about different consequences given the fate of DCVs would be very important for future research.

The model has a degree of universality because it can simulate DCV transport in type Ib and type III terminals, which have different morphology (Atwood et al. [12], Jia et al. [13], Menon et al. [14], Bulgari et al. [15]). The governing equations are the same; the only difference is how the model parameters are estimated.

It should be noted that it is not our goal to simulate the specific experiments reported in Wong et al. [1,5], Bulgari et al. [15], and Cavolo et al. [16]. Our goal is to develop a model capable of simulating the DCV circulation discovered by Levitan and colleagues and reported in the aforementioned references.

Predictions of our model can be tested and either proven or disproven. These predictions are important because it is difficult to experimentally investigate the fate of DCVs in boutons. In this paper, we compare the predicted DCV circulation under three different scenarios: (i) DCV capture in boutons is completely reversible and all captured DCVs eventually reenter the transiting pool (thus reentering the circulation); (ii) DCV capture in boutons is partly reversible and half of captured DCVs reenter the transiting pool while half are destroyed; and (iii) capture of DCVs in boutons is irreversible and all captured DCVs are eventually destroyed in boutons. The model published in Kuznetsov and Kuznetsov [10] considered only scenario (iii).

## Materials and Models

### Governing Equations.

We assumed a healthy axon (no disruptions of DCV transport) and we neglected any possible feedback effects on DCV transport. The latter is consistent with Wong et al. [5], who assumed that DCV delivery to boutons is controlled by DCV circulation and sporadic capture from this circulation rather than by any active address system for directing DCVs to a particular bouton. We simulate DCV transport in an axon terminal with two identical branches where each branch contains four boutons (Fig. 1). We used a multicompartment model (Jacquez [17]) to develop the governing equations for this problem. Our model includes five compartments: the axon and four boutons (Fig. 2). As in Wong et al. [5], we numbered the boutons from 4 (the most proximal bouton) to 1 (the most distal bouton). The model can be easily extended to terminals with different numbers of boutons and branches.

Our model accounts only for concentrations of those DCVs that are captured in boutons. DCVs transiting through the boutons and between the boutons are accounted for by simulating DCV fluxes. The variables utilized in our model are summarized in Table 1. We split model parameters in two groups, those whose values can be estimated from the literature (Table 2) and those that need to be estimated by considering balances for DCV concentrations and fluxes (Table 3). We estimated the latter parameters by considering DCV concentrations and fluxes in the beginning of the process (t = 0) and at steady-state.

Since the length of an axon is much larger than its diameter, we characterized the DCV concentration by its linear number density, expressed in vesicles per unit length of the axon. The concentration of DCVs in the resident state in a bouton can change due to DCV capture into the resident state or due to DCV destruction or reentering the transiting pool. Hereafter, we use “captured state” and “resident state” as synonyms. We also use “saturated concentration” and “steady-state concentration” as synonyms, and we use “DCV destruction” and “DCV degradation” as synonyms. In Fig. 2 DCV capture in boutons is shown by block arrows. Note that in boutons 4, 3, and 2 the capture rate of anterogradely moving DCVs may be different from the capture rate of retrogradely moving DCVs. We start by stating the conservation of captured DCVs in the most proximal bouton:
$L4dn4dt=min[h4a(nsat0,4−n4),jax→4]+H[t−t1]min[h4r(nsat0,4−n4),j3→4]−L4n4ln(2)T1/2$
(1)

The term on the left-hand side of Eq. (1) simulates DCV accumulation in the resident state.

The first two terms on the right-hand side of Eq. (1) simulate the rates of DCV capture into the resident state as DCVs pass bouton 4 in anterograde and retrograde directions, respectively.

We assumed that the rate of DCV capture is proportional to the difference between the DCV concentration in a bouton assuming infinite DCV half-life or infinite half-residence time in the bouton and the current DCV concentration in the bouton, $nsat0,i−ni$. The use of $nsat0,i$ is necessary so that DCV capture would occur even at steady-state to replace the loss of captured DCVs due to their destruction or reentering the circulation.

The number of anterogradely moving DCVs captured by bouton 4 cannot exceed the number of anterogradely transported DCVs entering bouton 4. This is simulated by the first term on the right-hand side of Eq. (1), which takes the lesser of these two values: (i) the capture rate that is determined by conditions in the bouton, provided there are enough DCVs transiting anterogradely through bouton 4 and (ii) the flux of DCVs entering bouton 4 from the axon.

The same is true for retrogradely transported DCVs; their capture rate is the lesser of (i) the DCV flux entering bouton 4 from bouton 3 and (ii) the rate of capture in bouton 4, provided unlimited supply of retrogradely moving DCVs. This is simulated by the second term on the right-hand side of Eq. (1). The multiplier $H[t−t1]$ in this term accounts for a 300 s delay that DCVs need to change anterograde to retrograde motors in bouton 1 (Levitan [11]).

The last term on the right-hand side of Eq. (1) simulates DCV removal from the resident state by DCV destruction or reentering the circulation. Note that for the purposes of stating the conservation of the number of DCVs in the resident state, it does not matter whether DCVs are destroyed or reenter the circulation (the particular mechanism by which DCVs leave the resident state is not important).

The equations representing the conservation of DCVs in boutons 3 and 2 take the form
$L3dn3dt=min[h3a(nsat0,3−n3),j4→3]+H[t−t1]min[h3r(nsat0,3−n3),j2→3]−L3n3ln(2)T1/2$
(2)

$L2dn2dt=min[h2a(nsat0,2−n2),j3→2]+H[t−t1]min[h2r(nsat0,2−n2),j1→2]−L2n2ln(2)T1/2$
(3)
The equation representing the conservation of DCVs in the most distal bouton (bouton 1) takes the form
$L1dn1dt=min[h1(nsat0,1−n1),j2→1]−L1n1ln(2)T1/2$
(4)

Unlike boutons 4, 3, and 2, DCVs pass bouton 1 only once. Therefore, there is no difference between anterograde and retrograde capture in bouton 1.

The equation stating the conservation of captured DCVs in the axon is
$Laxdnaxdt=jsoma+2j4→ax−2jax→4−Laxnaxln(2)T1/2,ax$
(5)

The parameter $jsoma$ in Eq. (5) should be viewed as the net rate of DCV production in the soma (the rate of DCV synthesis less the rate of DCV destruction in somatic lysosomes).

To achieve closure, Eqs. (1)(5) must be supplemented with equations simulating DCV fluxes between the axon and bouton 4 as well as DCV fluxes between the boutons (Fig. 2). The DCV flux from the axon to the most distal bouton is assumed to be proportional to the average DCV concentration in the axon
$jax→4=hinnax$
(6)

Here, we extended our previous model (Kuznetsov and Kuznetsov [10]) by assuming that captured DCVs can, after spending some time in a bouton, escape and re-enter the circulation. We used parameter $δ$ to characterize the portion of DCVs that escape from the captured state in boutons back into the circulation. $δ=0$ corresponds to the case when all captured DCVs are eventually destroyed in boutons. This case was investigated in Kuznetsov and Kuznetsov [10]. $δ=1$ corresponds to the case when all captured DCVs eventually reenter the circulation. It is also possible that those DCVs that escaped from the captured state in boutons are tagged for degradation in somatic lysosomes and return to the soma by retrograde transport (Rizzoli [18]) after they escape from boutons. For modeling purposes, it does not matter whether DCVs are destroyed in boutons or in the soma, as long as they are destroyed and do not re-enter the circulation. Therefore, this case can be simulated by $δ=0$. It should be noted that the assumption that the released DCVs return to the soma for degradation requires a nonzero retrograde flux from the terminal back to the axon, even in neurons with type III terminals.

In order to write the equations for the fluxes, it is necessary to postulate how the DCVs that escape from the captured state are split between the anterogradely and retrogradely moving populations. This should be done only for boutons 4, 3, and 2 because DCVs can only leave bouton 1 by retrograde flux. We assumed that in boutons 4, 3, and 2 DCVs are split in the following way: $ε$ of DCVs that escaped from the captured state join the anterograde flux and $(1−ε)$ join the retrograde flux. Equations for $j4→3$, $j3→2$, and $j2→1$ (Fig. 2) are then written utilizing the following DCV conservation requirement. The flux of anterogradely moving DCVs leaving a bouton is equal to the flux of anterogradely moving DCVs entering the bouton (the first term on the right-hand side) minus the rate of capture of anterogradely moving DCVs (the second term) plus the rate of DCV escape from the captured state into the anterogradely moving component of transiting DCVs (the third term)
$j4→3=jax→4−min[h4(nsat0,4−n4),jax→4]+δεL4n4ln(2)T1/2$
(7)

$j3→2=j4→3−min[h3(nsat0,3−n3),j4→3]+δεL3n3ln(2)T1/2$
(8)

$j2→1=j3→2−min[h2(nsat0,2−n2),j3→2]+δεL2n2ln(2)T1/2$
(9)
In writing equations for $j1→2$, $j2→3$, $j3→4$, and $j4→ax$ we assumed that the flux of retrogradely moving DCVs leaving a bouton is equal to the flux of retrogradely moving DCVs entering the bouton (the first term in the curly brackets) minus the rate of capture of retrogradely moving DCVs (the second term) plus the rate of DCV escape from the captured state into the retrogradely moving component of transiting DCVs (the last term on the right-hand side)
$j1→2=H[t−t1]{j2→1−min[h1(nsat0,1−n1),j2→1]}+δL1n1ln(2)T1/2$
(10)

$j2→3=H[t−t1]{j1→2−min[h2(nsat0,2−n2),j1→2]}+δ(1−ε)L2n2ln(2)T1/2$
(11)

$j3→4=H[t−t1]{j2→3−min[h3(nsat0,3−n3),j2→3]}+δ(1−ε)L3n3ln(2)T1/2$
(12)

$j4→ax=H[t−t1]{j3→4−min[h4(nsat0,4−n4),j3→4]}+δ(1−ε)L4n4ln(2)T1/2$
(13)
In Eqs. (6)(13) fluxes have units of vesicles/s. Equations (1)(5) are first-order ordinary differential equations and they require five initial conditions. It should be noted that in Drosophila, release by exocytosis is always incomplete. However, release by exocytosis, which is associated with molting behavior, can be massive (∼90% of content) (Levitan [11]). Therefore, for simplicity, we used the following initial conditions:
$n1(0)=0, n2(0)=0, n3(0)=0, n4(0)=0, nax(0)=nsat,ax$
(14)

Incomplete release by exocytosis can be easily incorporated into the model by modifying the initial conditions.

### Determining Parameter Values Using Published Experimental Data and Physical Estimates.

We estimated values of $jax→4$, $L1$, $L2$, $L3$, $L4$, $Lax$, $t1$, $T1/2$, $T1/2,ax$, $ϖ1$, $ϖ2$, $ϖ3$, and $ϖ4$ based on experimentally reported data and summarized the estimates in Table 2. Next, using DCV conservation and values given in Table 2, we estimated values of $nsat,ax$, $nsat,1$, $nsat,2$, $nsat,3$, $nsat,4$, $h1$, $h2a$, $h2r$, $h3a$, $h3r$, $h4a$, $h4r$, $hin$, $nsat0,1$, $nsat0,2$, $nsat0,3$, $nsat0,4$, $nsat,ax$, and $jsoma$ and summarized them in Table 3.

#### Saturated DCV Concentrations in Boutons, $nsat,i$ (i = 1,…,4).

In type Ib boutons at steady-state the number of DCVs captured in boutons is ∼200 per bouton, and this is the same in all of the boutons (Shakiryanova et al. [19]). Therefore, the saturated concentration of DCVs in type Ib boutons is estimated as $nsat,i=200/Li=40$ vesicles/μm (i = 1,…,4).

As noted in Bulgari et al. [15], the most distal bouton in type III terminals is always smaller than other boutons, while the other boutons are often comparable in size or one randomly located bouton is much larger. Here, as in Kuznetsov and Kuznetsov [10], we simulated a representative situation when bouton 3 is the largest, boutons 4 and 2 are of the same size, and bouton 1 is the smallest (the boutons are numbered as in Fig. 1). The relative size of boutons was quantified by the following equations:
$nsat,2=a2 nsat,1$
(15a)

$nsat,3=a3 nsat,1$
(15b)

$nsat,4=a4 nsat,1$
(15c)

where we used $a2=3$, $a3=9$, and $a4=3$, as in Kuznetsov and Kuznetsov [10].

According to Bulgari et al. [15], there is essentially no DCV circulation in type III terminals. If $δ=0$ (the situation studied in Kuznetsov and Kuznetsov [10]), all captured DCVs are eventually destroyed in boutons. This means that all DCVs that enter the terminal are eventually destroyed in boutons. The conservation of DCVs at steady-state then gives
$0=jax→4−(L1nsat,1+L2nsat,2+L3nsat,3+L4nsat,4)ln(2)T1/2$
(16)

We solved Eqs. (15) and (16), and obtained the following: $nsat,1=26.0$, $nsat,2=77.9$, $nsat,3=233.8$, and $nsat,4=77.9$ vesicles/μm. For consistency, we used the same saturated concentrations in type III boutons for $δ=1$ (all captured DCVs reenter the circulation).

#### Estimating the Saturated DCV Concentration in the Axon, $nsat,ax$⁠.

Following Kuznetsov and Kuznetsov [10], we estimated the DCV concentration in the axon as 10% of the saturated DCV concentrations in type Ib boutons. This means that $nsat,ax$ is estimated as $0.1×40$ vesicles/ μm = 4 vesicles/ μm. We used the same value for $nsat,ax$ for type Ib and type III terminals.

#### Estimating a Value of $hax$⁠.

We used the following balance equation for DCVs entering bouton 4, to estimate a value of $hax$
$haxnsat,ax=jax→4$
(17)

We obtained that $hax=0.0167$μm/s.

#### Mass Transfer Coefficients Characterizing DCV Capture into the Resident State in Boutons, $h1$⁠, $h2a$⁠, $h2r$⁠, $h3a$⁠, $h3r$⁠, $h4a$⁠, $h4r$⁠; and Saturated Concentrations of DCVs in Boutons at Infinite DCV Half-Life or Half-Residence Time, $nsat0,1$⁠, $nsat0,2$⁠, $nsat0,3$⁠, $nsat0,4$⁠.

For type Ib boutons, we used $ϖ4=ϖ3=ϖ2=0.1$ because the DCV flux drops by ∼10% after passing each bouton (Wong et al. [5]). For the most distal bouton, we used $ϖ1=0.4$. This was necessary for the model to predict the initial DCV accumulation in the most distal type Ib bouton, the effect that was reported in Wong et al. [5].

For type III boutons, the situation is more complicated. For anterogradely moving DCVs, the capture rate in boutons is high, but when passing boutons 2, 3, and 4 in the retrograde direction, DCVs are captured at the same reduced capture rate observed in type Ib boutons, $ϖIb=0.1$ (Levitan [11]). One possible hypothesis explaining the different capture rates of anterogradely and retrogradely moving DCVs in type III boutons could rely on different affinity to microtubules of kinesin and dynein motors that transport DCVs. If kinesin motors transporting anterogradely moving DCVs are more frequently detached from microtubules than dynein motors, anterogradely moving DCVs must be more frequently captured in boutons. Another possible explanation is that kinesin motors more easily release DCVs than dynein motors. The process of DCV release by motors may be differentially controlled (Shakiryanova et al. [19], Bulgari et al. [3], Cavolo et al. [16]).

We wrote equations for the reduction of the DCV flux after passing boutons 4, 3, and 2, respectively, at t = 0 (at the moment when there are no resident DCVs in the bouton). These equations are applicable to both type Ib and type III boutons. We obtained the following:
$ϖ4jax→4=h4a(nsat0,4−0)$
(18a)

$ϖ3jax→4(1−ϖ4)=h3a(nsat0,3−0)$
(18b)

$ϖ2jax→4(1−ϖ4)(1−ϖ3)=h2a(nsat0,2−0)$
(18c)
At steady-state, the balance of captured and destroyed (or reentering the circulation) DCVs gives the following:
$hia(nsat0,i−nsat,i)+hir(nsat0,i−nsat,i)=Linsat,iln(2)T1/2, i=2,…,4$
(19)
For type III boutons, Eqs. (18) and (19) must be supplemented with equations relating the mass transfer coefficients characterizing DCV capture rates when DCVs move anterogradely through boutons 4, 3, and 2 and when DCVs move retrogradely through the same boutons. Since the DCV capture rates in retrograde motion through boutons 2, 3, and 4 are presumably reduced to that of type Ib boutons (Levitan [11]), we set
$h2r=ϖIbϖ2h2a, h3r=ϖIbϖ3h3a, h4r=ϖIbϖ4h4a$
(20a)
to account for the reduced capture efficiency of retrogradely moving DCVs. For type Ib boutons Eq. (20a) gives that
$h2r=h2a, h3r=h3a, h4r=h4a.$
(20b)
The exact solution of Eqs. (18)–(20) is given by Eqs. (A1)(A6) in the Appendix. We solved Eqs. (18)–(20) with $ϖ4=ϖ3=ϖ2=0.1$ for type Ib boutons and obtained the following:
$h4a=h4r=8.65×10−5, h3a=h3r=6.99×10−5, h2a=h2r=5.48×10−5 μm/s$
(21)

$nsat0,4=77.1, nsat0,3=85.9, nsat0,2=98.5 vesicles/μm$
(22)

In type III terminals the boutons are of different sizes. The largest bouton is randomly located. We simulated a representative morphology where the largest bouton is bouton 3 (bouton numbering is shown in Fig. 1). To prevent too much DCV accumulation in the most proximal bouton, we used $ϖ4=0.2$. For other boutons, we used $ϖ3=ϖ2=ϖ1=0.65$ (65% initial capture efficiency is based on the report by Bulgari et al. [15]), see also footnote “e” after Table 2.

We solved Eqs. (18)–(20) with $ϖ4=0.2$, $ϖ3=0.65$, and $ϖ2=0.65$ for type III boutons and obtained the following:
$h4a=6.42×10−5, h3a=9.27×10−6, h2a=1.67×10−5 μm/s$
(23)

$h4r=3.21×10−5, h3r=1.43×10−6, h2r=2.57×10−6 μm/s$
(24)

$nsat0,4=208, nsat0,3=3740, nsat0,2=727 vesicles/μm$
(25)
Writing DCV balances for bouton 1 at t = 0 and at steady-state, we obtained the following equations that replace Eqs. (18) and (19):
$ϖ1jax→4(1−ϖ4)(1−ϖ3)(1−ϖ2)=h1(nsat0,1−0)$
(26)

$h1(nsat0,1−nsat,1)=L1nsat,1ln(2)T1/2$
(27)

Equation (27) accounts for the fact that DCVs pass bouton 1 only once. As shown in Fig. 2, DCVs pass all other boutons twice, first as they move anterogradely and then as they move retrogradely.

The exact solution of Eqs. (26) and (27) is given by Eqs. (A7) and (A8) in the Appendix. Equations (26) and (27) were solved for type Ib boutons with $ϖ4=ϖ3=ϖ2=0.1$ and $ϖ1=0.4$. We obtained the following:
$h1=3.26×10−4μm/s, nsat0,1=59.7 vesicles/μm$
(28)
We then solved Eqs. (26) and (27) for type III boutons with $ϖ4=0.2$, $ϖ3=0.65$, $ϖ2=0.65$, and $ϖ1=0.65$ (using a = 3, b = 9, and c = 3) and obtained the following:
$h1=3.08×10−6μm/s, nsat0,1=1380 vesicles/μm$
(29)

### Calculating the Rate of DCV Production in the Soma.

Since at steady-state the rate of DCV production in the soma must be equal to the rate of DCV destruction due to a finite half-life of DCVs, $jsoma$ is calculated from the following equation:
$jsoma=(1−δ)(2L1nsat,1+2L2nsat,2+2L3nsat,3+2L4nsat,4)ln(2)T1/2+Laxnsat,axln(2)T1/2,ax$
(30)

If $δ=1$, then DCVs can only be destroyed in the axon (no destruction in boutons, since $δ=1$ implies that all captured DCVs eventually reenter the transiting pool). In Eq. (30), $jsoma$ is the net rate of DCV production, the rate of DCV synthesis minus the rate of DCV destruction in somatic lysosomes. According to Eq. (30), $jsoma$ depends on the fate of DCVs in boutons. Indeed, if DCVs are destroyed in boutons ($δ=0$), more DCVs need to be produced in the soma to compensate for the DCV loss.

### Numerical Solution.

In order to solve Eqs. (1)(14) numerically, we wrote a Matlab script (Matlab R2016b, MathWorks, Natick, MA). For the numerical solution of differential Eqs. (1)(5) Matlab's solver ode45 was utilized. The error tolerance parameters, RelTol and AbsTol, were set to 10–6 and 10–8, respectively. We checked that further decreases of RelTol and AbsTol did not change the solution. Because governing equations (Eqs. (1)(13)) are linear with respect to DCV concentrations, the solutions for a given set of parameters (Tables 2 and 3) are unique.

## Results

Figures that display mass transfer coefficients which characterize the rates of DCV capture (Fig. S1), saturated concentrations of DCVs at infinite DCV residence time (or infinite half-life) (Fig. S2), and concentrations of DCVs at steady-state (Fig. S3) are available under the “Supplemental Materials” tab for this paper on the ASME Digital Collection.

### The Case When All Captured DCVs Eventually Escape from Boutons and Reenter the Circulation (⁠$δ=1$⁠).

Figure 3 shows how concentrations of DCVs in boutons increase to their steady-state values. Concentrations in type Ib boutons increase to the same steady-state value in all boutons (Fig. 3(a)), while concentrations in type III boutons increase to different steady-state values depending on the size of a bouton (see also footnote “b” after Table 2). It takes ∼25 h for the concentrations in type Ib boutons to reach their steady-state values (Fig. 3(a)) while in type III boutons it takes ∼120 h.

The DCV concentration in the axon initially decreases for ∼10 h, but then increases as the DCV circulation develops and DCVs start to return from bouton 4 to the axon. The concentration in the axon reaches its steady-state value in about ∼300 h (Figs. 4(a) and 4(b)).

The DCV fluxes at steady-state in type Ib terminals are the same between all boutons and are equal to the flux of DCVs that enter the terminal (Fig. 5(a)). This is because for $δ=1$ at steady-state, the rate of capture into the resident state in a bouton is equal to the rate of escape from that bouton into the circulation. Thus, for $δ=1$, the model predicts a uniform DCV circulation at steady-state in type Ib terminals. However, the situation is more complicated in type III boutons because the anterograde fluxes decrease from bouton to bouton and the retrograde fluxes increase from bouton to bouton (Fig. 5(b)). This is because in type III boutons the mass transfer coefficients controlling the rate of DCV capture into the resident state are larger for anterogradely moving DCVs than for retrogradely moving DCVs (Eq. (20a), see also Fig. S1(b)). However, when DCVs escape from the resident state, we assumed that 50% of such DCVs join anterograde transport and 50% join retrograde transport ($ε=0.5$). Thus, although on average the number of captured DCVs at steady-state equals to the number of DCVs that reenter the transiting pool, more anterogradely moving DCVs are captured than reenter the transiting pool and more retrogradely moving DCVs reenter the transiting pool than are captured. This explains the trend in Fig. 5(b).

When DCV fluxes start after the release by exocytosis in type Ib boutons, they initially drop from bouton to bouton, but in ∼20 h they converge to the same curve (Fig. 6(a)). However, in type III boutons, the fluxes do not converge to the same curve (Fig. 6(b)). This agrees with the trend in Fig. 5(b), which is discussed earlier. It is interesting that fluxes reach their steady-state values slower than DCV concentrations in boutons (compare Figs. 3 and 6). This suggests that the time it takes for the DCV circulation to develop is larger than the time it takes for the DCV concentrations in boutons to become steady-state. The kink on the curves in Fig. 6(b) is explained by the fact that at ∼100 h the DCV concentrations in boutons reach their steady-state values, but the DCV concentration in the axon continues to change (Fig. 4(b)). The circulation becomes steady-state only after $nax$ reaches its steady-state value.

The largest DCV flux in Figs. 5 and 6 is ∼0.07 vesicles/s. This result can be used to estimate the number of DCVs in the transiting state as $jmaxLi/v$, where $jmax$ is the maximum DCV flux, $Li$ is the length of a bouton (5 μm), and v is the average DCV velocity (v is estimated to be ∼1 μm/s, Wong et al. [5], Kwinter et al. [4], de Jong et al. [20], Park et al. [21]). This estimates the number of DCVs transiting through a bouton as ∼0.35 vesicles, which is much smaller than the number of vesicles in the resident state (∼200 DCVs reside in a type Ib bouton at steady-state, Shakiryanova et al. [19]).

Retrograde fluxes in type Ib and type III boutons start 300 s after anterograde fluxes begin; this delay simulates switching of kinesin to dynein motors in bouton 1 (Figs. 7(a) and 7(b)). In a type III terminal, retrograde fluxes are delayed even more because of a high anterograde capture efficiency of type III boutons. As a result, initially most DCVs are captured before they can turn around in the most distal bouton and start moving retrogradely (Fig. 7(b)).

### The Case When Half of Captured DCVs Escape From Boutons and Return Back into the Circulation and the Rest Are Destroyed in Boutons (⁠$δ=0.5$⁠).

The dynamics of change of the DCV concentration in type Ib boutons is not significantly affected by a decrease of $δ$ (compare Figs. 8(a) and 3(a)). This is explained by a small DCV capture rate in type Ib boutons. However, in type III boutons the DCV concentration increases faster, if $δ$ is decreased (compare Figs. 8(b) and 3(b)). This is because a decrease of $δ$ results in less DCVs escaping from boutons and more DCVs destroyed in boutons. Hence, to replace those DCVs that are being destroyed in boutons at steady-state, increased DCV synthesis in the soma is required. For example, in a neuron with a type III terminal $jsoma=0.0128$ vesicles/s for $δ=1$ and $jsoma=0.0795$ vesicles/s for $δ=0.5$ (Table 3). The larger supply of DCVs from the soma results in boutons being filled faster. Our model assumes that $jsoma$ remains constant as DCV stores in the terminal are being replenished after the release by exocytosis. It is possible that the rate of DCV synthesis in the soma is larger in the beginning of the process, but there is not enough experimental data at this point to quantify this process.

A larger rate of DCV synthesis in the soma results in a lesser depletion of the axon of DCVs. In the case of a type III terminal, for $δ=1$, the DCV concentration in the axon drops to ∼0.8 vesicles/ μm after ∼10 h (Fig. 4(b)), while for $δ=0.5$, it drops to only ∼2.4 vesicles/ μm (Fig. 9(b)).

At steady-state, the fluxes between the boutons in a type Ib terminal now decay from one bouton to the next (Fig. 10(a)), because for $δ=0.5$, there is some destruction of DCVs in boutons (compare with Fig. 5(a) for $δ=1$). Now, the fluxes between the boutons do not collapse to the same curve (Figs. 11(a) and 11(b)), unlike what happens for $δ=1$ (Figs. 6(a) and 6(b)). Due to a larger rate of DCV synthesis in the soma, the retrograde flux $j1→2$ is now slightly larger than for $δ=1$ (compare Fig. 12(b) with Fig. 7(b)).

### The Case When All Captured DCVs Are Eventually Destroyed in Boutons (⁠$δ=0$⁠).

This case was investigated in Kuznetsov and Kuznetsov [10] (although for different capture efficiencies in type III boutons); therefore, we present the results for this case in the Supplementary Material, which is available under the “Supplemental Materials” tab for this paper on the ASME Digital Collection. In summary, the dynamics of DCV concentration increase in type Ib boutons for $δ=0$ is not much different from that for $δ=0.5$ and $δ=1$ (compare Fig. S4a with Figs. 3(a) and 8(a)). This is due to a small capture rate in these boutons. However, in a type III terminal DCV concentrations in boutons increase faster, which is explained by a larger rate of DCV synthesis in the soma (Table 3). The depletion of the axon of DCVs, which occurs in type Ib boutons in the beginning of the process, is even smaller than in the previous cases (Fig. S5(a)). In type III boutons the concentration in the axon remains constant because in this case the DCV circulation in the terminal never develops (Fig. S5(b), see also Bulgari et al. [15]). DCV fluxes at steady-state now decrease faster from bouton to bouton because all captured DCVs are eventually destroyed in boutons (Figs. S6(a) and S6(b)). In type III terminals, no DCVs return from bouton 4 to the axon. All DCVs that enter the terminal are captured and destroyed in boutons (Fig. S7(b)). A larger rate of DCV synthesis in the soma results in larger retrograde fluxes (now not only $j1→2$ but also $j2→3$ take nonzero values immediately after 300 s into the process, Fig. S8(b).

Table 4 summarizes the percentage of DCVs that re-enter the axon at steady-state for three different assumptions concerning the fate of DCVs in boutons. Experimentally, it should be much easier to measure the percentage of DCVs that return to the axon than determine the fate of DCVs in boutons. The model can, thus, be used to obtain information about the fate of DCVs in boutons.

## Discussion and Future Directions

Utilizing the conservation of DCVs, we analyzed DCV transport in nerve terminals for three different assumptions concerning the fate of DCVs in boutons: (i) all captured DCVs eventually reenter the circulation, (ii) half of captured DCVs reenter the circulation and half are destroyed, (iii) all captured DCVs are eventually destroyed in boutons. If all DCVs reenter the circulation after spending some time in the resident state in boutons, at steady-state, the flux of DCVs entering the terminal must be equal to the flux of DCVs leaving the terminal and returning to the axon. This applies to both type Ib and type III terminals. This means that a DCV circulation would exist not only in type Ib terminals but also in type III terminals. The absence of a DCV circulation in type III boutons, reported in Bulgari et al. [15], could be explained by the assumption that type III boutons studied in Bulgari et al. [15] were not at steady-state but were rather accumulating DCVs between the molts for later release by exocytosis (Levitan [11]).

If half of DCVs are destroyed in boutons and half reenter the transiting pool, then at steady-state, there are DCV circulations in both type Ib and type III terminals, but the circulation in a type III terminal is weaker. If all captured DCVs are eventually destroyed in the boutons, then there is circulation in a type Ib terminal, but there is no circulation in a type III terminal, because all DCVs are captured before they leave the terminal.

The model can be easily extended by incorporating different values for the DCV half-life and half-residence time in boutons (currently they were both set to 6 h). Future research should address molecular mechanisms of different capture rates in type III boutons for anterogradely and retrogradely moving DCVs, reported in Levitan [11]. Future research should also address developing a two-concentration model, which would simulate not only the concentration of captured DCVs in boutons but also the concentration of DCVs transiting through the boutons. Statistical physics models for parameters ε and δ should also be developed. Feedback effects could be incorporated by assuming that the net rate of DCV production in the soma, $jsoma$, depends on signals received from boutons. Extending the model to simulate DCV transport under various neurodegenerative conditions associated with DCV transport failure would be worthwhile. Incorporating stochastic nature of DCV transport into future models would be a promising extension of the current approach.

## Funding Data

• Alexander von Humboldt-Stiftung (Humboldt Research Award).

• National Science Foundation, Division of Chemical, Bioengineering, Environmental, and Transport Systems (Grant No. CBET-1642262).

### Appendix

The exact solution of Eqs. (18)–(20) is
$h4a=jax→4ϖ4nsat,4−L4ϖ4ln(2)T1/2(ϖ4+ϖIb)$
(A1)

$nsat0,4=jax→4nsat,4T1/2(ϖ4+ϖIb)jax→4T1/2(ϖ4+ϖIb)−L4nsat,4ln(2)$
(A2)

$h3a=−jax→4ϖ3(ϖ4−1)nsat,3−L3ϖ3ln(2)T1/2(ϖ3+ϖIb)$
(A3)

$nsat0,3=jax→4nsat,3T1/2(ϖ4−1)(ϖ3+ϖIb)jax→4T1/2(ϖ4−1)(ϖ3+ϖIb)+L3nsat,3ln(2)$
(A4)

$h2a=jax→4ϖ2(ϖ3−1)(ϖ4−1)nsat,2−L2ϖ2ln(2)T1/2(ϖ2+ϖIb)$
(A5)

$nsat0,2=jax→4nsat,2T1/2(ϖ3−1)(ϖ4−1)(ϖ2+ϖIb)jax→4T1/2(ϖ3−1)(ϖ4−1)(ϖ3+ϖIb)−L2nsat,2ln(2)$
(A6)
The exact solution of Eqs. (26) and (27) is
$h1=−jax→4ϖ1(ϖ2−1)(ϖ3−1)(ϖ4−1)nsat,1−L1ln(2)T1/2$
(A7)

$nsat0,1=jax→4nsat,1T1/2ϖ1(ϖ2−1)(ϖ3−1)(ϖ4−1)jax→4T1/2ϖ1(ϖ2−1)(ϖ3−1)(ϖ4−1)+L1nsat,1ln(2).$
(A8)

## References

References
1.
Wong
,
M. Y.
,
Cavolo
,
S. L.
, and
Levitan
,
E. S.
,
2015
, “
Synaptic Neuropeptide Release by Dynamin-Dependent Partial Release From Circulating Vesicles
,”
Mol. Biol. Cell
,
26
(
13
), pp.
2466
2474
.
2.
Li
,
L.
,
Tian
,
X.
,
Zhu
,
M.
,
Bulgari
,
D.
,
Böhme
,
M. A.
,
Goettfert
,
F.
,
Wichmann
,
C.
,
Sigrist
,
S. J.
,
Levitan
,
E. S.
, and
Wu
,
C.
,
2014
, “
Drosophila Syd-1, Liprin-Alpha, and Protein Phosphatase 2A B' Subunit Wrd Function in a Linear Pathway to Prevent Ectopic Accumulation of Synaptic Materials in Distal Axons
,”
J. Neurosci.
,
34
(
25
), pp.
8474
8487
.
3.
Bulgari
,
D.
,
Deitcher
,
D. L.
, and
Levitan
,
E. S.
,
2017
, “
Loss of Huntingtin Stimulates Capture of Retrograde Dense-Core Vesicles to Increase Synaptic Neuropeptide Stores
,”
Eur. J. Cell Biol.
,
96
(
5
), pp.
402
406
.
4.
Kwinter
,
D. M.
,
Lo
,
K.
,
Mafi
,
P.
, and
Silverman
,
M. A.
,
2009
, “
Dynactin Regulates Bidirectional Transport of Dense-Core Vesicles in the Axon and Dendrites of Cultured Hippocampal Neurons
,”
Neurosci.
,
162
(
4
), pp.
1001
1010
.
5.
Wong
,
M. Y.
,
Zhou
,
C.
,
Shakiryanova
,
D.
,
Lloyd
,
T. E.
,
Deitcher
,
D. L.
, and
Levitan
,
E. S.
,
2012
, “
Neuropeptide Delivery to Synapses by Long-Range Vesicle Circulation and Sporadic Capture
,”
Cell
,
148
(
5
), pp.
1029
1038
.
6.
Bressloff
,
P. C.
, and
Levien
,
E.
,
2015
, “
Synaptic Democracy and Vesicular Transport in Axons
,”
Phys. Rev. Lett.
,
114
(
16
), p.
168101
.
7.
Bressloff
,
P. C.
, and
Karamched
,
B. R.
,
2016
, “
Model of Reversible Vesicular Transport With Exclusion
,”
J. Phys. A-Math. Theor.
,
49
(
34
), p.
345602
.
8.
Kuznetsov
,
I. A.
, and
Kuznetsov
,
A. V.
,
2013
, “
A Compartmental Model of Neuropeptide Circulation and Capture Between the Axon Soma and Nerve Terminals
,”
Int. J. Numer. Methods Biomed. Eng.
,
29
(
5
), pp.
574
585
.
9.
Kuznetsov
,
I. A.
, and
Kuznetsov
,
A. V.
,
2015
, “
Modeling Neuropeptide Transport in Various Types of Nerve Terminals Containing En Passant Boutons
,”
Math. Biosci.
,
261
, pp.
27
36
.
10.
Kuznetsov
,
I. A.
, and
Kuznetsov
,
A. V.
,
2015
, “
Can an Increase in Neuropeptide Production in the Soma Lead to DCV Circulation in Axon Terminals With Type III En Passant Boutons?
,”
Math. Biosci.
,
267
, pp.
61
78
.
11.
Levitan
,
E. S.
,
2017
, personal communication.
12.
Atwood
,
H. L.
,
Govind
,
C. K.
, and
Wu
,
C.-F.
,
1993
, “
Differential Ultrastructure of Synaptic Terminals on Ventral Longitudinal Abdominal Muscles in Drosophila Larvae
,”
J. Neurobiol.
,
24
(
8
), pp.
1008
1024
.
13.
Jia
,
X.-X.
,
Gorczyca
,
M.
, and
Budnik
,
V.
,
1993
, “
Ultrastructure of Neuromuscular-Junctions in Drosophila: Comparison of Wild-Type and Mutants With Increased Excitability
,”
J. Neurobiol.
,
24
(
8
), pp.
1025
1044
.
14.
Menon
,
K. P.
,
Carrillo
,
R. A.
, and
Zinn
,
K.
,
2013
, “
Development and Plasticity of the Drosophila Larval Neuromuscular Junction
,”
Wiley Interdiscip. Rev.-Dev. Biol.
,
2
(
5
), pp.
647
670
.
15.
Bulgari
,
D.
,
Zhou
,
C.
,
Hewes
,
R. S.
,
Deitcher
,
D. L.
, and
Levitan
,
E. S.
,
2014
, “
Vesicle Capture, Not Delivery, Scales Up Neuropeptide Storage in Neuroendocrine Terminals
,”
Proc. Natl. Acad. Sci. U. S. A.
,
111
(
9
), pp.
3597
3601
.
16.
Cavolo
,
S. L.
,
Bulgari
,
D.
,
Deitcher
,
D. L.
, and
Levitan
,
E. S.
,
2016
, “
Activity Induces Fmr1-Sensitive Synaptic Capture of Anterograde Circulating Neuropeptide Vesicles
,”
J. Neurosci.
,
36
(
46
), pp.
11781
11787
.
17.
Jacquez
,
J. A.
,
1985
,
Compartmental Analysis in Biology and Medicine
,
2nd ed.
,
The University of Michigan Press
, Ann Arbor, MI.
18.
Rizzoli
,
S. O.
,
2014
, “
Synaptic Vesicle Recycling: Steps and Principles
,”
EMBO J.
,
33
(
8
), pp.
788
822
.
19.
Shakiryanova
,
D.
,
Tully
,
A.
, and
Levitan
,
E. S.
,
2006
, “
Activity-Dependent Synaptic Capture of Transiting Peptidergic Vesicles
,”
Nat. Neurosci.
,
9
(
7
), pp.
896
900
.
20.
de Jong
,
E. K.
,
Vinet
,
J.
,
Stanulovic
,
V. S.
,
Meijer
,
M.
,
Wesseling
,
E.
,
Sjollema
,
K.
,
Boddeke
,
H. W. G. M.
, and
Biber
,
K.
,
2008
, “
Expression, Transport, and Axonal Sorting of Neuronal CCL21 in Large Dense-Core Vesicles
,”
FASEB J.
,
22
(
12
), pp.
4136
4145
.
21.
Park
,
J. J.
,
Cawley
,
N. X.
, and
Loh
,
Y. P.
,
2008
, “
A Bi-Directional Carboxypeptidase E-Driven Transport Mechanism Controls BDNF Vesicle Homeostasis in Hippocampal Neurons
,”
Mol. Cell. Neurosci.
,
39
(
1
), pp.
63
73
.