Detection of biological switches using the method of Gröebner bases.

BACKGROUND: Bistability and ability to switch between two stable states is the hallmark of cellular responses. Cellular signaling pathways often contain bistable switches that regulate the transmission of the extracellular information to the nucleus where important biological functions are executed. RESULTS: In this work we show how the method of Gröebner bases can be used to detect bistability and output switchability. The method of Gröebner bases can be seen as a multivariate, non-linear generalization of the Gaussian elimination for linear systems which conveniently seperates the variables and drastically simplifies the simultaneous solution of polynomial equations. A necessary condition for fixed-point state bistability is for the Gröbner basis to have three distinct solutions for the state. A sufficient condition is provided by the eigenvalues of the local Jacobians. We also introduce the concept of output switchability which is defined as the ability of an output of a bistable system to switch between two different stable steady-state values. It is shown that bistability does not necessarily guarantee switchability of every state variable of the system. We further show that, for a bistable system, the necessary conditions for output switchability can be derived using the Gröebner basis. The theoretical results are incorporated into an analysis procedure and applied to several systems including the AKT (Protein kinase B), RAS (Rat Sarcoma) and MAPK (Mitogen-activated protein kinase) signal transduction pathways. Results demonstrate that the Gröebner bases can be conveniently used to analyze biological switches by simultaneously detecting bistability and output switchability. CONCLUSION: The Gröebner bases provides a novel methodology to analyze bistability. Results clarify the distinction between bistability and output switchability which is lacking in the literature. We have shown that theoretically, it is possible to have an output subspace of an n-dimensional bistable system where certain variables cannot switch. It is possible to construct such systems as we have done with two reaction networks.

Dynamics and control of the ERK signaling pathway: Sensitivity, bistability, and oscillations.

Cell signaling is the process by which extracellular information is transmitted into the cell to perform useful biological functions. The ERK (extracellular-signal-regulated kinase) signaling controls several cellular processes such as cell growth, proliferation, differentiation and apoptosis. The ERK signaling pathway considered in this work starts with an extracellular stimulus and ends with activated (double phosphorylated) ERK which gets translocated into the nucleus. We model and analyze this complex pathway by decomposing it into three functional subsystems. The first subsystem spans the initial part of the pathway from the extracellular growth factor to the formation of the SOS complex, ShC-Grb2-SOS. The second subsystem includes the activation of Ras which is mediated by the SOS complex. This is followed by the MAPK subsystem (or the Raf-MEK-ERK pathway) which produces the double phosphorylated ERK upon being activated by Ras. Although separate models exist in the literature at the subsystems level, a comprehensive model for the complete system including the important regulatory feedback loops is missing. Our dynamic model combines the existing subsystem models and studies their steady-state and dynamic interactions under feedback. We establish conditions under which bistability and oscillations exist for this important pathway. In particular, we show how the negative and positive feedback loops affect the dynamic characteristics that determine the cellular outcome.

Optimum folding pathways for growing protein chains.

The folding of a protein is studied as it grows residue by residue from the N-terminus and enters an environment that stabilizes the folded state. This mode of folding of a growing chain is different from refolding where the full chain folds from a disordered initial configuration to the native state. We propose a sequential dynamic optimization method that computes the evolution of optimum folding pathways as amino acid residues are added to the peptide chain one by one. The dynamic optimization formulation is deterministic and uses Newton's equations of motion and a Go-type potential that establishes the native contacts and excluded volume effects. The method predicts the optimal energy-minimizing path among all the alternative feasible pathways. As two examples, the folding of the chicken villin headpiece, a 36-residue protein, and chymotrypsin inhibitor 2 (CI2), a 64-residue protein, are studied. Results on the villin headpiece show significant differences from the refolding of the same chain studied previously. Results on CI2 mostly agree with the results of refolding experiments and computational work.

Dynamic Modeling and Analysis of the Cross-Talk between Insulin/AKT and MAPK/ERK Signaling Pathways.

Feedback loops play a key role in the regulation of the complex interactions in signal transduction networks. By studying the network of interactions among the biomolecules present in signaling pathways at the systems level, it is possible to understand how the biological functions are regulated and how the diseases emerge from their deregulations. This paper identifies the key feedback loops involved in the cross-talk among the insulin-AKT and MAPK/ERK signaling pathways. We developed a mathematical model that can be used to study the steady-state and dynamic behavior of the interactions among these two important signaling pathways. Modeling analysis and simulation case studies identify the key interaction parameters and the feedback loops that determine the normal and disease phenotypes.

Folding dynamics of proteins from denatured to native state: principal component analysis.

Several trajectories starting from random configurations and ending in the native state for chymotrypsin inhibitor 2, CI2, are generated using a Go-type model where the backbone torsional angles execute random jumps on which a drift towards their native values is superposed. Bond lengths and bond angles are kept fixed, and the size of the backbone atoms and side groups are recognized. The large datasets obtained are analyzed using a particular type of principal component analysis known as Karhunen-Loeve expansion (KLE). Trajectories are decomposed separately into modes in residue space and time space. General features of different folding trajectories are compared in the modal space and relationships between the structure of CI2 and its folding dynamics are obtained. Dynamic scaling and order reduction of the folding trajectories are discussed. A continuous wavelet transform is used to decompose the nonstationary folding trajectories into windows exhibiting different features of folding dynamics. Analysis of correlations confirms the known two-state nature of folding of CI2. All of the conserved residues of the protein are shown to be stationary in the small modes of the residue space. The sequential nature of folding is shown by examining the slow modes of the trajectories. The present model of protein folding dynamics is compared with the simple Rouse model of polymer dynamics. Principal component analysis is shown to be a very effective tool for the characterization of the general folding features of proteins.

