SBML Level 3: an extensible format for the exchange and reuse of biological models.

Systems biology has experienced dramatic growth in the number, size, and complexity of computational models. To reproduce simulation results and reuse models, researchers must exchange unambiguous model descriptions. We review the latest edition of the Systems Biology Markup Language (SBML), a format designed for this purpose. A community of modelers and software authors developed SBML Level 3 over the past decade. Its modular form consists of a core suited to representing reaction-based models and packages that extend the core with features suited to other model types including constraint-based models, reaction-diffusion models, logical network models, and rule-based models. The format leverages two decades of SBML and a rich software ecosystem that transformed how systems biologists build and interact with models. More recently, the rise of multiscale models of whole cells and organs, and new data sources such as single-cell measurements and live imaging, has precipitated new ways of integrating data with models. We provide our perspectives on the challenges presented by these developments and how SBML Level 3 provides the foundation needed to support this evolution.

Analysis of cell division patterns in the Arabidopsis shoot apical meristem.

The stereotypic pattern of cell shapes in the Arabidopsis shoot apical meristem (SAM) suggests that strict rules govern the placement of new walls during cell division. When a cell in the SAM divides, a new wall is built that connects existing walls and divides the cytoplasm of the daughter cells. Because features that are determined by the placement of new walls such as cell size, shape, and number of neighbors are highly regular, rules must exist for maintaining such order. Here we present a quantitative model of these rules that incorporates different observed features of cell division. Each feature is incorporated into a "potential function" that contributes a single term to a total analog of potential energy. New cell walls are predicted to occur at locations where the potential function is minimized. Quantitative terms that represent the well-known historical rules of plant cell division, such as those given by Hofmeister, Errera, and Sachs are developed and evaluated against observed cell divisions in the epidermal layer (L1) of Arabidopsis thaliana SAM. The method is general enough to allow additional terms for nongeometric properties such as internal concentration gradients and mechanical tensile forces.

Pycellerator: an arrow-based reaction-like modelling language for biological simulations.

MOTIVATION: We introduce Pycellerator, a Python library for reading Cellerator arrow notation from standard text files, conversion to differential equations, generating stand-alone Python solvers, and optionally running and plotting the solutions. All of the original Cellerator arrows, which represent reactions ranging from mass action, Michales-Menten-Henri (MMH) and Gene-Regulation (GRN) to Monod-Wyman-Changeaux (MWC), user defined reactions and enzymatic expansions (KMech), were previously represented with the Mathematica extended character set. These are now typed as reaction-like commands in ASCII text files that are read by Pycellerator, which includes a Python command line interface (CLI), a Python application programming interface (API) and an iPython notebook interface. RESULTS: Cellerator reaction arrows are now input in text files. The arrows are parsed by Pycellerator and translated into differential equations in Python, and Python code is automatically generated to solve the system. Time courses are produced by executing the auto-generated Python code. Users have full freedom to modify the solver and utilize the complete set of standard Python tools. The new libraries are completely independent of the old Cellerator software and do not require Mathematica. AVAILABILITY AND IMPLEMENTATION: All software is available (GPL) from the github repository at https://github.com/biomathman/pycellerator/releases. Details, including installation instructions and a glossary of acronyms and terms, are given in the Supplementary information.

A Pycellerator Tutorial.

We present a tutorial on using Pycellerator for biomolecular simulations. Models are described in human readable (and editable) text files (UTF8 or ASCII) containing collections of reactions, assignments, initial conditions, function definitions, and rate constants. These models are then converted into a Python program that can optionally solve the system, e.g., as a system of differential equations using ODEINT, or be run by another program. The input language implements an extended version of the Cellerator arrow notation, including mass action, Hill functions, S-Systems, MWC, and reactions with user-defined kinetic laws. Simple flux balance analysis is also implemented. We will demonstrate the implementation and analysis of progressively more complex models, starting from simple mass action through indexed cascades. Pycellerator can be used as a library that is integrated into other programs, run as a command line program, or in iPython notebooks. It is implemented in Python 2.7 and available under an open source GPL license.

Application of a generalized MWC model for the mathematical simulation of metabolic pathways regulated by allosteric enzymes.

In our effort to elucidate the systems biology of the model organism, Escherichia coli, we have developed a mathematical model that simulates the allosteric regulation for threonine biosynthesis pathway starting from aspartate. To achieve this goal, we used kMech, a Cellerator language extension that describes enzyme mechanisms for the mathematical modeling of metabolic pathways. These mechanisms are converted by Cellerator into ordinary differential equations (ODEs) solvable by Mathematica. In this paper, we describe a more flexible model in Cellerator, which generalizes the Monod, Wyman, Changeux (MWC) model for enzyme allosteric regulation to allow for multiple substrate, activator and inhibitor binding sites. Furthermore, we have developed a model that describes the behavior of the bifunctional allosteric enzyme aspartate kinase I-homoserine dehydrogenase I (AKI-HDHI). This model predicts the partition of enzyme activities in the steady state which paves the way for a more generalized prediction of the behavior of bifunctional enzymes.

