Kinetic stability analysis of protein assembly on the center manifold around the critical point
 Tatsuaki Tsuruyama^{1}Email author
DOI: 10.1186/s1291801703917
© The Author(s). 2017
Received: 3 June 2016
Accepted: 5 January 2017
Published: 2 February 2017
Abstract
Background
Nonlinear kinetic analysis is a useful method for illustration of the dynamic behavior of cellular biological systems. To date, center manifold theory (CMT) has not been sufficiently applied for stability analysis of biological systems. The aim of this study is to demonstrate the application of CMT to kinetic analysis of protein assembly and disassembly, and to propose a novel framework for nonlinear multiparametric analysis. We propose a protein assembly model with nonlinear kinetics provided by the fluctuation in monomer concentrations during their diffusion.
Results
When the diffusion process of a monomer is selflimited to give kinetics nonlinearity, numerical simulations suggest the probability that the assembly and disassembly oscillate near the critical point. We applied CMT to kinetic analysis of the center manifold around the critical point in detail, and successfully demonstrated bifurcation around the critical point, which explained the observed oscillation.
Conclusions
The stability kinetics of the present model based on CMT illustrates a unique feature of protein assembly, namely nonlinear behavior. Our findings are expected to provide methodology for analysis of biological systems.
Keywords
Protein assembly Nonlinear kinetics FluctuationsBackground
Numerical simulation based upon multiparametric kinetic equations is the principal methodology for the analysis of the behavior of biological systems. Researchers often encounter a number of parameters in the governing equations of the system. Here, we introduce the center manifold theory (CMT) for simplification of the study of dynamic biological systems. CMT provides mathematical prescription for carrying out reduction of the number of parameters near the steady state, as well as information regarding the stability of the steady state. As a result, simulation is oriented to illustrate behavior around the critical point, at which system behavior drastically changes in the qualitative structure. The observable change is termed bifurcation, and the threshold values of the parameters are referred to as critical values or bifurcation values. The aim of this study was to provide a simple algorithm for the application of CMT to multiparametric kinetic equations, in order to clearly illustrate the behavior of the biological system. The CMT has been applied to the LotkaVolterra model of predator–prey system to provide important simulation results [1, 2]. In addition, several pioneering studies have applied CMT to neural network analysis [3]. Timedelay and diffusive effects play important roles in bifurcation phenomena [1, 4]. However, to date, there are few applications of the CMT to biochemical reaction models. We previously reported a model of cell signaling systems using nonlinear kinetics and demonstrated the phase transition phenomenon via a numerical simulation [5].
Pivotal proteinprotein interactions during cytoskeleton formation were selected as the application model for the present CMT method. Among the interactions between protein monomers, tubulin and actin polymerization are wellknown events that have been analyzed using the numerical method [6–10]. The physical robustness of the cytoskeleton is based on the biophysical properties of actin and tubulin. In particular, various mathematical models have been proposed to explain the kinetic behavior of tubulin assembly [6–11]. A theory of polymerization of macromolecules has been established on the basis of the kinetic model of aggregation [12, 13]. Oosawa and Asakura previously reported that polymerization is similar to micelle formation or crystallization, and that there is a critical monomer concentration above which monomers effectively polymerize. The authors additionally suggested that the nucleation step represents the ratelimiting step for polymerization. Nucleation and growth occur in parallel during the progression of polymerization. There is a gap in free energy change between initial nucleation and progression of linear polymerization [13]. The stable nucleus for polymerization consists of trimers or tetramers, and the growth of aggregates through elongation/ dissociation follows the formation of a thermodynamically unfavorable size of the nucleus. In the current study, we focused on the polymerization in the absence of de novo nucleation and the interaction between polymer and monomer (PM) interaction.
For stable growth, the lifespan of tubules is controlled by a guanidine triphosphate (GTP)cap that forms at their ends [14]. The structure and motility of growing tubules is influenced by intrapolymeric Brownian motion and fluctuation; this provides elasticity to the microtubules [15]. Polymerization/depolymerization is controlled by binding of adenosine triphosphate (ATP)/GTP, resulting in the assembly of monomeric proteins. The intermittent transition between slow growth and rapid shrinkage in polymeric assemblies of microtubules is termed dynamic instability [14]. Numerous models have been proposed to explain this instability; in particular, Zapperi and Mahadevan successfully identified two parameters: a structural mechanical parameter that characterizes the ratio of longitudinal to lateral interactions in an assembly, and a kinetic parameter that characterizes the ratio of timescales for growth and conformation change. These parameters serve to demarcate a region of uninterrupted growth from that of collapse [16].
Here, the first item on the right represents the diffusion rate. The second item, f (c _{ i }), denotes the function of kinetic rate of reactions other than the diffusion process. k _{ i } represents a coefficient.
Methods
Numerical simulation Numerical calculations were performed using Mathematica 8 (Wolfram Research, Inc., Champaign, IL).
Results
General formulation of an assembly
The model consists of several steps: (i) the monomer achieves an interactive state by binding a cofactor (ATP/GTP) that provides the monomer with the ability to interact; (ii) the monomer itself possesses the ability to hydrolyze the cofactor and lose assembly activity; (iii) the monomer has the ability to exchange the inactive hydrolyzed cofactor (ADP/GDP) with an active nonhydrolyzed one; and (iv) ATP/GTP are supplied continuously from the external environment. The second requirement indicates a selflimiting property of the monomer that causes dynamic instability during monomermonomer interaction. When examining protein interaction kinetics, analysis of the fluctuation in monomer concentrations was performed using Mathematica 9.
Protein interaction kinetics
In the above equations, k _{ 1 } and k _{ 2 } represent the kinetic coefficients for the addition of the monomer to the oligomer and the release of the monomer from the oligomer, respectively.
In the above approximation, we omitted D _{ X } D _{ Y } and k _{ 3 } as the diffusion coefficients and the direct conversion rate of X into Z is small.
Fluctuation of diffusion coefficient
In Eq. (14), the subscript ‘e’ signifies values at the steady state.
Therefore, the fluctuation y may be represented using −x−z. The fluctuation kinetics are thus provided by two parameters, namely x and z.
Eqs. (28) and (29) represent a master equation for the application of CMT.
Calculus simulation of concentration oscillations
For analysis of the behavior of the system, including multiparameters, the examination of the linearization of behavior of the system near a steady state provides insights into the qualitative behavior of the system in the neighborhood of the point. In particular, the eigenvalues of the linear part of the governing kinetic equations enable determination of the stability of the system behavior. CMT is a rigorous formulation of this observation that enables the reduction of a large number of parameters [22].
As a result, the fluctuations oscillate between decrease and increase in monomer concentrations, as shown in Fig. 2. When p <p _{ c }, the fluctuation was found tobe attenuated (Fig. 2d) and the monomer concentration reached a plateau. However, when p >p _{ c }, the fluctuation was found to diverge (Fig. 2f).
Evaluation of model stability using the center manifold around the equilibrium state
The eigenvalues of the Jacobian matrix, λ, in (33) are 0 and 2.9 × 10^{−4}. The given center manifold is an invariant manifold that is a tangent space of the center subspace, which is an eigenspace when the eigenvalue is equivalent to zero. The behavior of the fluctuation is complex when the real part of the eigenvalue is equivalent to zero. The above result in (36) shows that it is systematically necessary to analyze the behavior of the given system on the center manifold [22]. In order to analyze the behavior of the system, we investigated whether the change of the value in p around the critical value p _{c} gives u that satisfies du/dt = 0. When the two values of u are given, i.e., bifurcation of the system is shown, and oscillation and/or other interesting behaviors may be predicted.
As a result, as we described v and x had two amplitudes in (47) demonstrating the oscillation of the fluctuation by bifurcation in vε plane (Fig. 2). Thus, stability analysis enables prediction of the behavior of the fluctuation around the critical point of the protein assembly system.
Discussion
In this work, we presented a model for protein assembly kinetics and analyzed the stability around the critical point using CMT. The nonlinear kinetic equations include three parameters (X, Y, and Z); however, only two are independent. In the simulations, ATP/GTP or ADP/GDPbinding monomers periodically exhibit an oscillation between assembly and disassembly. This accurately reflects the microtubule kinetics showing unstable assembly [8].
To the best of our knowledge, this is one of the first reports on the application of CMT to the analysis of biological reaction systems [8]. The fluctuation of monomer concentrations was subjected to a perturbation expansion using a minimal increase in the supply of ATP/GTP near the concentration at the critical point. This mathematical method precisely treats nonlinear and multiparameter systems around the critical point. The fluctuation kinetics is expected to change from convergence to divergence of the concentration fluctuation of the monomer, i.e., from stable to unstable around the critical point, as shown in Fig. 2. Because of this high sensitivity to the concentration of ATP/GTP, protein assembly is dynamically regulated by minimal changes in the supply of ATP/GTP, which in turn is subject to metabolic control. Via modeling of microtubule growth at the mesoscopic scale, Zapperi et al. showed the time course of transition between slow growth and rapid shrinkage during microtubule polymerization [16]. The present simulation may explain microscopic tubulin oligomerization oscillations during the initial steps of microtubule assembly. In addition, the present model may explain the transition from microscopic oligomerization and aggregation to mesoscopic scale assembly. The quantitative evaluation of the theoretical basis of protein assembly requires further investigation through experimental studies.
The present center manifold analysis enables elucidation of detailed behavior around the steady state and oscillatory dynamics of protein monomer concentration. In the current study, we further developed the mathematical framework using CMT and aimed to describe Hopfbifurcation around the steady state, through the center manifold analysis, in a simple model. Coveney et al. have described a detailed model of protein assembly, including nucleation, its catalysis, and inhibition processes and performed a kinetic analysis of the initial nucleation process [23, 24]. The kinetic model of monomeroligomerization or nucleation requires multiple concentrations that describe variable oligomer and nucleation. As shown by Coveney et al., it was challenging to predict the behavior of the system using a multiparametric (dimensional) centermanifold on the model. In the current study, we utilized a monomeric parameter and showed bifurcation of the system around the critical point. Therefore, CMT in a simple model serves to reduce the dimensions of the system to signal dimensions, as shown in this study. We expect that the theoretical framework in the current study provides a general theory of protein assembly kinetics and signal transduction [5, 25].
The analysis of growth kinetics of polymerization, according to Oosawa’s model, has recently been reported by Michaels et al. [12]. The authors focused primarily on the dynamic phase of protein polymerization. As nucleation and polymerization to the nucleus proceeds in parallel, the analysis requires a detailed kinetic model of interaction between the nucleation and polymerization process [13, 14]. However, after the dynamic phase and before the plateau phase of polymerization, PM interactions are dominant during signal transduction. The present analysis illustrates the dynamics of cytoplasm in the stable state, and the corresponding influence on cell motility.
The present simulation was applied to such a quasistatistic state, and the results revealed a possibility that oscillation of monomer concentration may occur when the ATP/GTP concentration exceeds the critical concentration. The calculated critical concentration of ATP/GTP, based on Hopfbifurcation in (46) and amplitude of the fluctuation, coincided well with the amplitude obtained via the present simulation. The consistency in values in the simulation is important for verification. The periodic change in concentration may contribute to the coherently spatialperiodic viscosity and subsequently to contraction and elongation during cell movement. A recent study demonstrated the role of cytokeratin in determining keratinocyte motility and shape [26] and experimental method has greatly developed [27]. Structural components of cells determine nonlinear cellular structural behavior and the contribution of various cell components to stability in response to mechanical stimuli. The cytoskeleton plays key roles in determining cellular stiffness. Our model captures nonlinear structural behaviors including variable compliance along the cell surface and resistance to pullout force [28]. The role of the microtubules in dynamic behavior may be investigated from the viewpoint of cell geometries. Measurement of the oscillation and determination of the critical concentration of ATP/GTP may reveal physical properties such as elasticity and compressibility.
Conclusion
Our model is expected to be useful for computing biophysical behavior in response to minute changes in GTP/ATP concentration using fluorescence intensity meter in twodimensional cell geometries. In addition, the present model is expected to be suitable for use in algorithms for simulation of metabolic processes. Although further experimental studies are necessary for verification, our findings show that the current nonlinear model of dynamic instability analysis captures the nonlinear behaviors of cellular chemical and mechanical responses.
Abbreviations
 ATP:

