The fast and the slow axonal transport: a unified approach based on cargo and molecular motors coupled dynamics
TThe fast and the slow axonal transport: a uni(cid:133)ed approachbased on cargo and molecular motors coupled dynamics.
Alexandre Y. C. Cho, Victor R. C. M. Roque, Carla GoldmanUniversidade de Sao Paulo, Instituto de FisicaRua do Matao, 137105508-090 Sao Paulo, SP, BrazilMay 2020
Abstract
The origins of the large di⁄erences observed to the rates with which the diverse particlesare conveyed along axonal microtubules are still a matter of debate in the literature. Thereis evidence that certain neurodegenerative diseases may be triggered by disturbances to therelated transport processes. Motivated by this, we employ a model to investigate the mobilityproperties of certain cargoes which dynamics are coupled with that of molecular motors oncrowded microtubules. For certain initial and boundary conditions, we use the method ofcharacteristics to resolve perturbatively the pair of equations of the Burgers type resultingfrom a mean-(cid:133)eld approach to the original microscopic stochastic model. Extensions to thenon-perturbative limits are explored numerically. In this context, we were able to (cid:133)gure outconditions under which cargos average velocities may di⁄er up to orders of magnitude justby changing the number of motors on the considered track. We then discuss possibilities toconnect these theoretical predictions with available experimental data about axon transport. keywords - fast and slow axonal transport; molecular motors; tra¢ c jam; ASEP models. Introduction
The diverse types of particles that are usually transported along axons can be grouped into twomajor categories characterized by the average speed v at which this transport is put into e⁄ect.The group of fast-moving particles with v (cid:24) : (cid:0) (cid:22)m=s comprises membranous organelles as theGolgy-derived vesicles, mitochondria, endosomes and lysosomes, among others. The groups of slow-moving particles with v (cid:24) : (cid:0) : (cid:22)m=s for slow component- a (Sc a ) and v (cid:24) : (cid:0) : (cid:22)m=s for slow component- b (Sc b ) comprise nonmembranous neuro(cid:133)laments, cytoskeletal and cytosolicproteins, among others.The general biological interest of related studies relies on the fact that certain neurologicaldiseases are believed to be directly associated with failures of the system to keep the transport atthe right rates leading, eventually, to local particle accumulation [1]. The review of S. Roy [2]o⁄ers an historical account of experimental achievements in this (cid:133)eld and of the major hypothesesin the literature to explain the huge di⁄erences in the rates of particle transport along axons. Oneof the models used in these studies, the Dynamic Recruitment model proposed to explain the slowmovement, requires the presence of a common carrier structure [3],[4] or even fast moving cargovesicles [5] onto which cytosolic proteins would transiently get attached to. The slow rate withwhich these proteins are conveyed along the axon would then be attributed to the relative largeintervals of time along which they stay detached from the carrier [6]. A stochastic model accountingfor these ideas is considered in [7]. The other major hypothesis, the
Stop and Go model , attributesthe causes of slow movement to the ability of neuro(cid:133)laments, for example, to pause during transit.The duration of such pauses would be determinant to the resulting slow rates, although their originsand eventual mechanisms of control are not clari(cid:133)ed yet.Evidences indicating that individual Cytoskeletal polymers conveyed by slow transport in axonscan also move as fast as membranous organelles suggest that both fast and slow transport in axonsmay be due to a single mechanism. Actually, a uni(cid:133)ed view of axonal transport had already beenproposed by Ochs [8] according to which the only actual directed movement would be that of fastcargoes, the detected slow movement would be due to local and casual rearrangements of particles,not resulting in long-range transport. In line with this it is argued these days that such uni(cid:133)edmechanism is supported by the dynamics of molecular motors and on their ability to transport avariety of particles, referred generically as cargoes, as they move on structured "tracks", i.e. alongaxonal microtubules [9]. Yet, one still may question about the speci(cid:133)c motors that would drive2argoes at such slow rates as shown by Sc a and Sc b [1]. In view of this, we believe that althoughthe description based on molecular motors seems promising, it requires to be complemented bymechanisms that regulate cargo-motor interactions, able to explain such huge di⁄erences betweenthe rates of slow and fast transport in the absence of a putative carrier. Our intention here is tocontribute in this regard.To that end, we investigate speci(cid:133)c properties of a stochastic model proposed elsewhere [10],[11], [12]. It is conceived on the idea that the transport of cargoes by molecular motors is basedon a mechanism of cargo hopping among motors which, in turn, depends crucially on the quantityof motor available on the track. Here, we seek for quantitative predictions for the cargo averagevelocities under these conditions. We explore situations under which motor jamming on the trackmay promote cargo transport.We should remark that there are always two main aspects concerning the dependence of cargotransport on molecular motors in such systems. One of these is related to the microscopic mech-anisms of motor-cargo binding a¢ nity. The other concerns the dependence of the transport onthe number of available motors on track, not only on the type of motor directly bound to eachcargo [13], [14], [15]. In order to investigate the dependence on motor occupation, we base ourquantitative analysis on a model that describes the microscopic dynamics of motors and cargoescoupled through stochastic processes of ASEP type (Asymmetric Simple Exclusion Process). Themodel allows for two microscopic mechanisms for cargo movement. One of these describes themovement of cargo while bound to a motor; the other describes cargo hopping between pairs ofneighbor motors on a track. There is only one type of motor with bias in a de(cid:133)nite direction;without prejudice for the analysis, this is chosen anterograde. We have already studied possibilitiesto observe bidirectional movement in similar model systems [10],[11],[12]. The main concern of thepresent study is related to quantitative estimations for the average velocities developed by cargoessubject to these two mechanisms. In this context, we show the conditions under which the averagecargo speed may vary up to orders of magnitude simply adjusting the amount of motor available onaxons. This opens a new possibility to understand the observed di⁄erences in the rates of transportdepicted by real axons in quantitative terms and under a uni(cid:133)ed perspective.In Section 2 we derive the mean-(cid:133)eld equations for the model. The resulting non-linear coupledequations of Burgers type is studied analytically in Section 3 for certain initial and boundaryconditions using the Method of Characteristics. To our knowledge, such solutions to this particularsystem of equations represent an original contribution. As we shall show, these are achieved under3 perturbative approximation. The analytical expressions obtained in this way pointed out to therelevant regions of parameters, and particular initial and boundary values illustrating the role ofeach mechanism in producing such di⁄erences to the average cargo velocities, as observed. Formore general conditions, away from the perturbative limits, we base our analysis on a numericalstudy to resolve the particular system of non-linear equations of interest [16]. This is explained inSection 4. We discuss our results in Section 5. Some conclusions and additional remarks are inSection 6. The kind of ASEP model we shall consider to describe the stochastic dynamics of molecular motorsand cargoes has already been proposed elsewhere to study the so called bidirectional movement[10],[11],[12]. The model version analyzed here embraces both mechanisms of cargo transport atthe microscopic level namely, the isolated motion of a cargo attached to a motor and also themechanism of cargo hopping from motor to motor. The idea is to circumvent di¢ culties for long-range cargo transport due to motor jamming (congestion) on axons [17].Accordingly, an axonaltrack is represented by a one-dimensional lattice which sites may be occupied at each time eitherby an unloaded motor or by a loaded motor or else it can be empty. The model is de(cid:133)ned by thedynamics of occupation, as follows: ( a ) 10 ! with probability p ( b ) 20 ! with probability q ( c ) 12 ! with probability f ( d ) 21 ! with probability g (1)Label is assigned to an empty site. Label is assigned to a site that is occupied by a motor carryingno cargo - an unloaded motor. Label is assigned to a site occupied by a motor attached to a cargo -a loaded motor. The pairs shown in (1) indicate that each of the possible processes is represented bythe occupation of two neighboring sites. Parameters p; q; f; g are the probabilities for the occurrenceof each process per time interval ((cid:1) t ) : Process (a) represents an elementary movement of a biasedmotor from cell body towards axon·s terminal (anterograde direction). Process (b) represents anelementary process of the kind considered in (a) of a motor when carrying a cargo particle. Theremaining processes describe the exchanging of cargoes between neighboring motors, either to the4eft (c) or to the right (d). A condition for processes (a) and (b) to take place is that the site at theright of the motor remains empty during the time interval (cid:1) t . The two processes of exchangingdepend on the presence of a free motor at the right (c) or at the left (d) of the one attached to thecargo. In order to obtain the mean-(cid:133)eld equations corresponding to processes (1 a-d) for continuous-timeand space intervals, we follow the procedure described in [18] that nevertheless accounts only forone type of particle. In order to extend it to take into account the two types of particles present inour model, we de(cid:133)ne the vector: P t [ n ; n ; n ; :::n i ; n i +1 ; :::; n i + k ; :::; n N ] : (2)For any j = 1 ; ::N; n j = f ; ; g ; it represents the probability of any con(cid:133)guration of particles andholes occupying the sites of the one-dimensional lattice at time t . The marginal probabilities arede(cid:133)ned accordingly: P ti [ n i ; n i +1 ; :::; n i + k ] = X f n j g (1 (cid:20) j
The system (9) can be expressed in a vectorial form, @U ( z; t ) @t + A ( U ( z; t )) @U ( z; t ) @z = 0 : (11)where the vector U and the matrix A are de(cid:133)ned by: U = u ( z; t ) c ( z; t ) (12) A = A ( U ) = (cid:11) (1 (cid:0) u (cid:0) (cid:18)c ) (cid:0) (cid:11)(cid:18)u (cid:0) (cid:12)(cid:17)c (cid:12) (1 (cid:0) c (cid:0) (cid:17)u ) : (13)In order to proceed into a quantitative analysis of the equation (11), we examine the case f = g; for which (cid:18) = (cid:17) = 1 : In addition, we parameterize the constants (cid:11) and (cid:12) related to the hoppingprobabilities p and q , as (cid:11) = (cid:13) + " and (cid:12) = (cid:13) (cid:0) " (14)We may then express the eigenvalues (cid:21) and (cid:21) of A in terms of constants (cid:13) and " which become,for "=(cid:13) (cid:28) ; (cid:21) ’ (1 (cid:0) u (cid:0) c ) (cid:13) (cid:20) "(cid:13) (cid:18) u (cid:0) cu + c (cid:19)(cid:21) (cid:21) ’ (1 (cid:0) u (cid:0) c ) (cid:13) (cid:20) (cid:0) "(cid:13) (cid:18) u (cid:0) cu + c (cid:19)(cid:21) (15)Since (cid:21) and (cid:21) are distinct of each other and are both real, the system (9) is strictly hyperbolic. (cid:21) and (cid:21) are related to the velocities of the travelling wave solutions to Eq. (11). The correspondingleft eigenvectors (cid:0)! q and (cid:0)! q are calculated perturbatively up to (cid:133)rst order in the small parameter "=(cid:13); resulting: (cid:0)! q ’ (cid:18) (cid:0) ; (cid:20) (cid:0) (cid:0) "(cid:13) (cid:18) (cid:0) u (cid:0) cu + c (cid:19)(cid:21) (cid:19) (cid:0)! q ’ (cid:18) ; (cid:20) (cid:0) uc + 2 "(cid:13) uc (cid:18) (cid:0) u (cid:0) cu + c (cid:19)(cid:21) (cid:19) (16)The Riemann Invariants R ( U ) and R ( U ) are scalar quantities satisfying (cid:0)!r R i ( U ) = (cid:0)! q i ( U ) ;i = f ; g [19], [20]. These quantities are conserved by the dynamics along the set of corresponding7haracteristic curves z ( t ) and z ( t ) de(cid:133)ned by @ t z = (cid:21) and @ t z = (cid:21) . That is (cid:21) and (cid:21) are theinstantaneous velocities at each point of the characteristics on the ( z; t ) plane. In fact, q i (cid:1) (cid:18) @U@t + (cid:21) i @U@z (cid:19) = (cid:16) (cid:0)!r R i ( U ) (cid:17) (cid:1) (cid:18) @U@t + (cid:21) i @U@z (cid:19) = (cid:18) @R i @t + (cid:21) i @R i @z (cid:19) (cid:17) dR i d(cid:30) i = 0 (17)where we have used (11). With regard to expressions (15), one observes that if both u and c are constants along the characteristic curves, then both (cid:21) and (cid:21) are also constants so that thecharacteristics are straight lines. As indicated in the de(cid:133)nition above, the derivative dd(cid:30) i is takenalong each of the corresponding characteristic curves. Therefore, in order to (cid:133)nd R i ; i = f ; g theexpressions in (16) shall be inserted into the LHS of the (cid:133)rst equality in the expression above inorder to (cid:133)nd a general solution R i by integrating dR i d(cid:30) i = 0 : For q ; it results: (cid:0) (cid:18) @u@t + (cid:21) @u@z (cid:19) (cid:0) (cid:18) @c@t + (cid:21) @c@z (cid:19) (cid:0) "(cid:13) (cid:0) u (cid:0) cu + c (cid:18) @c@t + (cid:21) @c@z (cid:19) = 0 (18)which can be rewritten as dd(cid:30) (cid:20) (1 (cid:0) u (cid:0) c ) (cid:0) ln(1 (cid:0) u (cid:0) c ) + 2 "(cid:13) c (cid:21) = 0 (19)where the derivative is taken along the characteristics: dd(cid:30) = @@t + (cid:21) @@z : This means that R = (1 (cid:0) u (cid:0) c ) (cid:0) ln(1 (cid:0) u (cid:0) c ) + 2 "(cid:13) c (20)is the conserved quantity along the set of characteristic curves associated to q ( u; c ) .For q ; it results: c (cid:18) @ ( u=c ) @t (cid:0) (cid:21) @ ( u=c ) @z (cid:19) + 2 "(cid:13) uc (cid:0) u + c ) u + c (cid:18) @c@t + (cid:21) @c@z (cid:19) = 0 : (21)In order to write this in a form analogous to (19) we examine the region of densities for which: u ( z; t ) = 1 = x ( z; t ) c ( z; t ) = 1 = y ( z; t ) (22)with j x ( z; t ) j ; j y ( z; t ) j << = representing respectively, the excess (or depletion) of motors andcargoes at each point of the axon with respect to a considered density "background" set at = .With these, Eq. (21) can be approximated by dd(cid:30) (cid:18) x + 1 = y + 1 = (cid:19) ’ (23)8here dd(cid:30) (cid:17) @@t + (cid:21) @@z . This yields: R ’ (cid:18) x + 1 = y + 1 = (cid:19) : (24)For consistency, the expressions for the remaining quantities which are relevant to the analysis thatfollows must be reviewed in this region of densities: R ’ K + x + (1 + 2 "(cid:13) ) y (25)where K (cid:17) ln 2 + 12 + 14 (cid:18) "(cid:13) (cid:19) is a constant. For the eigenvalues, (cid:21) = (cid:21) ( x; y ) ’ (cid:0) x + y ) (cid:13) [1 + 2 "(cid:13) ( x (cid:0) y )] (cid:21) = (cid:21) ( x; y ) ’ ( (cid:0) ( x + y )) (cid:13) [1 (cid:0) "(cid:13) ( x (cid:0) y )] : (26)Observe that although the results (25) for R and (24) for R hold in limit "(cid:13) << only if both j x ( z; t ) j ; j y ( z; t ) j << = ; the expressions (26) for (cid:21) and (cid:21) hold valid for any x and y that keepthe densities u ( z; t ) and c ( z; t ) (22) within the interval [0 ; . We may now use (24), (25) and (26)to examine the behavior of both x ( z; t ) and y ( z; t ) through the analysis of the characteristic curvesfor certain initial and boundary conditions. We intend to examine the behavior of solutions c ( z; t ) and u ( z; t ) for all z (cid:21) and t (cid:21) satisfyingthe pair of equations (9), for initial conditions (IC) ( z > t = 0) : c ( z;
0) = 1 = u ( z;
0) = 1 = h (27)meaning that x ( z;
0) = h and y ( z;
0) = 0 ; and boundary conditions (BC) at z (cid:20) t (cid:21) ;c ( z (cid:20) ; t ) = 1 = y r i.e. y ( z (cid:20) ; t ) = y r u ( z (cid:20) ; t ) = 1 = x r i.e. x ( z (cid:20) ; t ) = x r : (28)9he constants x r ; y r ; h are such that j x r j ; j y r j ; j h j << = . Here, h represents the excess ( h > or depletion ( h < of motor density with respect to the value = along the axon at the initialtime t = 0 . The quantities x r and y r represent, respectively, the excess or depletion of motors andcargoes densities in the reservoir at z (cid:20) and t (cid:21) .In the following we present a quantitative analysis of the behavior of the two functions u ( z; t ) and c ( z; t ) with respect to time at each point of the domain ( z (cid:21) using the method of char-acteristics. We shall restrict this analysis, however, to certain initial and boundary values suchthat shocks and rarefaction of characteristic lines related to each of the eigenvalues (cid:21) or (cid:21) wouldnot be relevant up to a (cid:133)rst approximation in the small parameter "=(cid:13) to (cid:133)gure out solutionswithin considered time and space domains. As we shall see, this approximation is supported bynumerical data. Cases for which it does not hold but may be interesting concerning the trans-port phenomena of interest here are treated numerically and considered afterward for qualita-tive analysis as an extension of the perturbative results. A quantitative analysis of the equa-tions (9) for more general initial and boundary values is not on the scope of the present work. > At t = 0 and z (cid:21) ;(cid:21) ( x ( z; ; y ( z; (cid:21) ( h; ’ (cid:0) h(cid:13) (cid:20) "(cid:13) h (cid:21) < a ) (cid:21) ( x ( z; ; y ( z; (cid:21) ( h; ’ (cid:18) (cid:0) h (cid:19) (cid:13) (cid:20) (cid:0) "(cid:13) h (cid:21) > b ) (29)We shall consider "(cid:13) h << : This means that the set of characteristic lines de(cid:133)ned as C ( I ) alongwhich R in (25) is conserved have negative slope while the characteristics from the set C ( II ) alongwhich R in (24) is conserved have positive slope. Figure 1 suggests that these features may beextended consistently to the solutions x ( z; t ) and y ( z; t ) for all t (cid:21) : In fact, notice (cid:133)rst that thecharacteristic z ( t ) = (cid:21) ( h; t from the set C ( II ) emerging from z = 0 divides the plane into tworegions. Any point within the region z ( t ) > (cid:21) ( h; t referred to as Quiet Region (QR) , is reachedby characteristics that belong to both sets C ( I ) and C ( II ) emerging at t = 0 from all points z (cid:21) . The points within the region z ( t ) < (cid:21) t referred as Region 1 (Z1), are also reached from the setof characteristics C ( I ) emerging at t = 0 from the points z (cid:21) . However, the curves from the set10 ( II ) that reach the points within Z1 are those emerging from the boundary z = 0 at t > : Thus,while the quantities x ( z; t ) and y ( z; t ) are resolved within the QR by extending back along bothsets of characteristics C ( I ) and C ( II ) from any point ( z; t ) up to the initial conditions, it happensthat within Z1 the corresponding quantities x ( z; t ) and y ( z; t ) are resolved extending back fromboth the initial con(cid:133)guration through C ( I ) and from the boundary through C ( II ) : Consider then any point Q = ( z Q ; t Q ) at the crossing of two characteristics within the QR , onefrom the set C ( I ) and the other from the set C ( II ) . Since the initial conditions are such that bothdensities are constants for z > , one may write: R ( z;
0) = K + h = K + x Q + (1 + 2 "(cid:13) ) y Q R ( z;
0) = 1 + 4 h = 1 + 4 x Q y Q (30)where we have de(cid:133)ned x Q (cid:17) x ( z Q ; t Q ) and y Q (cid:17) y ( z Q ; t Q ) : From these two equations, one concludesthat at any point Q within QR , x Q = h and y Q = 0 : This means that the characteristics are straightlines within this region with slopes given in (29).Consider now the region Z1 de(cid:133)ned by the points ( z; t ) such that z (cid:20) (cid:21) ( h; t , t > . Figure1 shows two lines that belong to the set C ( I ) crossing a line from the set C ( II ) at points R and S within Z1 and at points P and Q within the QR . It follows that R ( z R ; t R ) = R ( z S ; t S ) R ( z P ; t P ) = R ( z R ; t R ) R ( z Q ; t Q ) = R ( z S ; t S ) (31)Because y Q = y P and x Q = x P ( QR ), one concludes that ( y R (cid:0) y S ) (1 + "(cid:13) ) = 2( x S y R (cid:0) x R y S ) (32)which has the trivial solutions y R = y S if x R = x S ; implying that the characteristics from C ( II ) are straight lines within the region Z1 . Now, in order to resolve for the quantities x (cid:3) + (cid:17) x ( z (cid:3) ; t (cid:3) ) and y (cid:3) + (cid:17) y ( z (cid:3) ; t (cid:3) ) at any point ( z (cid:3) ; t (cid:3) ) of Z1 , consider that11 (0 ; t ) = R ( z; R ( z (cid:3) ; t (cid:3) ) = R (0 ; t ) R ( z (cid:3) ; t (cid:3) ) = R ( z; (33)Using expressions (24) and (25), one (cid:133)nds that, for the chosen initial and boundary conditions andfor h > , the excess or depletion of cargoes y (cid:3) + (cid:17) y ( z (cid:3) ; t (cid:3) ) and motors x (cid:3) + (cid:17) x ( z (cid:3) ; t (cid:3) ) at any point ( z (cid:3) ; t (cid:3) ) within the Z1 are given by: y (cid:3) + = y r x (cid:3) + = h (cid:0) (cid:18) "(cid:13) (cid:19) y r (34)With regard to the characteristics C ( I ) within Z1 , one must examine the set of equations relatingboth quantities R ( z; t ) and R ( z; t ) evaluated at any two points R and T , as shown: R ( z R ; t R ) = R ( z T ; t T ) = R ( z; R ( z R ; t R ) = R ( z T ; t T ) = R (0 ; t ) (35)The last equality follows from the particular choice of constant boundary conditions as in (28). Fromthese, one may conclude that x R = x T and y R = y T meaning that both x and y are individuallyconserved along the characteristics C ( I ) : In turn, this implies that (cid:21) is constant and thus thecharacteristics of C ( I ) are also straight lines within Z1 . From (34), one concludes that for y r < h (cid:18) "(cid:13) (cid:19) (36)the condition x (cid:3) + + y (cid:3) + > holds assuring that (cid:21) ( x (cid:3) + ; y (cid:3) + ) < within Z1, and also that for y r > h (cid:18) "(cid:13) (cid:19) (cid:18) (cid:0) h (cid:19) (37)condition = (cid:0) ( x (cid:3) + + y (cid:3) + ) > holds assuring that (cid:21) ( x (cid:3) + ; y (cid:3) + ) > within Z1 so that the characteristicsfrom C ( II ) have positive slope in this region. These conditions give support to the picture outlinedin Figure 1. 12n the forthcoming examples we show that for the initial conditions (27) with h > one canchoose values for x r and y r for which both (cid:21) ( x; y ) and (cid:21) ( x; y ) are approximately conserved withinthe domain up to the order "(cid:13) : This means that the characteristic lines from each set C ( I ) and C ( II ) are nearly parallel to each other. Therefore, in looking for solutions to u ( z; t ) and c ( z; t ) ; shockformation as well as rarefaction regions can be neglected up to this order of approximation. On thebasis of the results discussed above one concludes that the system evolves in time exhibiting tworegions on the ( z; t ) plane presenting distinct values for the excess of cargoes and motors densities.Within QR , these quantities preserve their initial values whereas within Z1 , these quantities showalso a dependence on the excess of cargoes y r at the boundary z = 0 . The speed with which Z1 invades QZ is (cid:21) ( x (cid:3) + ; y (cid:3) + ) (29). In the present context this characterizes a wave carrying excess (ordepletion) of motors and cargoes that travels towards the positive direction (anterograde transport)along the axon at this speed. We shall see next how this picture may change for h < . < In this case one notices that at t = 0 , both (cid:21) ( h; and (cid:21) ( h; are positive. We shall show thatthese features are consistent with the solutions x ( z; t ) and y ( z; t ) for all times t > and space z (cid:21) domains, as suggested by Figure 2. In the absence of rarefactions, the characteristics from bothsets C ( I ) and C ( II ) emerging from z = 0 at t = 0 divide the plane z (cid:21) , t > into three regions,which will be referred to as the Quiet Region QR (upper plane z ( t ) > (cid:21) t ), intermediate region Z1 ( (cid:21) t < z ( t ) < (cid:21) t ) and lower region Z2 ( (cid:21) t > z ( t ) ).Likewise the case h > , the functions x ( z; t ) and y ( z; t ) can be resolved at any point ( z; t ) within the QR by extending back from ( z; t ) to the initial values at t = 0 along characteristics fromboth C ( I ) and C ( II ) leading to the same equations as in (30). This allows us to conclude thatboth C ( I ) and C ( II ) are straight lines and also that x Q (cid:17) x ( z Q ; t Q ) = h and y Q (cid:17) y ( z Q ; t Q ) = 0 ; at any point Q within this region.The analysis of the solutions within the intermediate region ( Z1 ), are performed by choosingany pair of points at the crossing of any two characteristics chosen from the sets C ( I ) and C ( II ) ,say points R and S; both at the same characteristic of C ( II ) ; as shown in Figure 2. Analogousequations (31) can be set for these two points from which we conclude that the C ( II ) are straightlines within Z1 . However, the (cid:133)rst of the equations formulated in (33), i.e. R (0 ; t ) = R ( z; doesnot hold in this case. In order to (cid:133)nd the solutions x (cid:3)(cid:0) and y (cid:3)(cid:0) to x and y for h < at any point ( z (cid:3) ; t (cid:3) ) of Z1 , we observe that the set (35) holds also for h < yielding x R = x T and y R = y T at13ny two points R and T of a C ( I ) . This implies that (cid:21) is also conserved within Z1 and thus thecharacteristics of C ( I ) are also straight lines within this region. The last two equations in (33) canthen be resolved, which yields in the present case, y (cid:3)(cid:0) = 14 (1 + 4 h ) (1 + 4 y r ) (cid:0) (1 + 4 x r ) (cid:18) "(cid:13) (cid:19) (1 + 4 y r ) + (1 + 4 x r ) x (cid:3)(cid:0) = h (cid:0) (cid:18) "(cid:13) (cid:19) y (cid:3)(cid:0) (38)The characteristics of C ( I ) have positive slope, i.e. (cid:21) ( x (cid:3)(cid:0) ; y (cid:3)(cid:0) ) > at all points within Z1 if x (cid:3)(cid:0) + y (cid:3)(cid:0) = h (cid:0) "(cid:13) y (cid:3)(cid:0) < . Using the result (38) for y (cid:3)(cid:0) , one shows that this condition is equivalentto: (1 + 4 y r ) > (cid:18) "(cid:13) + 4 h (cid:19)(cid:18) "(cid:13) (cid:0) h (cid:19) (1 + 4 x r ) : (39)Consider now any point A = ( z A ; t A ) within the Z2 , at the crossing of two characteristics, onefrom the set C ( I ) and the other from C ( II ) : The quantities x ( z A ; t A ) and y ( z A ; t A ) are resolvedwithin the Z2 by extending back from ( z A ; t A ) along both of these characteristics up to the boundaryat (0 ; t ) ; for any t (cid:21) : Since the excess densities x r and y r at the boundary z = 0 are kept at theconstant values as time passes these lead to: R ( z A ; t A ) = R (0 ; t ) R ( z A ; t A ) = R (0 ; t ) (40)It follows that x A = x r and y A = y r (41)at any point A = ( z A ; t A ) within Z2 . Both C ( I ) and C ( II ) are straight lines within Z2 . Moreover,for x r + y r < (42)the signs of both slopes do not change with regard to the other regions QZ and Z1 .From this analysis one concludes that for h < ; and under absence of shocks or rarefactionof characteristics, the system evolves in time exhibiting three regions on the ( z; t ) plane each ofwhich presenting, in the most general case, di⁄erent values for the excess of cargoes and motors.14ithin QZ , these quantities preserve their initial values. Within Z1 , these quantities show botha dependence on the initial values and on the values of excess of motors x r and that of cargoes y r at the boundary z = 0 , as shown in (40). Within Z2 , x and y coincide with their values at theboundaries (41). Under these conditions, the speed at which Z1 invades QZ is (cid:21) ( x (cid:3)(cid:0) ; y (cid:3)(cid:0) ) . Thespeed at which the region Z2 invades Z1 is (cid:21) ( x r ; y r ) : The above results for the behavior of both motor and cargo densities achieved analytically for certaininitial and boundary conditions are limited to situations for which shocks among characteristic linesand rarefaction regions are either absent or can be neglected within the considered perturbativelimit. The numerical simulations performed in parallel to the analysis presented above to studythe temporal behavior of u ( z; t ) and c ( z; t ) aims to extend the scope of applications attempting toapproach more closely the available data on axonal transport. In this Section we brie(cid:135)y discusssome features of the numerical code developed to this end.The non-linearity of the set of Eqs. (9) can create solutions with both discontinuities and smallscale smooth structures. These features together with the expectation that the system developsvelocities at very di⁄erent orders of magnitude, led us to choose, among many numerical methodsavailable, a (cid:133)nite-di⁄erence (FD) with a high-resolution shock-capturing (HRSC) scheme, widelyused in astrophysical and cosmological codes [21], [22], [23].The FD method has been chosen since it does not require the solutions of Riemann problemsat each interface between computational cells [16] reducing in this way the time costs if comparedto the conventional (cid:133)nite volume method (FV). The HRSC methods o⁄er high order of accuracy,sharp descriptions of discontinuities and convergence to the physically correct solution. It has theadvantage of treating discontinuous solutions consistently and automatically wherever they appearin the (cid:135)ow [24].In this work we have implemented our own numerical code using modern FORTRAN witha modular approach. The code uses up to (cid:133)fth-order reconstruction in the characteristic (cid:133)eldsand a local Lax-Friedrichs (cid:135)ux splitting [25]. For a reconstruction scheme, we implement the clas-sic Weighted Essentially Non-Oscillatory (WENO) schemes with third (WENO3) and (cid:133)fth order(WENO5) [16], [26] and two di⁄erent improved methods, a third order WENO3p [27] and a (cid:133)fthorder WENO5Z [28]. In addition to the aforementioned references, a brief and practical expla-15ation of these methods can be found in Ref. [23]. The time integration of the ODE·s obtainedfrom the discretization of Eqs. (9) has been made using a third order Strong Stability-PreservingRunge-Kutta (SSP-RK) scheme [29].In order to ensure accuracy in the implementation of the code, we (cid:133)rst wrote it for the classicalEuler equations to run some of the usual validation tests. We then performed a few minor changesto (cid:133)t the code for our purposes. Figure 3 depicts the results for the density curves obtainedfrom three di⁄erent tests and corresponding reconstruction schemes. In the upper left, we showthe Sod test [30] evolved up to the time t end = 0 : with cells. This test does not imposesigni(cid:133)cant computational di¢ culties. Analytical solutions to the equations allows one to check forthe accuracy with which the discontinuities are described by the implemented schemes. The blastwave considered in the second test (Figure 3, upper right) has been evolved up to t end = 0 : alsowith cells. This is a stronger test because the initial pro(cid:133)le presents a gradient of (cid:133)ve ordersof magnitude in pressure. The evident di⁄erences in accuracy reached by the di⁄erent schemesexpresses the superiority of (cid:133)fth order methods. Finally, the results using the iterative blast waves[31] are shown in Figure 3, lower left up to t end = 0 : , for cells. These results suggest that thebest accuracy among the considered methods is provided by the WENO5Z which solves slightlybetter for the peaks and valleys than the classic WENO5. Since this problem does not presentanalytical solutions, we compare the results for this cells with the outcomes of a more accuratecalculation, referred to as "exact" and obtained using WENO5 with cells . All the above tests suggest that the methods have been well implemented for the classicalone-dimensional Euler set composed of three equations for four unknowns (mass density, velocity,pressure·and energy density) which must be complemented with an equation of state (EoS) pro-viding in this way a relationship among the unknowns. The equations considered in the presentwork compose a closed set with two equations and two unknowns (the motor and cargo densities),making it possible to perform some simpli(cid:133)cations in respect to the application to the Euler prob-lem, in order to improve the computational costs. An additional test applied to the LinearizedGas Dynamics was then performed [32] . This corresponds to a classical Riemann Problem with a(cid:135)uid at rest and density equal to 1.0 on the left side and 0.5 at its right. We then investigate theevolution of the density and velocity pro(cid:133)les of an ideal gas with transmissive boundary conditions.In
Figure 3, lower right it is shown the results for density and velocity pro(cid:133)les obtained bothanalytically and also by means of a simulation run for 100 cells using a WENO5Z scheme. The twosolutions are in good agreement with each other, even near discontinuities.16
Discussion of the results: theory and simulations
The analytical expressions obtained in the preceeding section for motor and cargo densities arerestricted to the perturbative regime for which the parameter "=(cid:13) << : In addition, the initialand boundary values have been chosen as small di⁄erences x and y with respect to the densitiesat the de(cid:133)ned background corresponding to the value = in Eq. (22). Nonetheless, the resultsobtained for the two eigenvalues (cid:21) and (cid:21) in Eq. (15), are only restricted by the smallness of "=(cid:13): The simpli(cid:133)ed expressions (26) or (29) for the particular initial conditions expressed in terms ofthe parameter h; resulted from direct substitution of the de(cid:133)ned quantities x and y . This suggeststhat the speed of the interface separating the QZ from Z that in the absence of shocks is givenby v c = (cid:21) ( x (cid:3) ; y (cid:3) ) ’ (cid:21) ( h; , can be analyzed for any h within the interval (cid:0) = < h < = : Thisallows us to predict that, for j h j << the speed of the interface is close to its maximum value (cid:21) max2 = 0 : whereas for h (cid:24) : the speed (cid:21) ( h; of the referred interface shall attain very smallvalues. Therefore, although the expressions derived above for (cid:21) and (cid:21) in terms of the excesses x and y do not provide exact results for the speed of the interfaces, especially in cases for which shockand rarefaction of characteristics are determinants to the dynamics, these expressions function asimportant guides to direct our numerical experiments in order to scrutinize the parameter spaceand investigate the diverse behaviors of the quantities of interest. The time evolution of motorand cargo densities for small j h j are expected to display very di⁄erent pro(cid:133)les compared to thoseproduced for h (cid:24) : . Although one cannot predict based on of the above calculations what wouldbe the fates of the motor and cargo density pro(cid:133)les at these relatively high values of h; one canindeed envisage that the interface separating the QZ from its neighboring region, might travel atextremely low speeds, if compared to the corresponding speeds at small values of j h j : The following examples support these expectations. Cases (1)-(6) shown in Table 1 illustratethe behavior of the excess of cargoes and motors for j h j << at the indicated initial x = h ; y = 0 and boundary values x r = x (0 ; t ) and y r = y (0 ; t ) . It is also depicted the speed of the correspondinginterface between QZ and Z : For these, the motor hopping rates are (cid:133)xed at p = 0 : (unloaded)and q = 0 : (loaded). Although other possibilities could be validated, these choices keep theexpansion parameter at small values "(cid:13) = 0 : . The cargo hopping rates between neighbormotors are (cid:133)xed at f = g = 0 : . Fixing f = g leads to (cid:17) = (cid:18) = 1 (10) that is consistent with thedevelopment in Section . The choices for h and boundary values x r and y r in each example andthe corresponding numerical v num and analytical v c values obtained using the results of precedingSection are compiled in Table1. v c = (cid:21) ( x (cid:3) ; y (cid:3) ) characterizes the velocity of sharp interfaces between17he regions QR and Z . Table1 -
Analytical and numerical results for density wave pro(cid:133)les andaccompanying velocities within the perturbative regime.The results for x (cid:3) and y (cid:3) evaluated with the help of expressions (34) for h > or (38) for h < ; and the corresponding x (cid:3) num and y (cid:3) num obtained by numerical simulation for the excess of motorsand cargoes within region Z are also shown in Table1. In Figure 4 there are represented the wavedensity pro(cid:133)les for cases (2) (Fig. 4a) and (3) (Fig. 4b) at di⁄erent instant of times, as indicated.These are chosen from examples at typical boundary values x r and y r for which the system presentsat initial times either slight depletion of motors h < or a slight excess h > with respect to thebackground at u = 1 = . The corresponding di⁄erences in the average cargo velocities re(cid:135)ect thee⁄ects due to motor jamming. The analytical results achieved for the velocities and density pro(cid:133)lesare in fact very well represented by their corresponding numerical results. This supports the ideaof using these expressions to guide the choice of parameters in order to investigate the behavior ofthe quantities of interest for the axon transport problem within regions away from the perturbativelimits.Results in Table 2 for speed of the wave front v pq at which the QZ is invaded as time passeshave been obtained by numerical simulation, as described above. These are indexed by the valuesfor the parameters p and q as shown in each column. We use ( p = 0 : , q = 0 : ; ( p = 0 : , q = 0 : p = 0 : , q = 0 : and ( p = 0 : , q = 0 : . These results for v pq cannot be reproducedby the analytical expressions derived above because the initial and boundary values have been(cid:133)xed, respectively, at h = 0 : , x r = (cid:0) : and y r = 0 : which do not fully satisfy the perturbativeconditions. Nonetheless, we expect in these cases that the pulse (cid:133)xed by the boundary valuesadvancing as a wave pro(cid:133)le over Z in each case do that at very small speeds, as suggested by theexpression for (cid:21) in (26). The results below illustrating the dependence of v pq on the microscopicparameters f and g corroborate these predictions: for large h; i.e. h = 0 : the numerical valuesdepicted for v pq are typically orders of magnitude less than the values v pq obtained for j h j = 0 : : Figure 5 illustrates the corresponding cargo and motor wave density pro(cid:133)les at di⁄erent instantsof time obtained numerically for two situations, case (5) Fig.5(a), and case (6) Fig.5(b), takenfrom Table 2. For chosen values of parameters f and g; we compare the pro(cid:133)les of h = 0 : h = (cid:0) : at (cid:133)xed p = 0 : and q = 0 : . These results suggest that to the relatively largedi⁄erences depicted by cargo velocities there are associated large di⁄erences on the wave pro(cid:133)le, asit is observed in experimental data [2]. Table 2 -
Cargo average velocities v pq for motor hopping parameters p (unloaded), q (loaded) and excess motorson the track h at initial times as indicated, within the non-perturbative regime.Table 2 depicts typical behaviors of the average speed v pq at which cargoes from the boundaryat z = 0 advance over the QZ . We observe that:(i) v pq decreases as the excess of motors h on the axon increases. As h reaches : ; keepingthe remaining parameters at (cid:133)xed values v pq may decrease to such small values up to 3 orders ofmagnitude smaller than those reached at j h j = 0 : . This is due to an increase of motor tra¢ c jamthat hampers loaded-motor motion through process q . In this case, the important mechanism fordriving cargo movement is cargo hopping regulated by rates g (forwards) and f (backwards). Inturn, this may explain the existence of variations in the velocity of slow components that travel onaxons known as slow component-a (Sca) and slow component-b (Scb), whose speed may di⁄er fromeach other up to one order of magnitude. These variations are compatible with results in Table 2as parameters g and f change at high values of h: (ii) If f = g and for h not too large, p and q become the dominant processes to drive cargoes.The relatively high values of v pq reached at j h j = 0 : (Table 2) re(cid:135)ect the success of mechanisms p and q to drive the cargoes in the absence of heavy tra¢ c jam (congestion).(iii) In general, v pq increases with parameters p and q , as expected. However, the e⁄ects ofincreasing p (unloaded motors) and q (loaded motors) in a situation of high motor density shouldnot result in signi(cid:133)cant increasing of velocities since tra¢ c jam is the determinant e⁄ect in thiscase, as mentioned. Consistently with this, we expect that the e⁄ect of increasing v pq with p and/or q shall be noticeable at low motor densities, as observed in the results depicted above.(iv) The e⁄ects on v pq of changing cargo hopping rates from g < f to g > f are such that v pq increases at high motor excess h = 0 : but decreases at small j h j = 0 : : For discussing thise⁄ect one should notice that the dominant process for driving cargoes at high densities is hoppingcontrolled by the rates f and g . We should then expect that v pq increases by increasing the forwardhopping rate g with respect to f . On the other hand, one observes in Table 2 that at low motordensities j h j = 0 : the forward speed v pq does not increase with g . On the contrary, it increases19s the backward cargo hopping rate f increases with respect to g: This can be understood sinceat low densities, the clusters of motors eventually assembled are likely to appear isolated from theothers so that for g > f the loads would end up stacked in front of each cluster. In this case, andfor relative small q; each cluster would take a long time to reach the back end of the next cluster.On the other hand, for g < f the cargoes would be more prone hopping backward to accumulateat the back end of each cluster until another cluster, led by its front motors through process p; would join the former providing in this way extra home for process f to continue. Accordingly,the accumulation of cargoes at the back end of the new cluster of motors assembled in this way bythe union of the previous two would be the source of shock waves driving cargoes to the forwarddirection.In this context, we want to argue about a possible connection between the results achievedhere and those reported from experiments performed to investigate the e⁄ects of Kinesin-BindingProteins (KBP) on the transport of particles along microtubules [13]. In that study, one investigatesthe e⁄ects on particle mobility under variations of the amount of KBP in the considered systems.E⁄ects of increasing KBP have been correlated to a decrease of KIF1A motors attached to themicrotubules, and to an increase in the average velocities of the remaining moving KIF1A motorsafter the introduction of KBP. On the contrary, it is observed in another experimental series thataverage velocities of Kinesin based mobile Rab vesicles (cargoes) decrease at increasing KBP (i.e.,decreasing the amount of motor attached).Figure 6 depicts the predictions of the model for the behavior of average velocities, undervariations of the excess h for certain choices of parameters, as indicated. Both velocities presenta unique global maximum at h = h max , that depends on the choices for f and g . This behaviorcan be understood as resulting from two e⁄ects: increasing the excess of motors at low densities ( h < h max ) would contribute to drive the cargoes that are released from the reservoir. On the otherhand, e⁄ects due to tra¢ c jam would prevail at large motor densities ( h > h max ) , as discussedabove. We may then suggest that the data in [13] might be understood as resulting from a balancebetween these two e⁄ects under the condition that the related experiments explore di⁄erent regionsof motor densities. 20 Conclusions and additional remarks
Our model describes axon transport of di⁄erent particles in a uni(cid:133)ed mode as these are conveyeddirectly by molecular motors, in the absence of intermediates i.e. vesicles or any other a priori moving structure [4], [5]. The idea is that the cargoes are able to detach from molecular motorsavailable on the axons and attach directly to neighbor motors at certain probability rates. Thishopping mechanism is favored under situations of a heavy tra¢ c jam at which the motors assembleinto large clusters. We propose that such clusters may supply the need for the alleged intermediatecarriers. The model accounts also for processes leading to direct movement of cargoes attached tomotors if allowed by local jamming conditions. These features are able to predict quantitativelythe observed di⁄erences between the velocities of the particles conveyed in axons at slow and fastmodes, essentially as a consequence of changing the density of motors available because it is thedensity of motors that regulate the tra¢ c jam. Moreover, the model suggests that the stochasticmechanism of cargo hopping from motor to motor can also explain the di⁄erences between the twocomponents at slow rates (Sca and Scb) by tuning the corresponding hopping parameters f and g .In order to be able to extract these properties from the model and describe phenomenologicalaspects of these systems concerning cargo velocities, we have based our analysis on particularsolutions to the pair of coupled equations of the Burgers type describing the dynamics of de(cid:133)neddensities of motors and cargoes. These equations resulted from a mean-(cid:133)eld approach of theoriginal stochastic description. We use a perturbative scheme for certain initial and boundaryvalues considering the particular choice for which f = g . To our knowledge, the analytical solutionsachieved in this way to the system of quasi-linear PDE·s are new in the literature. In turn, thesesolutions pointed out to the regions of parameters, extending far beyond the perturbative limitsto be investigated numerically in order to obtain information about the system that may relate tothe existing data. Such interplay between numerical and analytical approaches was shown to becrucial in exploring the properties of interest.It is worth mentioning that the description of interacting particle system considered here allowsfor studying the e⁄ects produced by one kind of carrier motor and one kind of cargo at a time.Di⁄erent motors and cargoes should be taken into account by a separate dynamics characterizedby some speci(cid:133)c set of parameters. We have presented results for motors that have been chosenpossessing anterograde (plus-ended) bias. Our predictions suggest that depending on the initialconditions, cargoes normally recognized as slow components (either Sc a or Sc b ) may also be con-21eyed in fast-moving mode as already observed, for example, in studies exploring the dynamicsof synapsin proteins [33]. Accordingly, what would determine if a set of cargoes move in a fastor slow mode is a combination of factors comprising their ability to interact with motors and thequantity of motor available on the axon. In a forthcoming study, we shall consider this same kindof description attempting to unravel the corresponding pulse-like wave pro(cid:133)les observed in axontransport since these present distinct and typical features for slow and fast components. We wish to thank the faculty members of the Laborat(cid:243)rio de Automa(cid:231)ªo e Controle (LAC) fromEscola PolitØcnica da Universidade de Sªo Paulo (Poli-USP) for their criticisms and helpful sug-gestions. VRM acknowledges the (cid:133)nancial support provided by Universidade de Sªo Paulo (USP).
References [1] S. Millecamps, J-P Julien, Axonal transport de(cid:133)cits and neurodegenerative diseases, NatureReviews Neuroscience 14: 161-176 (2013).[2] S. Roy, Seeing the unseen: the hidden world of slow axonal transport, The Neuroscientist :20, 71-81 (2014).[3] M. Tytell, M. M. Black, J. A. Garner, R. J. Lasek, Axonal transport: each major rate com-ponent re(cid:135)ects the movement of distinct macromolecular complexes, Science 214: 179-181(1981).[4] J. A. Garner, R. J. Lasek, Cohesive axonal transport of the slow component b complex ofpolypeptiedes. J. Neurosci 2:1824-1835 (1982).225] Y. Tang, D. Scott, U. Das, D. Gitler, A Ganguly, S. Roy, Fast vesicle transport is required forslow axonal transport of synapsin, J. Neurosci: 33(39): 15362-15375 (2013).[6] See, for example, F. Berger, C. Keller, M. J. I. M(cid:252)ller, S. Klumpp, and R. Lipowsky, Co-operative transport of molecular motors, Biochem. Soc. Trans. 39: 1211-1215 (2011).[7] J. J. Blum, M. C. Reed, A model for fast axonal transport, Cell Motility 5:507-527 (1985).[8] S. Ochs, Characterization of fast orthograde transport, Neuroscience Res. Prog. Bull. :20, 19- 31 (1980).[9] A. Brown, Axonal transport of membranous and nonmembranous cargoes: a uni(cid:133)ed perspec-tive, The J. Cell Biol.:160, 817 - 821 (2003).[10] C. Goldman, E. T. Sena, The dynamics of cargo driven by molecular motors in the context ofan asymmetric simple exclusion process, Physica A 388: 3455-3464 (2009);[11] C. Goldman C., A Hopping Mechanism for Cargo Transport by Molecular Motors on CrowdedMicrotubules, J. Stat. Phys. 140: 1167-1181 (2010).[12] L. W. Rossi, P. K. Radke, C. Goldman, Long-range cargo transport on croweded microtubules:the motor jamming mechanism, Physica A:401, 319-329 (2014).[13] J. T. Kevenaar, S. Bianchi, M. van Spronsen, A. Ahkmanova, M. O. Steinmetz, C. C. Hoogen-raad , Kinesin-Binding Protein Controls Microtubule Dynamics and Cargo Tra¢ cking by Reg-ulating Kinesin Motor Activity, Current Biology, 26: 849 - 861 (2016).[14] A. Kunwar, S.K. Tripathy, J. Xu, M.K. Mattson, P. Anand, R Sigua, M. Vershinin, R. J.McKenney, C.C. Yu, A. Mogilner, S.P. Gross, Mechanical Stochastic Tug-of-war models cannotexplain bidirectional lipid-droplet transport, Proc. Natl. Acad. Sci. USA 18960 - 18965 (2011).[15] W. M. Saxton, P.H. Hollenbeck, The axonal transport of mitochondria, J Cell Science, 125:2095 - 2104 (2012).[16] C.-W. Shu, High-order Finite Di⁄erence and Finite Volume WENO Schemes and DiscontinuousGalerkin Methods for CFD, Int. J. Comput. Fluid D. 17, 107 (2003).[17] L. Conway, D. Wood, E. T(cid:252)zel, J. L. Ross, Motor transport of self-assembled cargoes incrowded environments, Proc. Natl. Acad. Sci. 109: 20814-19 (2012).2318] A. Pelizzola, M. Pretti, Cluster approximations for the TASEP: stationary state and dynamicaltransition, Eur. Phys. J. B. 90: 183-190 (2017).[19] L. C. Evans,
Partial Di⁄ erential Equations 2nd Ed., (cid:15)
Figure 1 - Characteristic lines for the system (9) with initial and boundary conditions asspeci(cid:133)ed in (27) and (28), h > : (cid:15) Figure 2 - Characteristic lines for the system (9) with initial and boundary conditions asspeci(cid:133)ed in (27) and (28), h < : (cid:15) Figure 3 - Density pro(cid:133)les for numerical tests validation.
Upper left : Sod(cid:146)s problem.
UpperRight : Blast Wave Problem.
Lower left : Interactive Blast Waves.
Lower right : LinearizedGas. (cid:15)
Figure 4 - Time evolution of motor (gray) and cargo (black) density pro(cid:133)les for illustrativeexamples taken from Table 1. (a) - case (2) ; (b) - case (3), for the initial cargo occupationat c = 0 : and reservoir excesses (28) at x r = (cid:0) : and y r = 0 : : For the scale factorwe use (cid:24) = 1 what sets (cid:11) = p and (cid:12) = q (10). (cid:15) Figure 5 - Time evolution of motor (gray) and cargo (black) density pro(cid:133)les for illustrativeexamples taken from Table 2, for p = 0 : and q = 0 : . (a) h = 0 : case (5) ; (b) h = (cid:0) : case (6), for the initial and reservoir occupations, and scale factor set as in Figure 4. (cid:15) Figure 6 - Variation of the average (a) cargo and (b) motor velocities with the excess motordensity h; at c = 0 : p = 0 : and q = 0 : ; for f = g = 0 : (cross); f = 0 : ; g = 0 : (plus); f = 0 : ; g = 0 : (dot). 25 igure igure 2 igure 3 igure 4(a) igure 4(b) igure 5(a) igure 5(b) igure 6(a) igure 6(b) able 1able 1