Effect of Strain Rate on Single Tau, Dimerized Tau and Tau-Microtubule Interface: A Molecular Dynamics Simulation Study

Microtubule-associated protein (MAP) tau is a cross-linking molecule that provides structural stability to axonal microtubules (MT). It is considered a potential biomarker for Alzheimer’s disease (AD), dementia, and other neurological disorders. It is also a signature protein for Traumatic Brain Injury (TBI) assessment. In the case of TBI, extreme dynamic mechanical energies can be felt by the axonal cytoskeletal members. As such, fundamental understandings of the responses of single tau protein, polymerized tau protein, and tau-microtubule interfaces under high-rate mechanical forces are important. This study attempts to determine the high-strain rate mechanical behavior of single tau, dimerized tau, and tau-MT interface using molecular dynamics (MD) simulation. The results show that a single tau protein is a highly stretchable soft polymer. During deformation, first, it significantly unfolds against van der Waals and electrostatic bonds. Then it stretches against strong covalent bonds. We found that tau acts as a viscoelastic material, and its stiffness increases with the strain rate. The unfolding stiffness can be ~50–500 MPa, while pure stretching stiffness can be >2 GPa. The dimerized tau model exhibits similar behavior under similar strain rates, and tau sliding from another tau is not observed until it is stretched to >7 times of original length, depending on the strain rate. The tau-MT interface simulations show that very high strain and strain rates are required to separate tau from MT suggesting Tau-MT bonding is stronger than MT subunit bonding between themselves. The dimerized tau-MT interface simulations suggest that tau-tau bonding is stronger than tau-MT bonding. In summary, this study focuses on the structural response of individual cytoskeletal components, namely microtubule (MT) and tau protein. Furthermore, we consider not only the individual response of a component, but also their interaction with each other (such as tau with tau or tau with MT). This study will eventually pave the way to build a bottom-up multiscale brain model and analyze TBI more comprehensively.


Introduction
Blast-induced TBI has been considered highly significant on the battlefield and in sports [1,2]. Numerous studies have been performed to address this injury mechanism and consequences from different perspectives, such as post-traumatic stress disorder, aggregation of neurological disorders, damage threshold of cytoskeletal components, etc. When the head is impacted by blast-like mechanical force, the resulting pressure wave transmits to the interior of the brain and causes macroscale, tissue-level, and cellular level damage [3][4][5][6][7][8][9][10][11]. A major part of the brain's interior is built on over 100 billion neuron cells and its surrounding extra-cellular matrices. Structurally, a single neuron is composed of a soft cell body (soma), fibrous dendrite branches, and an axon fiber. Both axons and dendrites are supported by the cytoskeleton, which is primarily composed of neurofilaments logical conditions around 7 of them are actually phosphorylated [24]. The abnormally phosphorylated tau proteins are 3-4 times more phosphorylated than the normal ones, and~2 sites per mole of tau protein can be phosphorylated in normal conditions [25]. Recently, cryo-EM technology has been able to obtain 3-4 Å resolution images of paired helical filaments (PHF) in AD-affected brains [26], which is directly related to abnormal phosphorylation in tau, a clinically and pharmacologically relevant variation [27]. From the structural biology perspective, hyperphosphorylation has been perceived as one of the critical parameters which induce self-assembly of tau into paired helical filaments (PHF) or straight filaments [28]. In this regard, recent experimental and computational studies have provided critical insights into the structure and interaction behavior of tau depending on its phosphorylation state [29,30].
Not only phosphorylation and dephosphorylation are relevant, but also the strain rate dependent structural response of tau regarding this study. In some simulation studies of microtubule bundles cross-linked with tau protein, the response under different levels of uniaxial tension for different distributions of cross-links and discontinuities have been observed, and the estimated Young's modulus used for calculation was 5 MPa [31], although they have admitted that it is decidedly approximate, and tau-tau, tau-MT, etc. interactions are important candidates to determine the mechanical behavior. Their model showed that bundle failure occurred due to the failure of cross-links. Tau proteins, through complementary dimerization with other tau proteins, form bridges to nearby microtubules to form bundles [32]. The tau proteins have an edge-to-edge distance of 20 nm and hexagonal packing [33]. Furthermore, computational models have been built to study the behavior of cytoskeletal filaments [34,35] and cross-linked networks have also been investigated [36][37][38]. Earlier computational and theoretical models have shown that shear resistance provided by the cross-linked network greatly increases MT bending stiffness [39], and found bundle stiffness regime [40]. Discrete bead-spring models have been widely used to build filament and network structures [41,42]. In another finite element analysis, Young's Modulus of tau protein was assumed to be 5 MPa [43], but elevated to 62.5 MPa to reduce the flexural behavior so that the spring constant can be the same as Peter et al. [31]. The tensile test results reasonably agree with that of Peter et al., and it also analyzed the behavior of bundle under torsion. Prior continuum and computational axon modeling give us useful insights on tau protein, such as viscoelastic shear lag model [44], the failure mechanism of axon study [45], and so on. The computational continuum models and rate-dependent tests can determine the properties of the cytoskeletal components such as microtubules, and therefore are relevant to the current study.
The structure and phosphorylation studies cannot depict the complete structural and functional scenario of tau, because another critical aspect of tau protein is the MT binding region, which can work as a backbone of the structure and from which the projection domain can spread around [46]. Some recent studies have further clarified the position and structure of MT binding sites in tau, such as Rosenberg et al. [19]. Each of the four binding regions of tau is around 18 amino acids (AA) long, separated by inter-repeats of 13-14 AA length. The positively charged repeat or inter-repeat region is believed to facilitate electrostatic interactions between tau and negatively charged surface of MT [16,18,37]. In the N-terminal of the MT binding region, there is a positively charged and proline-rich region, which has phosphorylation sites responsible for the regulation of MT association to taus [47,48]. The other portions of tau protein have significant functionalities as well, such as the N-terminal region, which can contain zero to two negatively charged inserts (each 29 AA long). The N terminus and proline-rich region make up the projection domain of tau protein, which extends outward from the MT surface and determines the inter-MT distance in MT bundles [33]. Other studies show that MT binding with tau is very fast, does not depend on the MT location in the axonal shaft, and varies with MT curvature [49]. A study of intrinsically disordered tau protein folding on MT shows that tau locally folds into a stable structure upon binding [47].
Aside from studying standalone tau protein response, our interest is to obtain insight into tau-MT interaction. Tau protein organization on MT can be described as a coating or surface decoration, as suggested by some works [50]. The mode of tau protein binding with MT has been studied, and it suggests that tau protein binds along as well as across protofilaments [51,52]. Now, while considering a feasible approach to computationally determine single tau, dimerized tau, and tau-MT behavior, we must consider the earlier studies that have been performed on intrinsically disordered proteins (IDPs). Therefore, by focusing on computational studies on IDPs (especially on tau proteins), we have found that recent years have observed significant improvement in computational studies by using predicted structures. Notable examples are electrostatic study [53], structure-based molecular docking and dynamics to find inhibitors and drug targets of tau [54,55], and general aggregation behavior study on tau [56]. However, studies focusing on the mechanical response of tau are still non-existent. A recently published comprehensive review study has pointed out this current ongoing limitation [57], and an MD simulation study has shown the effect of phosphorylation on tau [58]. Therefore, this study is an attempt to address one of the significant ongoing limitations in the current literature by applying strictly mechanical loading.
Based on the existing literature as discussed above, there are certain parameters yet undetermined, such as mechanical behavior single tau filament due to strict mechanical loading (deformation at high strain rate), required stretching for separation of dimerized tau proteins, tau-MT binding strength, and whether MT subunit binding between themselves is stronger than the tau-MT bond, etc. This study attempts to address these specific questions using an atomistic computational method called molecular dynamics (MD). The rationale behind choosing this tool are twofold: (1) reliable experimentation at this length scale is rare and (2) MD is an appropriate tool to capture the molecular details of the elements involved. Specifically, for the single tau, we have obtained the secondary structure from i-TASSER predictor software [20]. For dimerized tau, we have made the model with overlapped projection domain as depicted by Rosenberg et al. [19]. For tau-MT interaction, we have used Chau et al. proposals of interaction, where tau can interact with one or both subunits of MT [16]. The effects of high strain rate on the deformation mechanism of single tau, dimerized tau and tau-MT interface are studied.

