Mechanical Behavior of Axonal Microtubules; the Effect of Fluid on the Rupture of Axonal Microtubules

Axonal microtubules are dynamically instable bundles in the interior part of the axon. The dynamics of these bundles are of vital importance in the behavior of axon such as their degeneration. Each axon typically contains 10~100 microtubule bundles with average length of 4μm. These bundles are coated with cytoplasm and are cross linked with random number of tau proteins. In some circumstances such as acceleration or deceleration of head in space or during the strike, they are placed in tension which may cause rupture of these bundles or disconnection of tau protein cross links. Mechanical behavior and rupture modality of microtubule bundles are becoming more and more important recently. In our model, viscoelastic microtubule bundles constituted from several discrete masses connected to the neighboring mass with a standard linear solid (SLS), a spring damper model. In addition we take into account the effect of cytoplasm by Dissipative Particle Dynamic (DPD) to investigate the rupture nature and mechanical behavior of these bundles and the effect of cytoplasm on their mechanical behavior. We obtain these results for various amounts of suddenly applied end forces to the group of axonal microtubule bundles.

A b s t r a c t -Axonal microtubules are dynamically instable bundles in the interior part of the axon. The dynamics of these bundles are of vital importance in the behavior of axon such as their degeneration. Each axon typically contains 10~100 microtubule bundles with average length of 4µm. These bundles are coated with cytoplasm and are cross linked with random number of tau proteins. In some circumstances such as acceleration or deceleration of head in space or during the strike, they are placed in tension which may cause rupture of these bundles or disconnection of tau protein cross links. Mechanical behavior and rupture modality of microtubule bundles are becoming more and more important recently. In our model, viscoelastic microtubule bundles constituted from several discrete masses connected to the neighboring mass with a standard linear solid (SLS), a spring damper model. In addition we take into account the effect of cytoplasm by Dissipative Particle Dynamic (DPD) to investigate the rupture nature and mechanical behavior of these bundles and the effect of cytoplasm on their mechanical behavior. We obtain these results for various amounts of suddenly applied end forces to the group of axonal microtubule bundles.