Adenosine triphosphate ATP
 CMT:

Center manifold theory
 GTP:

Guanidine triphosphate
 PM:

Polymer and monomer (PM)
Declarations
Acknowledgements
We are very thankful for Prof. Kenichi Yoshikawa, Dr. Masatoshi Ichikawa, and Prof. Masayuki Imai.
Funding
This work was supported by the Ministry of Education, Culture, Sport, Science and Technology, Japan, under the project name “Synergy of Fluctuation and Structure”. Project number 2502. The funders had no role in the study design, data collection and analysis, decision to publish, or preparation of the manuscript.
Availability of data and materials
The dataset supporting the conclusions of this article is included within the article.
Author’ contributions
TT designed the study, implemented the final model, and wrote the manuscripts.
Competing interests
The authors declare that they have no competing interests.
Consent for publication
Not applicable.
Ethics approval and consent to participate
Not applicable.
Open AccessThis article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated.
Authors’ Affiliations
References
 Chang X, Wei J. Stability and Hopf bifurcation in a diffusive predatorprey system incorporating a prey refuge. Math Biosci Eng. 2013;10:979–96.View ArticlePubMedGoogle Scholar
 Zhang X, Zhao H. Bifurcation and optimal harvesting of a diffusive predatorprey system with delays and interval biological parameters. J Theor Biol. 2014;363:390–403.View ArticlePubMedGoogle Scholar
 Xiao M, Zheng WX, Cao J. Hopf bifurcation of an (n + 1) neuron bidirectional associative memory neural network model with delays. IEEE Trans Neural Netw Learn Syst. 2013;24:118–32.View ArticlePubMedGoogle Scholar
 Yamaguchi I, Ogawa Y, Jimbo Y, Nakao H, Kotani K. Reduction theories elucidate the origins of complex biological rhythms generated by interacting delayinduced oscillations. PLoS One. 2011;6:e26497.View ArticlePubMedPubMed CentralGoogle Scholar
 Tsuruyama T. A model of cell biological signaling predicts a phase transition of signaling and provides mathematical formulae. PLoS One. 2014; (in press).
 Hazra P, Inoue K, Laan W, Hellingwerf KJ, Terazima M. Tetramer formation kinetics in the signaling state of AppA monitored by timeresolved diffusion. Biophys J. 2006;91:654–61.View ArticlePubMedPubMed CentralGoogle Scholar
 Wu Z, Wang HW, Mu W, Ouyang Z, Nogales E, Xing J. Simulations of tubulin sheet polymers as possible structural intermediates in microtubule assembly. PLoS One. 2009;4:e7291.View ArticlePubMedPubMed CentralGoogle Scholar
 VanBuren V, Cassimeris L, Odde DJ. Mechanochemical model of microtubule structure and selfassembly kinetics. Biophys J. 2005;89:2911–26.View ArticlePubMedPubMed CentralGoogle Scholar
 Symmons MF, Martin SR, Bayley PM. Dynamic properties of nucleated microtubules: GTP utilisation in the subcritical concentration regime. J Cell Sci. 1996;109:2755–66.PubMedGoogle Scholar
 Voter WA, Erickson HP. The kinetics of microtubule assembly. Evidence for a twostage nucleation mechanism. J Biol Chem. 1984;25:10430–8.Google Scholar
 Zilberman M, Sofer M. A mathematical model for predicting controlled release of bioactive agents from composite fiber structures. J Biomed Mater Res A. 2007;80:679–86.View ArticlePubMedGoogle Scholar
 Oosawa F, Kasai M. A theory of linear and helical aggregations of macromolecules. J Mol Biol. 1962;4:10–21.View ArticlePubMedGoogle Scholar
 Michaels TC, Garcia GA, Knowles TP. Asymptotic solutions of the Oosawa model for the length distribution of biofilaments. J Chem Phys. 2014;140:194906.View ArticlePubMedGoogle Scholar
 Chretien D, Jainosi I, Taveau JC, Flyvbjerg H. Microtubule’s conformational cap. Cell Struct Funct. 1999;24:299–303.View ArticlePubMedGoogle Scholar
 Oosawa F, Asakura S. Thermodynamics of the Polymerisation of Proteins. New York and London: Acdemic Press; 1975. p. 204.
 Zapperi S, Mahadevan L. Dynamic instability of a growing adsorbed polymorphic filament. Biophys J. 2011;101(2):267–75.View ArticlePubMedPubMed CentralGoogle Scholar
 Wustner D, Solanko LM, Lund FW, Sage D, Schroll HJ, Lomholt MA. Quantitative fluorescence loss in photobleaching for analysis of protein transport and aggregation. BMC Bioinformatics. 2012;13:296.View ArticlePubMedPubMed CentralGoogle Scholar
 Dorsaz N, De Michele C, Piazza F, Foffi G. Inertial effects in diffusionlimited reactions. J Phys Condens Matter. 2010;22:104116.View ArticlePubMedGoogle Scholar
 Kasche V, de Boer M, Lazo C, Gad M. Direct observation of intraparticle equilibration and the ratelimiting step in adsorption of proteins in chromatographic adsorbents with confocal laser scanning microscopy. J Chromatogr B Analyt Technol Biomed Life Sci. 2003;790:115–29.View ArticlePubMedGoogle Scholar
 O’Leary TJ. Concentration dependence of protein diffusion. Biophys J. 1987;52:137–9.View ArticlePubMedPubMed CentralGoogle Scholar
 Kenneth H. A diffusion model with a concentrationdependent diffusion coefficient for describing water movement in legumes during soaking. J Food Sci. 1983;48:618–23.View ArticleGoogle Scholar
 Guckenheimer J, Holmes PJ. Nonlinear oscillations, dynamical systems, and bifurcations of vector fields. 1st ed. New York: Springer; 1983. p. 1–459.Google Scholar
 Wattis JAD, Coveney PV. Analysis of a generalized beckerdoring model of selfreproducing micelles. Proc T Soc Lond A. 1996;452:2079–102.View ArticleGoogle Scholar
 Wattis JAD, Coveney PV. Mesoscopic models of nucleation and growth processes : a challenge to experiment. Phys Chem Chem Phys. 1999;1:2163–76.View ArticleGoogle Scholar
 Babu CVS, Song EJ, Yoo YS. Modeling and simulation in signal transduction pathways: a systems biology approach. Biochimie. 2006;88:277–83.View ArticleGoogle Scholar
 Nakata T, Okimura C, Mizuno T, Iwadate Y. The role of stress fibers in the shape determination mechanism of fish keratocytes. Biophys J. 2016;110:481–92.View ArticlePubMedPubMed CentralGoogle Scholar
 McGarry JG, Prendergast PJ. A threedimensional finite element model of an adherent eukaryotic cell. Eur Cell Mater. 2004;7:27–33.View ArticlePubMedGoogle Scholar
 Burk AS, Monzel C, Yoshikawa HY, Wuchter P, Saffrich R, Eckstein V, et al. Quantifying adhesion mechanisms and dynamics of human hematopoietic stem and progenitor cells. Sci Rep. 2015;5:9370.View ArticlePubMedGoogle Scholar