Method
The tau protein structure is obtained from the i-TASSER predictor software [20] (isoform F, which is also known as Tau-4, 2N4R is used [59]. For the amino acid sequence and definition of the regions such as projection domain, MT binding domain, and tail domain, please refer to Supplementary Material Part 4). We have used the model with the C score of 1.06 (the C score is determined based on the significance of threading template alignments and the convergence parameters of the structure assembly simulations) which we have assumed is satisfactory for an IDP. For single tau protein, we have solvated the obtained structure with TIP3P water using CHARMM-GUI [60] solvator and quick-MD simulator modules, which have been substantiated in recent MD simulation works on axonal cytoskeletal components [61,62]. Required numbers of 0.15 M KCl ions were added to obtain charge neutralization for explicit solvent simulations. The amount of 0.15 M KCl is determined at the solvation step of the "Quick MD Simulator/Solution Builder" module, where the total charge due to the presence of the protein in the system is calculated, and then 0.15 M KCl is added accordingly to neutralize the charge. The 0.15 M KCl is part of the protocol to achieve the charge neutrality of the system. The system achieves this by adding the necessary number of cations and anions as a part of the solvation. Although adding similar salt(s) such as NaCl is also acceptable, currently the module only supports adding 0.15 M KCl, which is a general protocol for charge neutralization in MD simulation of biomolecules.
The dimerized model was created by using UCSF Chimera [63], in which we have overlapped the projection domains of two identical tau proteins. We have used an implicit sol-vent technique of CHARMM [64] in LAMMPS [65] (pair_style lj/charmm/coul/charmm/ implicit command) for dimerized tau, which facilitates a bigger box size, faster calculation, and convenient observation of two tau system. The implicit solvent computes with a modification of adding an extra r −1 term for Coulombic energy calculation and skipping long-range Coulombic energy calculations.
The MT structure is obtained from the existing model built by Wells et al. [14], each subunit of which has one GTP or a GDP, along with one Mg 2+ ion in the junction. By repeating the helically arranged 13 subunits periodically along the length direction, we have created a virtually infinite MT (similar methodology has been followed in the work of Wu et al. for MT study [8]). The tau-MT interaction system was created by placing the tau binding sites in close proximity to the MT binding sites as proposed by Chau et al. [16]. We have used the implicit solvent technique for tau-MT as described above for dimerized tau.
With periodic boundary conditions in all three directions, we have minimized and equilibrated the structures for all cases (single tau, dimerized tau, and tau-MT) for 100 ps to minimize the potential energy at a targeted temperature of 310 K. Single tau, dimerized tau and tau-MT contained 6424, 12,848 and 33,292 atoms respectively. The LJ potentials are used with an inner and outer cutoff of 10 Å and 12 Å, respectively. Long-range coulombic interactions are computed by pppm style, facilitating a particle-particle particle-mesh solver with a 3D mesh.
We have used CHARMM36 [64,66,67] potential parameters with appropriate CMAP corrections [64] for all the simulations. Potentials for GDP and GTP are taken from that of ADP and ATP respectively [14]. The equilibration was performed in NVT canonical ensemble, with the temperature damping parameter of 100 fs.
In single tau protein tensile tests, we have fixed the MT binding domain and pulled the projection domain at different strain rates (10 8 -2 × 10 9 s −1 ) along the x-axis. For dimerized tau, we have fixed the MT binding domain of one protein and pulled away from another binding domain along the x-axis (10 9 s −1 and 2 × 10 9 s −1 ). For tau-MT, While keeping the upper and lower few layers of atoms fixed of the MT, the tau protein was pulled away by few atoms of the projection domain (strain rate: 10 9 s −1 and 2 × 10 9 s −1 ). For the relevance of the simulations with explicit water molecules as solvent vs. implicit water molecules, readers may refer to the Supplementary Materials of Wu et al. [8]. Although biomolecule simulations are more realistic with explicit solvation, we have adopted implicit techniques for dimerized tau and tau-MT due to the high box size required for the high strain, and for convenient observation of unfolding, stretching, and separation.
The stress-strain plots are obtained by the per-atom stress calculation and summation in LAMMPS. However, as the output is in (pressure × volume) unit, we must divide the obtained stress value by the volume of the protein (or a certain portion of the protein). The general formulation used by stress per atom command is P = (P xx + P yy + P zz )/(3 × V), where P xx , P yy , and P zz are the summations of stress/atom values for all atoms in the x, y, and z-direction respectively, and V is the summation of the volume of the atoms of the protein being considered. The approximated volume was obtained by Voronoi cell approximation, adapted from LAMMPS voro++ package [68]. The strain is simply obtained by the displacement of the atoms from the initial position. All the tensile tests are performed in the NVT ensemble, with 100 fs temperature damping parameters. The visualizations of the tensile tests are carried out by OVITO software [69].
All of the simulations were carried out by the STAMPEDE2 supercomputer of Texas Advanced Computing Center (TACC).