Reproducibility, fidelity, and discriminant validity of mRNA amplification for microarray analysis from primary hematopoietic cells.

Analysis of gene expression in clinical samples poses special challenges, including limited RNA availability and poor RNA quality. Quantitative information regarding reliability of RNA amplification methodologies applied to primary cells and representativeness of resulting gene expression profiles is limited. We evaluated four protocols for RNA amplification from peripheral blood mononuclear cells. Results obtained with 100 ng or 10 ng of RNA amplified using two rounds of cDNA synthesis and in vitro transcription were compared with control 2.5-microg RNA samples processed using a single round of in vitro transcription. Samples were hybridized to Affymetrix HG-U133A arrays. Considerable differences in results were obtained with different protocols. The optimal protocol resulted in highly reproducible gene expression profiles from amplified samples (r = 0.98) and good correlation between amplified and control samples (r = 0.94). Using the optimal protocol dissimilarities of gene expression between mononuclear cells from a normal individual and a patient with myelodysplastic syndrome were primarily maintained after amplification compared with controls. We conclude that small variations in methodology introduce considerable distortion of gene expression profiles obtained after RNA amplification from clinical samples and too strong a focus on a very small number of genes picked from an array analysis could be unduly influenced by seemingly acceptable methodologies. However, it is possible to obtain reproducible and representative results using optimized protocols.

Using cellzilla for plant growth simulations at the cellular level.

Cellzilla is a two-dimensional tissue simulation platform for plant modeling utilizing Cellerator arrows. Cellerator describes biochemical interactions with a simplified arrow-based notation; all interactions are input as reactions and are automatically translated to the appropriate differential equations using a computer algebra system. Cells are represented by a polygonal mesh of well-mixed compartments. Cell constituents can interact intercellularly via Cellerator reactions utilizing diffusion, transport, and action at a distance, as well as amongst themselves within a cell. The mesh data structure consists of vertices, edges (vertex pairs), and cells (and optional intercellular wall compartments) as ordered collections of edges. Simulations may be either static, in which cell constituents change with time but cell size and shape remain fixed; or dynamic, where cells can also grow. Growth is controlled by Hookean springs associated with each mesh edge and an outward pointing pressure force. Spring rest length grows at a rate proportional to the extension beyond equilibrium. Cell division occurs when a specified constituent (or cell mass) passes a (random, normally distributed) threshold. The orientation of new cell walls is determined either by Errera's rule, or by a potential model that weighs contributions due to equalizing daughter areas, minimizing wall length, alignment perpendicular to cell extension, and alignment perpendicular to actual growth direction.

An auxin-driven polarized transport model for phyllotaxis.

Recent studies show that plant organ positioning may be mediated by localized concentrations of the plant hormone auxin. Auxin patterning in the shoot apical meristem is in turn brought about by the subcellular polar distribution of the putative auxin efflux mediator, PIN1. However, the question of what signals determine PIN1 polarization and how this gives rise to regular patterns of auxin concentration remains unknown. Here we address these questions by using mathematical modeling combined with confocal imaging. We propose a model that is based on the assumption that auxin influences the polarization of its own efflux within the meristem epidermis. We show that such a model is sufficient to create regular spatial patterns of auxin concentration on systems with static and dynamic cellular connectivities, the latter governed by a mechanical model. We also optimize parameter values for the PIN1 dynamics by using a detailed auxin transport model, for which parameter values are taken from experimental estimates, together with a template consisting of cell and wall compartments as well as PIN1 concentrations quantitatively extracted from confocal data. The model shows how polarized transport can drive the formation of regular patterns.

Application of a generalized MWC model for the mathematical simulation of metabolic pathways regulated by allosteric enzymes.

In our effort to elucidate the systems biology of the model organism, Escherichia coli, we have developed a mathematical model that simulates the allosteric regulation for threonine biosynthesis pathway starting from aspartate. To achieve this goal, we used kMech, a Cellerator language extension that describes enzyme mechanisms for the mathematical modeling of metabolic pathways. These mechanisms are converted by Cellerator into ordinary differential equations (ODEs) solvable by Mathematica. In this paper, we describe a more flexible model in Cellerator, which generalizes the Monod, Wyman, Changeux (MWC) model for enzyme allosteric regulation to allow for multiple substrate, activator and inhibitor binding sites. Furthermore, we have developed a model that describes the behavior of the bifunctional allosteric enzyme aspartate Kinase I-Homoserine Dehydrogenase I (AKI-HDHI). This model predicts the partition of enzyme activities in the steady state which paves a way for a more generalized prediction of the behavior of bifunctional enzymes.

An enzyme mechanism language for the mathematical modeling of metabolic pathways.

MOTIVATION: As a first step toward the elucidation of the systems biology of complex biological systems, it was our goal to mathematically model common enzyme catalytic and regulatory mechanisms that repeatedly appear in biological processes such as signal transduction and metabolic pathways. RESULTS: We describe kMech, a Cellerator language extension that describes a suite of enzyme mechanisms. Each enzyme mechanism is parsed by kMech into a set of fundamental association-dissociation reactions that are translated by Cellerator into ordinary differential equations that are numerically solved by Mathematica. In addition, we present methods that use commonly available kinetic measurements to estimate rate constants required to solve these differential equations.