Regulatory networks and complex interactions between the insulin and angiotensin II signalling systems: models and implications for hypertension and diabetes.

The cross-talk between insulin and angiotensin II signalling pathways plays a significant role in the co-occurrence of diabetes and hypertension. We developed a mathematical model of the system of interactions among the biomolecules that are involved in the cross-talk between the insulin and angiotensin II signalling pathways. We have identified several feedback structures that regulate the dynamic behavior of the individual signalling pathways and their interactions. Different scenarios are simulated and dominant steady-state, dynamic and stability characteristics are revealed. The proposed mechanistic model describes how angiotensin II inhibits the actions of insulin and impairs the insulin-mediated vasodilation. The model also predicts that poor glycaemic control induced by diabetes contributes to hypertension by activating the renin angiotensin aystem.

Combining optimal control theory and molecular dynamics for protein folding.

A new method to develop low-energy folding routes for proteins is presented. The novel aspect of the proposed approach is the synergistic use of optimal control theory with Molecular Dynamics (MD). In the first step of the method, optimal control theory is employed to compute the force field and the optimal folding trajectory for the Cα atoms of a Coarse-Grained (CG) protein model. The solution of this CG optimization provides an harmonic approximation of the true potential energy surface around the native state. In the next step CG optimization guides the MD simulation by specifying the optimal target positions for the Cα atoms. In turn, MD simulation provides an all-atom conformation whose Cα positions match closely the reference target positions determined by CG optimization. This is accomplished by Targeted Molecular Dynamics (TMD) which uses a bias potential or harmonic restraint in addition to the usual MD potential. Folding is a dynamical process and as such residues make different contacts during the course of folding. Therefore CG optimization has to be reinitialized and repeated over time to accomodate these important changes. At each sampled folding time, the active contacts among the residues are recalculated based on the all-atom conformation obtained from MD. Using the new set of contacts, the CG potential is updated and the CG optimal trajectory for the Cα atoms is recomputed. This is followed by MD. Implementation of this repetitive CG optimization-MD simulation cycle generates the folding trajectory. Simulations on a model protein Villin demonstrate the utility of the method. Since the method is founded on the general tools of optimal control theory and MD without any restrictions, it is widely applicable to other systems. It can be easily implemented with available MD software packages.

Prediction of optimal folding routes of proteins that satisfy the principle of lowest entropy loss: dynamic contact maps and optimal control.

An optimization model is introduced in which proteins try to evade high energy regions of the folding landscape, and prefer low entropy loss routes during folding. We make use of the framework of optimal control whose convenient solution provides practical and useful insight into the sequence of events during folding. We assume that the native state is available. As the protein folds, it makes different set of contacts at different folding steps. The dynamic contact map is constructed from these contacts. The topology of the dynamic contact map changes during the course of folding and this information is utilized in the dynamic optimization model. The solution is obtained using the optimal control theory. We show that the optimal solution can be cast into the form of a Gaussian Network that governs the optimal folding dynamics. Simulation results on three examples (CI2, Sso7d and Villin) show that folding starts by the formation of local clusters. Non-local clusters generally require the formation of several local clusters. Non-local clusters form cooperatively and not sequentially. We also observe that the optimal controller prefers "zipping" or small loop closure steps during folding. The folding routes predicted by the proposed method bear strong resemblance to the results in the literature.

Determination of pair-wise inter-residue interaction forces from folding pathways and their implementation in coarse-grained folding prediction.

Interaction forces among residue pairs are determined from optimum folding pathways along which a protein represented as a coarse-grained chain of alpha-carbons goes from different initial configurations to a known native state. A dynamic optimization approach is employed that uses the coarse-grained model to compute the optimal folding pathways. The pair-wise interaction forces obtained in this manner are incorporated into the coarse-grained model which is then simulated to fold the protein from a new set of initial configurations in a predictive way. We show that the folding pathways predicted in this manner are near-optimal. We applied the technique to the secondary structures: helix and beta-sheet.

Optimum folding pathways of proteins: their determination and properties.

We develop a dynamic optimization technique for determining optimum folding pathways of proteins starting from different initial configurations. A coarse-grained Go model is used. Forces acting on each bead are (i) the friction force, (ii) forces from bond length constraints, (iii) excluded volume constraints, and (iv) attractive forces between residue pairs that are in contact in the native state. An objective function is defined as the total attractive energy between nonbonded residues, which are neighbors in the native state. The objective function is minimized over all feasible paths, satisfying bond length and excluded volume constraints. The optimization problem is nonconvex and contains a large number of constraints. An augmented Lagrangian method with a penalty barrier function was used to solve the problem. The method is applied to a 36-residue protein, chicken villin headpiece. Sequences of events during folding of the protein are determined for various pathways and analyzed. The relative time scales are compared and scaled according to experimentally measured events. Formation times of the helices, turn, and the loop agree with experimental data. We obtain the overall folding time of the protein in the range of 600 ns-1.2 micros that is smaller than the experimental result of 4-5 micros, showing that the optimal folding times that we obtain may be possible lower bounds. Time dependent variables during folding and energies associated with short- and long-range interactions between secondary structures are analyzed in modal space using Karhunen-Loeve expansion.