Single Tau Deformation
For the single tau model, we have performed the tensile tests at four different strain rates: 10 8 , 5 × 10 8 , 10 9 , and 2 × 10 9 s −1 . The MT binding region atoms are fixed, and the first few atoms of the projection domain are pulled at the −x direction. The calculated stress-strain graphs are shown in Figure 1.  (a) Simplified schematic for tau protein unfolding and stretching (εxx = tensile strain). Initially, there are multiple folding (which are distinct from each other). When loading is applied, the van der Waals force and electrostatic force between the folded portions are broken, and therefore the structure unfolds. At extreme strain, the structure becomes free of all the folds, and stretches to It can be observed that the stress-strain response of a single tau can be divided into two regions-the "unfolding" and "stretching" regions. When the mechanical load is applied, first the folds disappear one after another because of the breaking of the electrostatic and van der Waals forces between the folded portions of the filament (this phenomenon will be called "unfolding" onwards in this manuscript). However, the mentioned forces are relatively weaker than chemical bonds, as they are merely the result of proximity and positioning (3D confirmation of the structural portions in space). Therefore, during the unfolding region (see Figure 1) no significant rise of stiffness is observed. However, after the filament is free of the folded regions, the structurally linear filament will be deformed against the relatively stronger covalent bonds (this phenomenon will be called "stretching" onwards in this manuscript). Therefore, in this region there will be higher stress developed and a significant rise in stiffness will be observed, as shown in the "stretching region" in Figure 1. The unfolding and stretching phenomena are clarified by a simplified schematic in Figure 2 as well.
In Figure 1, multiple distinct slopes can be observed that reflect the stiffnesses of tau protein at different strain states. Up to~200% tensile strain, the slope is very small implying a relatively low unfolding stiffness of tau. After the unfolding is complete, we can observe pure stretching of the covalent bonds. Unfolding stiffness is obtained from the strain of 0-150%, while stretching stiffness is obtained from the strain of 200-300% where another distinct change of slope is observed. Figure 2 shows the unfolding and stretching snapshots. The unfolding stiffness and pure stretching stiffness values are shown in Table 1. According to Table 1 and Figure 1, it can be inferred that a single tau exhibits viscoelastic behavior to some extent and that its stiffness is strongly dependent on strain rate. With the increment of strain rate, tau acts like a stiffer material in both the unfolding and stretching zone, which is expected for viscoelastic soft biomaterial. The slope is in the Mega Pascal range: 0.5 GPa or~500 MPa. From the existent literature, we have not found any validated estimate of the stiffness of the projection domain of tau, but we have detail information regarding the stiffness of other cytoskeletal components, such as MT and microfilaments. Deformation studies and viscoelastic shear lag models of MT show that the MT stiffness can be in the Giga pascal range [31,44]. Studies on microfilament show that their stiffness can vary from few megapascals to >2 GPa [70][71][72][73][74]. Therefore, we can assume that our results obtained for the tau protein stiffness, both in the unfolding and stretching region, are reasonable (dependent on strain rate). From various studies on the axonal cytoskeletal components such as MT, actin, etc. we find that their mechanical behavior is strain rate dependent, and their stiffness is high (from megapascal to gigapascal range) [57]. As overall axon behavior is also dependent upon MT and tau-MT interaction, we expect that the crosslink is also sufficiently stiff to support the structure. Figure 1. Stress vs. strain plot of single tau projection domain. Up to ~200% strain, the protein keeps unfolding, and after that a sharp rise in the slope is observed, suggesting the pure stretching of covalent bonds. Unfolding and stretching zones are shown in rectangular boxes.