A mathematical model for the branched chain amino acid biosynthetic pathways of Escherichia coli K12.

As a first step toward the elucidation of the systems biology of the model organism Escherichia coli, it was our goal to mathematically model a metabolic system of intermediate complexity, namely the well studied end product-regulated pathways for the biosynthesis of the branched chain amino acids L-isoleucine, L-valine, and L-leucine. This has been accomplished with the use of kMech (Yang, C.-R., Shapiro, B. E., Mjolsness, E. D., and Hatfield, G. W. (2005) Bioinformatics 21, in press), a Cellerator (Shapiro, B. E., Levchenko, A., Meyerowitz, E. M., Wold, B. J., and Mjolsness, E. D. (2003) Bioinformatics 19, 677-678) language extension that describes a suite of enzyme reaction mechanisms. Each enzyme mechanism is parsed by kMech into a set of fundamental association-dissociation reactions that are translated by Cellerator into ordinary differential equations. These ordinary differential equations are numerically solved by Mathematica. Any metabolic pathway can be simulated by stringing together appropriate kMech models and providing the physical and kinetic parameters for each enzyme in the pathway. Writing differential equations is not required. The mathematical model of branched chain amino acid biosynthesis in E. coli K12 presented here incorporates all of the forward and reverse enzyme reactions and regulatory circuits of the branched chain amino acid biosynthetic pathways, including single and multiple substrate (Ping Pong and Bi Bi) enzyme kinetic reactions, feedback inhibition (allosteric, competitive, and non-competitive) mechanisms, the channeling of metabolic flow through isozymes, the channeling of metabolic flow via transamination reactions, and active transport mechanisms. This model simulates the results of experimental measurements.

Modeling the organization of the WUSCHEL expression domain in the shoot apical meristem.

MOTIVATION: The above-ground tissues of higher plants are generated from a small region of cells situated at the plant apex called the shoot apical meristem. An important genetic control circuit modulating the size of the Arabidopsis thaliana meristem is a feed-back network between the CLAVATA3 and WUSCHEL genes. Although the expression patterns for these genes do not overlap, WUSCHEL activity is both necessary and sufficient (when expressed ectopically) for the induction of CLAVATA3 expression. However, upregulation of CLAVATA3 in conjunction with the receptor kinase CLAVATA1 results in the downregulation of WUSCHEL. Despite much work, experimental data for this network are incomplete and additional hypotheses are needed to explain the spatial locations and dynamics of these expression domains. Predictive mathematical models describing the system should provide a useful tool for investigating and discriminating among possible hypotheses, by determining which hypotheses best explain observed gene expression dynamics. RESULTS: We are developing a method using in vivo live confocal microscopy to capture quantitative gene expression data and create templates for computational models. We present two models accounting for the organization of the WUSCHEL expression domain. Our preferred model uses a reaction-diffusion mechanism in which an activator induces WUSCHEL expression. This model is able to organize the WUSCHEL expression domain. In addition, the model predicts the dynamical reorganization seen in experiments where cells, including the WUSCHEL domain, are ablated, and it also predicts the spatial expansion of the WUSCHEL domain resulting from removal of the CLAVATA3 signal. AVAILABILITY: An extended description of the model framework and image processing algorithms can be found at http://www.computableplant.org, together with additional results and simulation movies. SUPPLEMENTARY INFORMATION: http://www.computableplant.org/ and alternatively for a direct link to the page, http://computableplant.ics.uci.edu/bti1036 can be accessed.

MathSBML: a package for manipulating SBML-based biological models.

UNLABELLED: MathSBML is a Mathematica package designed for manipulating Systems Biology Markup Language (SBML) models. It converts SBML models into Mathematica data structures and provides a platform for manipulating and evaluating these models. Once a model is read by MathSBML, it is fully compatible with standard Mathematica functions such as NDSolve (a differential-algebraic equations solver). MathSBML also provides an application programming interface for viewing, manipulating, running numerical simulations; exporting SBML models; and converting SBML models in to other formats, such as XPP, HTML and FORTRAN. By accessing the full breadth of Mathematica functionality, MathSBML is fully extensible to SBML models of any size or complexity. AVAILABILITY: Open Source (LGPL) at http://www.sbml.org and http://www.sf.net/projects/sbml

Cellerator: extending a computer algebra system to include biochemical arrows for signal transduction simulations.

Cellerator describes single and multi-cellular signal transduction networks (STN) with a compact, optionally palette-driven, arrow-based notation to represent biochemical reactions and transcriptional activation. Multi-compartment systems are represented as graphs with STNs embedded in each node. Interactions include mass-action, enzymatic, allosteric and connectionist models. Reactions are translated into differential equations and can be solved numerically to generate predictive time courses or output as systems of equations that can be read by other programs. Cellerator simulations are fully extensible and portable to any operating system that supports Mathematica, and can be indefinitely nested within larger data structures to produce highly scaleable models.