INTRODUCTION
Axonal microtubules are dynamically instable and microtubules cycle through the length of axons.
Dynamic instability of microtubules causes the microtubules to occupy greater volume in the interior of axons (Buxton et al. 2010). Microtubules' dynamics are impressive in the function of axons and their degeneration (Buxton et al. 2010). Microtubular geometry is formed by microtubule-associated (Tau) proteins and tubulin polymerization (Shamloo et al. 2016).
Microtubular function can be affected by acetylation, glutamylation, glycation, phosphorylation and tyrosination (Perez et al. 1999, Shamloo and and dysfunctions of axonal microtubules can result in Alzheimer's disease (AD) neurons (Majid T. et al. 2015). In AD neuron, microtubules reduce in number and length, significantly (Buxton et al. 2010).
Typically each axon contains ~10-100 microtubule/cross sectional view (Fadic et al 1985) and the average length of each microtubule bundle is ~4µm in the axon (Yu et al. 1994) and they have polar orientations. Axonal microtubules are cross linked to their neighboring microtubule bundles with tau proteins and these tau proteins are considered as one of the important elements for axonal strength in response to the mechanical force (Prechsel et al. 1992).
In some circumstances, such as in acceleration or deceleration of head in the space, microtubule bundles are places under tension force (Tang-Schomer et al. 2010). The tension applied to the nervous system cause Traumatic Axonal Injuries. Axonal injury is one of the most common nervous system injuries (Smith et al. 2003, Sanjith 2011, Kirkcaldie et al. 2016) and patients who are suffering from axonal injuries are commonly experience numerous functional deficits (Adam et al. 1989, Sanjith 2011, Blumberg et al. 1989, Williams et al. 1990).
Microtubule bundles can resist mechanical stress (Peter and Mofrad 2012). Since the injuries of microtubules are shown to be one of the main consequences of the traumatic brain injury these bundles have been studied widely from mechanical point of view (Margulies et al. 1990.
Most studies also show that the rupture of microtubules causes severe damages in axonal brain injuries (Tang-Schomer et al. 2010, Mehrafrooz and. There are numerous models which study the mechanical injuries of microtubules computationally and experimentally (Kilinc et al. 2008and 2009, Meaney et al. 1994, Dixton et all 1988and 1991, Feeney et al. 1991, Mishra et al. 2005. Mathematical and computational modeling has been proved to be the easiest way to study the mechanics of nervous system (Erickson andO'Brein 1992, Flyvbjerg et al, Shamloo et al. 2015 B and. In axonal computational modeling, microtubule bundles are assumed to be composed of discerete masses connected with mechanical elements (Buxton et al. 2010).
Microtubule bundles have been recently modeled both as discrete masses and continuum elements (Shahinnejad et al. 2013, Teixeira et al. 2014, Ansari-Benam et al. 2015, Adnan et al. 2018). Shahinnejad et al. 2013 used finite volume to model axonal microtubules under tension and they also assumed microtubules as discrete masses. In 2012, Peter and Mofrad suggested a model for axonal microtubules. They endeavored to represent an intact model for axonal microtubules. In their model microtubules were assumed to be lots of discrete masses connected to their neighboring microtubule bundle with number of tau proteins (  al. he ior re al of microtubule comparing to their behavior in the absence of cytoplasm. In our model, the effect of fluid around the microtubules on the behavior of microtubules under the action of varying stresses is studied in a three dimensional dimensional (2D) platform. We determined the steady state behavior of axon for different force rates and different force magnitudes under the effect of cytoplasm and compared this behavior with one which is not under the effect of cytoplasm. In Fig. 1(a), a neuronal cell is shown with its connected axon, dendrites and microtubules. As it is shown in this figure, the axonal microtubules are positioned within the interior part of axon ( Fig. 1(b)).  Fig. 1(d). The space between each two microtubule bundle is 0.045 µm which is 4.5 times the distance between each two neighboring mass on an individual bundle which is filled with cytoplasm.

b. Fluid Modeling
Fluid around the microtubule bundles is modeled using Dissipative Particle Dynamic (DPD) method (Liu et al. 2014). Fluid is assumed to be a large number of DPD particles which filled the whole surrounding of axonal microtubules. The allocation of DPD particles is face-centered cubic lattice which is dense enough to model fluid particles. To model the random position of these particles, we changed the initial position of each particle randomly by a very small ratio (maximum of 0.1 of distance between each 2 particles). In addition, Maxwell-Boltzmann distribution has been used in this study to specify the initial velocity of each DPD particle (Walstad 2013). Fig. 3 illustrates 3 dimensional axonal microtubules surrounded with fluid particles.

c. Forces
Dissipative Particle Dynamic (DPD) follows the basic principles of Navier-Stokes equations with certain alteration (Liu et al. 2006). The computational cost of this method is much less than other similar methods such as Molecular Dynamic (MD) (Keaveny et al. 2005). DPD forces are considered all possible interparticle forces between fluid-fluid particles and fluid-Microtubule bundle particles There are three DPD forces exerted on particle 'i' by surrounding particle 'j' (Eq.
Conservative Repulsive Force is given by Eq.2. In this equation 'a ij ' is the maximum repulsion force between two particles 'i' and 'j', 'r ij ' is the distance between particle 'i' and particle 'j',

‫ݎ‬
is the unit vector that determines the direction from particle 'j' to 'i'. The Conservative repulsive force depends on the cut-off radius (r c ) which is the maximum range of F C .

‫ܨ‬
In addition, the dissipative and random forces are given by the equations 3 and 4, respectively, where 'γ' and 'σ' are characteristic strength of these forces and their relation are given by equation 5.
In these equations 'V ij ' is the relative velocity of particle 'i' related to particle 'j', θ ij is a random value with zero mean and unit variance, It can be considered that each particle can be affected by 7 effective forces: all possible external forces (F ext ij ); force between 2 neighboring mass in a same microtubule bundle (F SLS ij ); force that is acted upon tau proteins by a neighboring mass from the neighboring microtubule bundle (F Tau ij ); the force which is due to the effect of bending moment of a sole microtubule bundle (F b ij ); force between a DPD particle and a discrete mass of microtubule bundle which is known as particle-microtubule bundle force (F P-M ij ); and force between two DPD particles known as particle-particle force (F P-P ij ) and as it is noted earlier 'F P-M ij ' and 'F P-P ij ' are themselves divided into 3 separately part: Conservative repulsive force, dissipative force, and random force. SLS forces also contain 2 parts: elastic force which is because of the springs in a SLS unit, 'F el ij ', and dissipative force of damper in the SLS system, 'F P-P ij ' ‫ܨ(‬ ௌ ௌ ൌ ‫ܨ‬ ‫ܨ‬ ௗ ).
We have written a MATLAB code to solve our system. In our code each particle and microtubular discrete mass has an index which make particle known during various iterations. In each time step and for each particle we measure the total force from other particles as in equation 7. When the force of a particle is measured, we are able to determine the acceleration of the particle; as a result the velocity and location of the particle in the next time step can be define using numerical integrating and derivation methods. In equation 7, 'a', 'b', 'c' and 'd' are indicators of whether particle 'i' and particle 'j' are discrete masses in the same microtubule; discrete masses in 2 various neighboring microtubule; a fluid particle and a discrete mass of a microtubule bundle; or 2 fluid particle. These indicators are logical numbers and assume only binary values of unity or zero. To illustrate it more, if the particle 'i' is a discrete mass of a microtubule bundle and particle 'j' is a DPD particle, there would be exist only ‫ܨ'‬ , ି ெ ' and 'd' would be a unit value and 'a', 'b' and 'c' would be zero numbers. We let the system goes until it reach its steady state. In various time steps, if the strain of each two neighboring mass of microtubule bundles which are connected together with SLS unit become more than 50%, the SLS unit is fractured and when the fractures increase, the system breaks down and the rupture occurs.

RESULT
We have modeled axonal microtubule bundles computationally and we are interested in finding the mechanical behavior of them under the effect of suddenly applied end forces in the presence of the fluid around them and compare our result with their behavior in the absence of fluid to find the effect of fluid on their mechanical behavior.
In the beginning, we make all of our equations non-dimensional. For finding the appropriate time step we used 10 discrete masses connected to each other with SLS units (Fig. 2). Using the system we found the appropriate time step will be equal to ~ 2.245×10 (-2) ps. Table 1 shows the material parameters used in this study (Peter and Mofrad 2012). If the microtubule bundles shown in Fig. 1(b) are assumed to be attached to a constant wall in one side (Cell body side) and act upon on other side (Axonal terminal side) by an external force, microtubule bundles will show time-varying deformation. When bundles start deforming, the overall strain of the system in each time step can be defined by Eq. (8). The first and the second summation in this equation are the over the rightmost and leftmost points. 'n' and 'L * ' are the number of microtubules and the initial non-dimensional length of the system.
According to Janmey et al (1991), we assumed that when the strain between two neighboring point mass of microtubule bundles become more than 50%, the SLS unit between them breaks down. During the real traumatic brain injury, the mechanical force on the neurons can be considered to grow rapidly in the time interval of [0, τ ] ('τ' is in the nanoseconds scale.) from zero to its final value (F) and then become a constant value in the further time steps (Fig.4). We investigated the response   It can be inferred from the Fig. 5 that when the effect of discontinuities is modeled, the steady strain of the system is increased. We could point out that the discontinuities cause growth of the system's strain.
When the surrounding fluid is modeled, the steady state strain of the system is decreased significantly and this shows that fluid decrease the stored energy of the system. According to this figure and our modeling we conclude that the fluid around the microtubules cause decreasing of final strain by 12% and this causes postponing the rupture of axonal system. which is coated with fluid particles, when the force is applied to the system, the deformation starts until the system reaches its steady state. As we increase the applied force, the steady state strain is increased and when the final non-dimensional applied force reaches its critical value of 0.275 (~ 13 µN) the rupture occurs and the connection between neuron cell and axonal terminals breaks down. We studied the nature of rupture in this study, including how rupture starts, expands and failure occurs.  We studied the effect of different parameters on the critical area of axonal microtubule bundles, where bundles are more likely to fail. Based on our observation over various forces and different force-growing times (τ in Fig. 4), critical area is always nearer to the axonal terminal compared to its distance from cell body (Fig.8). When we increase τ (force grows slowly) while force is over-critical, the critical area moves to the middle of the axonal microtubules and when τ is decreased, critical area becomes nearer to the axonal terminals.  Fig. 4) decreases, critical area shows tendency to become nearer to axonal terminals

CONCLUSION
Knowing the exact behavior of axonal microtubules is very helpful in investigating the damage on central nervous system. We modeled axonal microtubules computationally using Spring-Damper methods. We assumed microtubule bundles as a series of discrete masses connected to their neighbors with SLS units while each bundle has a point of discontinuity. Surrounded cytoplasm of bundles is modeled with Dissipative Particle Dynamic (DPD) method and interaction of microtubule bundles with fluid particles are considered to simulate the mechanical response of axonal microtubules to a specific external forces.
These forces were acted to the end of each microtubule bundle. Now, we are able to predict the nature of axonal rupture and failure of axonal microtubules and we identified the critical damage areas of microtubules. We modeled microtubules which are relatively close to the real behavior of axonal microtubules, comparing with other models. We demonstrate that the failure of axonal microtubules occur when the applied force become about 13µN.
The effect of cytoplasm on axonal microtubules is studied in our model. It is indicated that fluid around the microtubules dissipate the energy of system and cause the strain to decrease by 12% and as a result decrease the feasibility of mechanical failure of axonal microtubules which are responsible to transfer the information from cell body to axonal terminals.
Computational modeling of axonal microtubules is very important in knowing the behavior of microtubules due to the fact that it is not only time-efficient, but also easy and cost-effecting. Mechanical damage of neurons is one of the most common central nervous system injuries and microtubular failure is one of the main reasons of traumatic brain injury. Using the exact information of computational modeling of axonal microtubules assists us to develop new methods to prevent the damage on nervous system and culture microtubule bundles and replace the damaged part of real microtubules with cultured ones. In this model we improved other models by taking into account the effect of cytoplasm on the behavior and mechanical response of axonal microtubules.