Figure 2. (a)
Simplified schematic for tau protein unfolding and stretching (εxx = tensile strain). Initially, there are multiple folding (which are distinct from each other). When loading is applied, the van der Waals force and electrostatic force between the folded portions are broken, and therefore the structure unfolds. At extreme strain, the structure becomes free of all the folds, and stretches to Initially, there are multiple folding (which are distinct from each other). When loading is applied, the van der Waals force and electrostatic force between the folded portions are broken, and therefore the structure unfolds. At extreme strain, the structure becomes free of all the folds, and stretches to a relatively linear filament. (b) Initial single tau structure (strain = 0%), (c) tau protein being unfolded due to pulling at 10 9 s −1 (strain = 36%), (d) tau protein being stretched (strain = 206%, only the projection domain is shown for the convenience of visualization). Color legends: green: projection domain, red: MT binding region, blue: N terminus or tail, white, and the rest: inter-repeats between the MT binding regions. The enlarged snapshots are for the first~1100 atoms for convenient visualization. Water molecules are not shown.

Dimerized Tau Deformation
The dimerized tau model was similar to the single tau, except that we have fixed the MT binding site of one tau and pulled away from the binding site of another tau to observe the developed stress and possible sliding out of the tau projection domain at extreme strain. We have plotted the stress-strain curves (calculation procedure was similar to single tau) for the projection domain of the tau protein that has been pulled at the x-direction, at the strain rate of 10 9 s −1 and 2 × 10 9 s −1 . The implicit-solvent dimerized tau model shows several stages of tension during the test. We are referring to the protein with fixed MT binding region as protein 1, and the protein being pulled as protein 2. The stages observed are (for the strain rate of 2 × 10 9 s −1 ): (i) Unfolding of protein 2 (up to 163% strain), (ii) Stretching of protein 2 (up to 257% strain), (iii) Unfolding of protein 1 (up to 334% strain), (iv) Stretching of protein 1 (up to 395% strain), (v) Disentanglement of the overlapped projection domains of the tau proteins (up to 721% strain), (vi) Sliding out or projection domain along with stretching (very fast, occurs at around 722-758% strain region), and (vii) Separation of proteins (~758% strain). For a strain rate of 10 9 s −1 , the separation occurs just above the stretch of~800%. Therefore, the observations imply that the separation of dimer also depends on the strain rate, and at a certain high strain rate, the separation occurs significantly early, although the stretch following the unfolding seems the same for both strain rates. It is to be noted that the separation stretch that is mentioned here are the percentage of strain at which the two dimers are visibly separated from each other, although the pulling away of one tau from another starts earlier (around 30% of strain earlier), and before complete separation, there are sub-stages of untangling of the overlapped projection domains of the two tau proteins. Figure 3a shows the stress-strain curves of the projection domain of protein 2 and Figure 4 shows snapshots of the stages observed. The single tau study on this paper already shows that tau is highly stretchable. Therefore, when we are considering a dimerized protein model, it is expected for both proteins to be highly stretched before the observation of sliding out of one projection domain out of another. This set of simulations provides new insight into the shear mechanism and sliding threshold of dimerized tau, which is yet non-existent in literature. Figure 3b shows the potential energy trend, and how the potential energy decreases drastically at the separation, irrespective of applied strain rate.

Tau-MT Interaction
To determine the nature of single tau-MT interaction, one N-system MT adopted from Wells et al. [14] was used. One tau protein MT binding site was attached to the surface of the MT. The tensile test results were quite straightforward: the stretching is followed by unfolding in the projection domain of the tau structure, and even at these extreme strain rates, tau had to be stretched for >11 times of its initial length before getting completely separated from MT. The unfolding and stretching manners are similar to the single and dimerized tau tests, that is, the projection domain gets unfolded for a long time, then gets significantly stretched, leading to the eventual separation from the MT surface. As we have seen in the case of the dimerized tau models, we observe earlier separation at a higher strain rate and vice versa for tau-MT models. We can conclude from the two strain rate results that at a certain strain rate range, the tau can be stretched significantly while still being attached to the MT surface. On the other hand, at a higher strain rate, it is not allowed to be stretched in that manner, depicting the effect of strain rate. In the dimerized tau, we have observed the development of higher stress in tau before final separation for lower strain rate. However, in tau-MT interaction, we have not observed any significant difference in the stress-strain trend for the two different strain rates, rather only in the separation stretch. The untangling sub-stages also take place before separation, as we have observed in dimerized tau. Figure 5a shows stress-strain graphs for the applied strain rates, and Figure 5b shows the potential energy graph, which suggests that potential energy decreases significantly at separation, for both cases. Figure 6 shows various significant stages during the pulling of single tau away from the MT surface.
Lastly, we are also interested to find out the strain rate-dependent nature of dimerized tau-MT interaction. We have attached a dimerized tau on the surface of MT in a similar way to the single tau-MT model and applied a high strain rate. This simulation was to compare the relative strength of the tau-MT bond to the tau-tau bond, and as the simulation shows, the tau-tau bond is much stronger than the tau-MT bond, because although the pulled tau was stretched significantly at the higher strain rate (~360% before separation from MT surface), it did not disentangle from the other tau. Rather, the entire dimerized structure got separated from the MT surface, suggesting that the tau-tau bond is stronger than the tau-MT bond. For the lower strain rate, the dimerized tau subunits got significant sliding over each other but eventually got separated from the MT surface before getting entirely separated from each other (~825% before separation). Figure 7 shows that potential energy significantly reduces at separation, and Figure 8 shows different significant stages up to and at separation.

