1 INTRODUCTION
Over the past several decades, mathematical modeling has become an important tool in biological research. Due to the rapid expansion of biological knowledge, kinetic modeling has become a realistic goal, particularly for experimentally well-characterized systems.
Kinetic models have been essential for the progress of the fields of systems and synthetic biology, where scientists have successfully applied these models to explain many biological phenomena. For example, researchers recently constructed a whole cell kinetic model of the bacterium
Mycoplasma genitalium, which was then used to predict phenotypic features of the bacterium [
1]; kinetic modeling and synthetic circuits were recently used to successfully examine the role of stochastic pulses in regulating the stress response in bacteria [
2], quantitatively understand quorum sensing in bacteria [
3] and create genetic oscillators [
4].
A kinetic model essentially represents a mathematical integration of existing data and mechanisms on a particular system, and is useful in a number of ways. A kinetic model can be used to test the consistency in the experimental data or mechanisms [
5], provide mechanistic explanations for counter-intuitive observations [
6], to facilitate the formulation of experimentally testable hypotheses [
7] or to provide insight into emergent properties, such as robustness [
3,
8,
9], which may be otherwise difficult to grasp intuitively. Kinetic models are often simulated using various numerical algorithms. Because this process is computationally intensive, model simulations using those algorithms are usually done using modeling software.
Many of the latest annotated kinetic models submitted to the BioModels Database [
10] use either COPASI, MATLAB, or both for modeling or simulation needs. Although these pieces of software are powerful for model creation and simulation, they have a high learning threshold. This high learning threshold makes it difficult for beginning modelers with little or no programming experience to start modeling gene regulatory networks. To alleviate this problem and facilitate modeling, many GUI programs have been created, such as OpenAlea [
11], CellDesigner [
12], VCell [
13] and E-Cell [
14]. However these modeling programs have some limitations. For example, OpenAlea does not cater to general modeling needs and is specialized for solving plant modeling problems. CellDesigner allows GUI modeling of biochemical networks, but it requires an external mathematical package (COPASI) for simulating networks or fitting parameters. Both E-Cell and VCell allow simulation, GUI modeling, and provide very high functionality, but they have a high learning threshold, making them hard to intuitively use for beginners. Many of the other programs are generally available only as plugins and require some other main package or software for either simulating or modeling a network. There are also proprietary tools available such as Monolix [
15], which allow both modeling and simulation capabilities in a single user friendly package, but these are not openly accessible, hampering the facilitation of teaching and research.
Therefore, there is a need for an open, user-friendly, and GUI-based software package with a low learning threshold for both simulating and modeling synthetic and natural gene regulatory networks. Dynetica 2.0 caters to this need. To increase the user friendliness and utility of Dynetica 2.0 over that of version 1.0 [
16], we have incorporated many new features, such as an improved GUI, expression variables, modular organization, support for import and export of models in the form of MATLAB code and Systems Biology Markup Language (SBML) models, improved sensitivity analysis, parameter searching, and many bug fixes.
To allow greater adoption and usage of Dynetica 2.0 by the synthetic and systems biology community, we have released the source code of Dynetica 2.0 at under the GPL Version 2 license.
2 RESULTS
A general overview of the design paradigm of Dynetica is shown in Figure 1. Models in Dynetica 2.0 are represented by the substances and reactions that make up the system, as well as the parameters and expressions that govern them. Models are drawn in a graph layout and can be compartmentalized into individual modules. These modules can be copied, saved, and imported independently, thus facilitating the construction of large systems from basic subcomponents.
Simulations may be run using any of the following algorithms: Runge-Kutta (Fourth Order, Fixed Time Step), Runge-Kutta-Fehlberg (Fourth Order, Variable Time Step), Gillespie algorithms (both First Reaction and Direct Method implementations), Gibson algorithm, and the Euler-Maruyama method for stochastic differential equations.
Additionally, Dynetica 2.0 provides a number of useful tools for analyzing the behavior of models. Sensitivity analysis allows the user to examine how the time courses or characteristic metrics of a system change with perturbations in the system’s parameters or initial conditions. Repeated stochastic simulations allow the user to visualize the distribution of a system’s substance values across many trials. Parameter searching perturbs the system’s parameter values in order to match the system’s output to some user-specified behaviors.
Finally, Dynetica 2.0 supports the import and export of models in SBML format, and provides tools for exporting models to either deterministic or stochastic MATLAB scripts, allowing the user to perform any analyses not covered by Dynetica.
In comparison to the previous version, the areas of improvement in Dynetica 2.0 can be summarized in four segments: interface enhancements, sensitivity analysis support, genetic parameter searching, and import/export functions.
2.1 Interface
The three primary interface enhancements in Dynetica 2.0 are expression variable support, modular system construction, and import/export support.
2.1.1 Expression Variables
An expression variable is an entity in Dynetica 2.0 that can be represented by a mathematical or Boolean function of substance concentrations, parameters, and time. Expression variables can be used in reaction kinetics equations just like substance concentrations or parameter values. Because of this flexibility, they can be used as time-dependent system inputs, intermediate expressions, and output variables or indicators.
2.1.2 Modules
When building large or complex systems in Dynetica, it often becomes difficult to keep track of the various components and how they interact. One way to mitigate this problem is by taking advantage of the modular nature of many of these large networks. Although the total number of components may be large, most components only interact with a small subset of all of the components. The large networks tend to form small sub-networks, which then interact with one another via a handful of interconnects.
Dynetica 2.0 allows users to take advantage of this phenomenon by creating such sub-networks, called “modules.” Figure 2, displaying a model for programmed altruistic death, exemplifies how Dynetica modules can be harnessed to enhance conceptual simplicity. The substances comprising the cheater and cooperator modules, on the bottom right and bottom left corners of the network graph, respectively, have each been isolated into their own module. Another advantage of such functionally defined modules is that they can be easily reused across or within models, speeding up new model construction. Further details regarding the use of modular systems in Dynetica can be found in section 1 of the Supplemental Materials.
2.1.3 Import and export
A major improvement in Dynetica 2.0 over the previous version is its ability to export models in SBML or MATLAB. The earlier version of Dynetica only allowed models to be stored in the native DYN format used by Dynetica. This made it difficult for models created in Dynetica to be further extended using more advanced modeling software, as they would have to be recreated in the new modeling platform. This new feature allows models created in Dynetica to be imported by other software packages that recognize the SBML format or MATLAB code.
In addition, Dynetica 2.0 allows the user to merge existing models with each other or add one model to another as a module. These features greatly facilitates creation of new models because it allows existing models to be reused or extended. Moreover, these new features allow the user to quickly create and test synthetic circuits by combining circuit elements from a library. Another advantage of using this new feature is that the initial MATLAB or SBML code for a model can be generated much more quickly and robustly, since the user will not have to rewrite the entire code when the model is modified or extended. Moreover, parameter values and constants are also preserved during the merger.
2.2 Sensitivity Analysis
2.2.1 Metrics
Metrics are the key units used to analyze the behavior of a system in Dynetica 2.0. A metric is a value that describes some characteristic of a single time course. These are used in both sensitivity analysis and parameter searching.
The metrics implemented in Dynetica 2.0 include:
• Final, min, and max value;
• Range;
• Maximum rate;
• Time to reach a user-specified fraction of steady state value;
• Area under curve;
• Peak frequency in the power spectrum;
• Correlation coefficient when compared to another substance or expression variable in the same system.
2.2.2 Simulation types
Dynetica 2.0 provides the user the ability to explore how the behavior of the system varies with changes in either parameter values or initial substance concentrations, through its sensitivity analysis tool. The user selects a parameter or substance from the model to act as the independent variable, and defines any number of metrics to evaluate. The simulation is run multiple times over a range of values for the independent variable. To visualize the change in behavior as the independent variable is perturbed, the value of each metric can be plotted as a function of the parameter or initial concentration value, or alternatively, the time courses for each trial can be plotted on the same axes.
2.2.3 Deterministic simulations
One type of sensitivity analysis supported by Dynetica uses deterministic simulations. An example analysis, using the ePop system [
17], is shown in Figure 3A. Since the simulations are deterministic, there only needs to be one for each combination of parameter values.
The time courses produced by a deterministic sensitivity analysis can be programmatically analyzed to produce a more succinct plot. In Figure 3B, the cell count (n) time courses shown in Figure 3A are analyzed for their peak frequencies, which is plotted as a function of the parameter w.
2.2.4 Stochastic simulations
The sensitivity analysis tools may be used with stochastic algorithms as well. However, the results may prove less meaningful because the change in system behavior may be due to random noise rather than the deviations in the independent variable. To examine how stochastic systems respond to changes in parameters or initial concentrations, Dynetica 2.0 provides a tool for performing repeated stochastic simulations for each value of the independent variable, as demonstrated with the Rb-E2F system analysis displayed in Figure 3C. Both the time courses and the distribution of values can be viewed for user-specified time intervals, as shown in Figure 3D.
Dynetica 2.0 also provides the ability to perform sensitivity analyses using stochastic algorithms and generate distributions of values, such as the Rb-E2F system analysis seen in Figure 3C.
2.2.5 Summary of Sensitivity Analysis Types
In summary, the following types of sensitivity analysis are supported in Dynetica 2.0:
• Deterministic
o Basic (time courses);
o Processed according to one of the metrics in the section of metrics.
• Stochastic (repeated to obtain a distribution)
2.3 Parameter Search
Dynetica 2.0 provides tools for searching the parameter space of a system in order to optimize user-defined objectives. This is useful for when the user knows the desired behavior of a system and wants to find a set of parameters for which the system closely matches this behavior. The search may be performed over all or any subset of the system’s parameters.
The user defines the desired behavior of the system by specifying a number of target objectives, which the algorithm will attempt to match. The user selects a metric and specifies a goal for that metric’s value, which can be to minimize, maximize, or target its value to a specific number. Target objectives can also be given different weights, with higher weighted objectives having a greater influence on the search than lower ones.
The parameter search is done using a genetic algorithm, a computational technique that imitates natural selection to heuristically find an optimal solution given a set of constraints. The initial population is made up of individual parameter vectors. In order to calculate the evolutionary fitness of each individual parameter vector, its parameters are applied to the system and the simulation is run. The vector then receives a fitness score between 0 and 100 based on how well the simulation matches the user-defined metrics. At each generation, high scoring parameter vectors survive and reproduce, while low scoring vectors die off, maintaining a constant population size. The details of the implementation of this algorithm can be found in the supplementary materials.
To run a parameter search, the user must specify which parameters will be perturbed. All others will remain fixed at their user-defined values. The user then must specify the target objectives for the algorithm to optimize, each of which consists of a metric, a goal for that metric, and a relative weight. The user must also specify the simulation duration, population size, and maximum number of generations for which the algorithm will run. The user can additionally define a stopping threshold, which is a fraction specifying what score is needed to cause the search to complete. For example, the default value of 0.05 indicates that a fitness score of 95 or higher (out of 100) is sufficient to stop the search. Setting this value to 0 means the search will not stop unless the target objectives are matched exactly or the maximum number of generations is reached.
By default, each parameter vector is created by random normal perturbations around the user’s defined value. This allows the algorithm to use the user-given values as an initial guess. However, the user may also chose to create the initial population by selecting uniformly random numbers between the parameter’s specified minimum and maximum values. This may produce better results by exploring new areas of the parameter space. Note, however, that if this option is selected, it is important to specify realistic constraints for the parameter’s minimum and maximum values. Figure 4 shows a typical example of using the parameter search tool to identify a new parameter set for a specific target function.
Dynetica 2.0 also features a tool for searching the parameter space of a system when the desired behavior is not a feature of a single simulation, but of the dose response or sensitivity analysis curves. In addition to the inputs above, the user must specify a parameter or initial substance concentration to use as the independent variable, as well as a range of values to scan over. The same genetic search algorithm is used; however, rather than using a single metric value, the fitness score is calculated as a function of the shape of curves when the specified metrics are plotted as against the independent variable.
For each metric selected, there are six target functions that can be used to calculate the fitness score. The range of the sensitivity plot can be minimized or maximized, resulting in the metric being minimally or maximally responsive to perturbations in the independent variable, respectively. The user may also select to maximize or minimize the difference between the first two local extrema, or the first local extremum and the steady state value. The user may also select to maximize the difference from the initial value to the first local extremum. Lastly, the user may choose for the response curve to be maximally linear, in which case a linear regression will be case the score will be based on the coefficient of determination. Figure 5 illustrates the use of the search algorithm to identify a parameter set to generate a strong biphasic response in a synthetic gene circuit.
Further details on Dynetica 2.0’s Parameter Search feature can be found in section 2 of the Supplemental Materials.
3 DISCUSSION
Significant interface and experience improvements have been made to Dynetica, in the form of the expression variable/module features, as well as import/export functionality. These enhancements ease the process of constructing and evaluating models in Dynetica 2.0, as well as the collaboration process when working with researchers who utilize other pieces of dynamic network simulation software.
Dynetica 2.0 has added support for much more sophisticated analyses of stochastic simulations. Stochastic simulations can be repeated with a fixed parameter set to generate outcome distributions, or with a varying parameter set to accomplish sensitivity analyses. For small systems, stochastic analyses can provide realistic insight into laboratory phenomena. For those systems that are too large for stochastic algorithms to be practical, sensitivity analyses can also be conducted using deterministic algorithms. In addition, a genetic search algorithm can be used to determine a set of parameter values that would produce the desired value of a specified metric.
Some of the limitations of Dynetica 2.0 are its inability to directly model spatial gradients and multiple compartments. However, many basic modeling tasks do not require the added complexity of multiple compartments, computationally intensive simulations, and spatial modeling. Despite these limitations, we believe that Dynetica 2.0 serves its purpose by providing a robust and easy to use modeling and simulation environment for non-programmers.
4 ACCESSING THE SOFTWARE
Dynetica is freely available to the research community.
(i) General users can access the binary code at: http://www.genome.duke.edu/labs/YouLab/software/dynetica/index.php.
(ii) Developers can access the source code at:
https://github.com/youlab/dynetica.
Higher Education Press and Springer-Verlag Berlin Heidelberg
PDF
(1172KB)