Tau-MT Interaction
To determine the nature of single tau-MT interaction, one N-system MT adopted from Wells et al. [14] was used. One tau protein MT binding site was attached to the surface of the MT. The tensile test results were quite straightforward: the stretching is followed by unfolding in the projection domain of the tau structure, and even at these extreme strain rates, tau had to be stretched for >11 times of its initial length before getting completely separated from MT. The unfolding and stretching manners are similar to the single and dimerized tau tests, that is, the projection domain gets unfolded for a long time, then gets significantly stretched, leading to the eventual separation from the MT surface. As we have seen in the case of the dimerized tau models, we observe earlier separation at a higher strain rate and vice versa for tau-MT models. We can conclude from the two strain rate results that at a certain strain rate range, the tau can be stretched significantly while still being attached to the MT surface. On the other hand, at a higher strain rate, it is not allowed to be stretched in that manner, depicting the effect of strain rate. In the dimerized tau, we have observed the development of higher stress in tau before final separation for lower strain rate. However, in tau-MT interaction, we have not observed any significant difference in the stress-strain trend for the two different strain rates, rather only in the separation stretch. The untangling sub-stages also take place before separation, as we have observed in dimerized tau. Figure 5a shows stress-strain graphs for the applied strain rates, and Figure 5b shows the potential energy graph, which suggests that potential energy decreases significantly at separation, for both cases. Figure 6 shows various significant stages during the pulling of single tau away from the MT surface.
Lastly, we are also interested to find out the strain rate-dependent nature of dimerized tau-MT interaction. We have attached a dimerized tau on the surface of MT in a similar way to the single tau-MT model and applied a high strain rate. This simulation was to compare the relative strength of the tau-MT bond to the tau-tau bond, and as the simulation shows, the tau-tau bond is much stronger than the tau-MT bond, because although the pulled tau was stretched significantly at the higher strain rate (~360% before separation from MT surface), it did not disentangle from the other tau. Rather, the entire dimerized structure got separated from the MT surface, suggesting that the tau-tau bond is stronger than the tau-MT bond. For the lower strain rate, the dimerized tau subunits got significant sliding over each other but eventually got separated from the MT surface before getting entirely separated from each other (~825% before separation). Figure 7 shows that potential energy significantly reduces at separation, and Figure 8 shows different significant stages up to and at separation.      . Figure 7. (a) Stress vs. Strain graph for the projection domain of tau during the pulling at different strain rates. The separation occurs early for higher strain rate (~360%), and much later for lower strain rate (~825%). (b) Potential energy vs. time plot for the dimerized tau-MT system. As expected, potential energy drastically reduces at separation (~360% strain for 2 × 10 9 s −1 , ~825% strain for 1 × 10 9 s −1 ). For both cases, the tau-tau bond is stronger than the tau-MT bond. time plot for the dimerized tau-MT system. As expected, potential energy drastically reduces at separation (~360% strain for 2 × 10 9 s −1 ,~825% strain for 1 × 10 9 s −1 ). For both cases, the tau-tau bond is stronger than the tau-MT bond.

Discussion
In this study, we have analyzed the response of tau protein and tau-MT interaction from a strictly mechanical point of view. We have performed tensile tests on a predicted structure of tau protein to determine the single tau projection domain stiffness, dimerized tau separation stretch, tau-MT separation stretch, and comparison between tau-tau vs. tau-MT bond. For a disordered protein, the confidence score (C-score) of −0.03 has been assumed as reliable in our simulation. The detail of quantification of the reliability for a protein structure predicted by i-TASSER [20] is discussed in the supplementary material (part 1) of this manuscript.
Admittedly, the strain rates that we have applied fall into high to very high range, but it facilitates us to obtain an insight of sub-axonal level response of this particular neural cytoskeletal component. In reality, the tissue level loading might be lower than the cellular level stress, meaning moderate level blow on the head may lead to high-level tissue deformation, which eventually leads to an extreme level of stress and failure in subaxonal level components, supported by recent finite element method (FEM) studies on axon [75,76]. These studies show that axonal level anisotropy and cellular level heterogeneity might play an instrumental role to determine failure criteria of the components and injury level. This level of strain rate is justified in the scenario of cavitation bubble collapse or blast wave exposure, which leads to intensely high stress in the sub-axonal component [8]. This particularly suggests that the cavitation implosion scenario might arise at the sub-

Discussion
In this study, we have analyzed the response of tau protein and tau-MT interaction from a strictly mechanical point of view. We have performed tensile tests on a predicted structure of tau protein to determine the single tau projection domain stiffness, dimerized tau separation stretch, tau-MT separation stretch, and comparison between tau-tau vs. tau-MT bond. For a disordered protein, the confidence score (C-score) of −0.03 has been assumed as reliable in our simulation. The detail of quantification of the reliability for a protein structure predicted by i-TASSER [20] is discussed in the supplementary material (part 1) of this manuscript.
Admittedly, the strain rates that we have applied fall into high to very high range, but it facilitates us to obtain an insight of sub-axonal level response of this particular neural cytoskeletal component. In reality, the tissue level loading might be lower than the cellular level stress, meaning moderate level blow on the head may lead to high-level tissue deformation, which eventually leads to an extreme level of stress and failure in sub-axonal level components, supported by recent finite element method (FEM) studies on axon [75,76]. These studies show that axonal level anisotropy and cellular level heterogeneity might play an instrumental role to determine failure criteria of the components and injury level. This level of strain rate is justified in the scenario of cavitation bubble collapse or blast wave exposure, which leads to intensely high stress in the sub-axonal component [8]. This particularly suggests that the cavitation implosion scenario might arise at the sub-axonal level. The choice of the range of the applied strain rates in this study is in light of this recently published insight, facilitating us capturing high strain rate response of tau and tau-MT interaction.
The model we have implemented attempts to capture the subaxonal level damage when susceptible to a high strain rate. As per the dimerization scheme suggested by Rosenberg et al. [19] and the role of tau protein in manipulating the deformation criteria of microtubules as suggested by Ahmadzadeh et al. [44], we find that the representative model that captures the major deformation aspects consists of a microtubule interacting with a dimerized tau protein. Therefore, in this study, we attempt to not only capture the single tau deformation but also the tau-MT interaction. In the simulation, the model uses the CHARMM force field [60,67], which considers both bonded and non-bonded (chargedependent and positional) interactions. In other words, not only the covalent bonds, angles, and dihedrals are considered, but also the electrostatic and van der Waals interactions are reflected in the interaction, making the force field a very reliable one to extract the deformation characteristics while submerged in explicit water molecules. Furthermore, the tau-MT interaction is dictated by the differently charged portions of tau and MT structure. When the extreme strain rate is applied, the position-dependent, non-bonded interaction force, i.e., van der Waals force also becomes important to determine the separation criteria, and other parameters such as hydrophobicity do not become primary driving attributes to characterize the behavior of the protein deformation [77]. Therefore, this study particularly focuses on electrostatic and van der Waals (or non-bonded) interactions.
Earlier studies that modeled MT response under mechanical loading did not incorporate tau protein with detail mechanical properties, rather studied MT protofibril response where tau protein is considered as a viscoelastic spring [44], the properties of which were adapted from the earlier characterization of pro-and anti-aggregant conformations of tau protein obtained by single-molecule force spectroscopy (SMFS) [78]. However, this study characterizes the response of only mutant conformation of the "repeat domain of tau", not the projection domain which is susceptible to unfolding and stretching under the application of a high strain rate. Therefore, the current study is an effective extension of the previous SMFS study with more comprehensive insight. In general, SMFS studies show the stretch of a molecule from force vs. displacement perspective focusing on the detachment peak force at the contour length of the molecule used, but in the TBI scenario, a more relevant representation is using directional stress vs. strain which can perfectly differentiate between the unfolding and stretching region of a protein under pure mechanical loading.
For dimerized tau, we have re-generated the dimerized tau structure with overlapped projection domains depicted in the study of Rosenberg et al. [19]. This study particularly established the importance of the projection domain which determines the inter-MT distance in an axonal bundle, as the tau-tau interaction is dependent on the length of the projection domain, and interaction with a surface can be adhesive to repulsive (or a combination of both) based on the configuration. This study highlighted the tau-tau interaction and tau-mica interaction and substantiated that the interaction force is a function of projection domain length. The current study undertakes the missing aspect of the study: the effect of shear that the dimerized conformation is susceptible to when undergoing high mechanical stress. Essentially our study shows that the projection domains of dimerized tau proteins are strong against shear and sliding force and that the scenario it undergoes in strictly mechanical loading is highly dynamic, consisting of both unfolding and stretching of both proteins. It further fortifies the observation that the proline-rich region and N terminus combination, which is the projection domain, is highly stretchable while in the dimerized conformation, although the individual tau proteins might reach the failure region earlier. This specific observation shows that mere electrostatic and van der Waals bonds between the negatively charged N terminal region of one tau protein and positively charged proline-rich region of another tau protein are sufficient to withstand mechanical loading to a significant extent. This mechanical behavior of tau is also highlighted in the MT modeling work performed earlier, where one of the key predictions is that tau protein may elongate differently, and according to the position along the MT bundle [44,79]. This high ability of tau protein to stretch aligns with another prediction of the study, that is, tau elongation facilitates the sliding of MT. The increased rigidity of the proline-rich region can also be held responsible for the particular behavior of the tau projection domain as suggested by NMR spectroscopy studies [80].
The tau-MT model in the current study is developed as per the proposal of Chau et al., which shows specific MT binding sites on tau can facilitate bonding with tau binding sites on MT [15,16]. The proposal suggests that tau-MT interaction is highly dynamic and that one MT binding region can interact with one or both subunits of MT (α and β). In our model, the tau protein has been attached to the specific location of the helically arranged protofibril to ensure that tau interacts with the C-terminal site by only R1 or R1-R2 inter-repeat, and with the internal site by the rest of the MT binding sites (R2-R4), which was the scheme suggested by Chau et al. This is a relevant representation of a TBI scenario, where the sub-axonal level stretch on the cytoskeletal component is high for a significant timespan, but the tau-MT bond is sufficiently strong to withstand mechanical load as long as the entire projection domain is not stretched. As stretching does not occur before the unfolding of all the conformed portions of the projection domain, the MT structure instability is not invoked before the directional stress is developed in the MT binding region. The "jaw interaction" between the flanking domain of tau protein and acidic outside of MT surface is also responsible for such strong affinity, as seen between the projection domain and MT (partially), tail domain and MT (fully), suggesting that tau-MT bond is strong irrespective of the intervention of MT binding sites, which is proposed by earlier NMR studies [80]. However, in our case, both the intervention of flanking region and MT binding region are present, which facilitates a stronger bond, as suggested by a separate study [81]. Similar to the dimerized tau system, this set of simulations shows that the interaction between charged portions (MT binding sites of tau and C-terminal as well as the last 1/3rd portion of the C-terminal region excluding the terminal 12 amino acids) provides significant mechanical strength, even when susceptible to high to very high strain rate.
Finally, the dimerized tau-MT interaction system was a representation of a more comprehensive scenario of TBI, in which there is a competition between the mechanical strength of tau-tau and tau-MT bonds. As all of the involved regions, in this case, are charged (proline-rich region of tau, negatively charged N-terminus region of tau, negatively charged outside surface of MT, etc.), we can assume that electrostatic interaction is more important in this case than the van der Waals bonds. The dimerized tau-MT bond shows that it requires the development of more stress in the tau-tau interaction region than required in the tau-MT region, which suggests that proline-rich region interaction with negatively charged projection domain is stronger than negatively charged MT surface interaction with positively charged MT binding sites, and under high strain rate, we can expect tau separation from MT surface. In the case of multiple occasions of tau separation from MT surface in a single injury phenomenon, it may lead to the MT system collapse before the tau system collapse, as tau actually folds into stable conformation upon binding with MT [82].
For high strain rate scenarios, molecular dynamics simulation with implicit solvation technique has been used for axonal cytoskeletal components in recent literature [8]. This is a reasonable compensation when modeling a large representative set (tau-tau dimer or tau-microtubule interaction model susceptible to high strain rate). However, in terms of capturing the mechanical behavior (stress vs. strain), we have not sacrificed accuracy, as in the Supplementary Material Part 2, we have presented a detailed benchmark study (explicit vs. implicit) to substantiate that the stretching fashion is similar after significant strain (~57%). Furthermore, we have also provided the scheme to determine the difference in terms of stress and strain between explicit and implicit solvation techniques. About the possible weakening of the interaction between tau and microtubule, at high strain rates used in this study, the most prominent effect is depicted in the stretching of individual proteins, not the nonbonded interactions between the respective proteins. In LAMMPS implicit solvation, the comparison remains reasonable, as implicit solvation does not render the system to be merely simulated in an empty box, rather it creates a system where the coulombic interaction with a modified formula to incorporate simplified implicit solvation with additional screening [67], which addresses the issue of possible weakening of protein-protein interaction in explicit solvation up to a reasonable extent. Therefore, it can be argued that explicit vs. implicit will have a similar effect when susceptible to a high strain rate. Moreover, the structural viewpoint has already been highlighted in the limitation paragraph of the discussion section. In short, the system we have presented is a reasonable compensation to capture the mechanical behavior of a large representative system of mechanical deformation of axonal cytoskeletal components. Furthermore, steered molecular dynamics (SMD) is also a viable methodology to capture the deformation of biomolecules. In the preliminary model considerations, we have explored the viability of steered molecular dynamics in LAMMPS. However, the current scheme with standard molecular dynamics simulation with the calculation of per atomic stress for the group atoms being pulled provides similar convenience as the steered molecular dynamics option with constant velocity pulling ("cvel" option in LAMMPS). In other words, the standard molecular dynamics is an equally viable candidate for simulation design in our case, as depicted by successful implementation in our works on axonal cytoskeletal components [58,61].
In the unfolding and stretching of the dimerized tau and tau-MT interaction model, there is a collaborative effect of interaction between the charged portions of the structures, and structural coordination at a given time point [2,17]. At the molecular level, it will be interesting to study the residue-specific interactions, which can be achieved by molecular docking or protein-protein interaction-specific servers. In this study, the applied strain rate allows us to only capture the significant points in terms of massive structural changes, such as initiation of unfolding, initiation of stretching, or separation.
As tau protein contains a proline-rich region, it is highly relevant to study its characteristics and assesses its behavior from a neurodegenerative disease perspective. The importance of such proline-rich proteins has been evident in separate neurodegenerative disorder studies [83,84]. These proteins are heavily implicated in neurodegenerative diseases and traumatic brain injury, in which tau is also involved. The PRR is a speculative binding site in proteins, so future studies could include protein biochemistry experiments focused on protein-protein interactions.
Evident from the strength of the tau-MT bond demonstrated in this paper, despite having fast-binding kinetics, we suggest that the intrinsic disorder of tau facilitates this phenomenon. Tau regularly alters its conformation, so its inherent flexibility is a likely source of the protein's ability to remain bound despite the increased strain rate. We suggest tau absorbs the strain throughout its length and relies on the strong capability of the projection domain to remain bound to the microtubule during events where the brain undergoes significant strain and stretch. The strength of the tau-MT bond is particularly important for the field of deformation of axonal cytoskeletal components, where strain and axonal stretch is thought to be a primary mechanism of injury [85]. Future studies might include studying comparing the IDP-substrate bond of other IDPs to see if the high bond strength is a common feature across IDPs as compared with ordered proteins.
Lastly, it is also evident from our MD simulation study that incorporation of physical chemistry perspective (such as posttranslational modification like domain focused or residue focused phosphorylation) along with the mechanical viewpoint is important to obtain comprehensive insight on tau protein behavior. In all the tensile tests, tau protein has shown the dependence on the applied strain rate as single tau behaves as a stiffer material at higher strain rate in both unfolding and stretching region and dimerized tau separation and tau-MT separation stretches have shown strong dependency on strain rate, which suggests the importance of a separate study of viscoelastic characterization of tau protein.
As a computational work, this study has certain limitations, which are highlighted here. The first limitation is the reliability of the i-TASSER predicted structure. At the beginning of the discussion section, it is substantiated that such structure prediction for tau protein is appropriate. However, readers should exercise caution while considering such modeling for implementation in an entirely structure-focused study. As an IDP, tau protein does not possess a defined structure. Therefore, capturing its mechanical behavior when subjected to a high strain rate is a challenge. In this regard, experimental techniques have not been able to provide sufficient insight. Using a computational approach, the reasonable approach is compensating between the disadvantage of implementing predicted structure and limitations in the time scale of equilibrating biomolecular structure in atomistic simulation, and the advantage of being able to capture both single tau, dimerized tau, and tau-MT interaction at atomistic detail. In this regard, we have discussed in supplementary material part 1 that using the i-TASSER predicted structure is an attempt to generate a representative structural conformation of tau protein, which due to the specific parameterization of i-TASSER provides structure containing compaction, and while it can capture the mechanical behavior of the projection domain by appropriate design of tensile tests, it should not be considered as the only viable repertoire to obtain structural insight. In particular, to obtain initial stress versus strain response, this model may not provide an ideal prediction. Secondly, we use the implicit solvation technique in some of our simulations, and therefore, the resulting initial structure of tau protein shows significant compaction. In Supplementary Material Part 2, we provide detailed snapshots of a benchmark comparison between explicit and implicit systems and substantiate that such limitation is overcome at sufficiently stretched tau. Lastly, in supplementary material part 3, it is shown that independently predicted structures show similar mechanical behavior in terms of deformation. In short, there are several limitations in the study from generating tau structure and solvation technique. However, proceeding with such limitations and making reasonable compensation facilitates us providing critical insights into the mechanical deformation of a crucial axonal cytoskeletal component, rendering the study design viable. In addition, due to the intrinsically disordered nature of tau structure, little mechanical insight from the experimental approaches is existent in the literature. Most of the recent experimental maneuvers on tau are limited to obtain biochemical phenomena such as antioxidant properties [86], posttranslational characteristics [87], etc. Even very recent studies use all-atom simulations to determine the structural attributes due to recent advancements of force field development and computational power [88]. Recent near-atomic resolution cryo-EM study, however, provides us inklings regarding tau-MT interaction, which suggests that such interaction may not be limited to one tubulin unit, rather one tau may elongate along the MT protofibril to interact with different tubulin units [52]. In the interest of capturing the representative interaction simplistically, such modes of interactions are not captured in the study, meaning that there is scope to more comprehensively consider different modes of interaction between tau and MT in future studies.
As an additional insight, from an entirely structural point of view, shuffling the sequence of tau will alter the three-dimensional coordination significantly, leading to a different mechanical property. Furthermore, as a biological member of the axonal cytoskeleton, the functionality of tau is also dependent upon the structure and sequence [89], as the sites for post-translational modification regions (such as candidate residues for phosphorylation) will change in position [58], and differently accessible by the potential activator(s). However, this is a highly interesting perspective for drug developers [90].

Conclusions
In this paper, we have computationally determined the stiffness of the projection domain of single tau protein and dimerized tau proteins, and the strength of the tau-MT interface. The necessity of this study has been depicted by a recent review study that pointed out that insight on the high strain rate behavior of tau protein is currently absent in the literature. Due to the length scale associated with the axonal cytoskeleton, the MD simulation approach has been utilized as a viable maneuver. From our MD simulations, the major findings can be summarized as below:

1.
Single tau protein is highly stretchable. It shows unfolding to a great extent before being purely stretched. The unfolding stiffness range is between 50 MPa and 500 MPa while stretching stiffness can be >2 GPa. The stiffness in both regions increases with the increase in the strain rate.

2.
Dimerized tau-tau bond is strong, and the dimer structure does not dissociate before being stretched at >7.5 times of the initial length.

3.
Tau protein can be separated from MT only at a very high strain (>11 times the initial length of tau), and the tau-MT bond is stronger than the MT subunit bond. From the dimerized tau-MT simulation, we have obtained that the tau-tau bond is stronger than the tau-MT bond. We can hypothesize about mechano-chemical events which can trigger MT-tau separation, which has the timescale(s) that MD simulation can capture.

4.
Strain rate affects the separation stretch and developed stress in tau for the dimerized tau model, while it only affects the separation stretch significantly for the tau-MT model. A higher strain rate causes early separation and vice versa.
Bottom-up computational modeling of axons requires insight on the mechanical behavior of individual components, and therefore, this study provides required insight on the strain rate-dependent mechanical behavior of individual tau protein as well as tau-MT interface interaction. In the injury biomechanics area and especially in multiscale axonal cytoskeletal component deformation studies, these findings will play an instrumental role to determine damage criteria at the sub-axonal level. This MD simulation study particularly finds out the sub-axonal level response of axonal cytoskeletal components, which are relevant to the TBI scenario, where nanoscale injury propagates (axonal damage, MT instability, tau unfolding, and stretching, tau-MT separation) due to macroscale impact (head injury).