Efficient parameter estimation for spatio-temporal models of pattern formation: case study of Drosophila melanogaster

ISRES+: an improved evolutionary strategy for function minimization to estimate the free parameters of systems biology models

MOTIVATION

Mathematical models in systems biology help generate hypotheses, guide experimental design, and infer the dynamics of gene regulatory networks. These models are characterized by phenomenological or mechanistic parameters, which are typically hard to measure. Therefore, efficient parameter estimation is central to model development. Global optimization techniques, such as evolutionary algorithms (EAs), are applied to estimate model parameters by inverse modeling, i.e. calibrating models by minimizing a function that evaluates a measure of the error between model predictions and experimental data. EAs estimate model parameters "fittest individuals" by generating a large population of individuals using strategies like recombination and mutation over multiple "generations." Typically, only a few individuals from each generation are used to create new individuals in the next generation. Improved Evolutionary Strategy by Stochastic Ranking (ISRES), proposed by Runnarson and Yao, is one such EA that is widely used in systems biology to estimate parameters. ISRES uses information at most from a pair of individuals in any generation to create a new population to minimize the error. In this article, we propose an efficient evolutionary strategy, ISRES+, which builds on ISRES by combining information from all individuals across the population and across all generations to develop a better understanding of the fitness landscape.

RESULTS

ISRES+ uses the additional information generated by the algorithm during evolution to approximate the local neighborhood around the best-fit individual using linear least squares fits in one and two dimensions, enabling efficient parameter estimation. ISRES+ outperforms ISRES and results in fitter individuals with a tighter distribution over multiple runs, such that a typical run of ISRES+ estimates parameters with a higher goodness-of-fit compared with ISRES.

AVAILABILITY AND IMPLEMENTATION

Algorithm and implementation: Github-https://github.com/gtreeves/isres-plus-bandodkar-2022.

Research on Multi-Time-Delay Gene Regulation Network Based on Fuzzy Label Propagation

In view of time delay existing in gene regulation, by using the analysis idea and methods of complex network, this paper proposes a multi-time-delay gene regulation network analysis method based on the fuzzy label propagation. The algorithm takes the relative change trend coefficient, the correlation coefficient, and the mutual information as the similarity measurement indexes for the gene pair, fully reflecting the correlation of gene pairs and simultaneously obtaining the gene regulation relationship and the time delay through the fuzzy label propagation algorithm of the semisupervised learning. Experimental results of the cell cycle-regulated genes of yeast show that the proposed construction method of GRN can not only correctly select potential regulation genes but also provide details about the gene regulator model, thereby more accurately constructing gene regulation network.

Global optimization using Gaussian processes to estimate biological parameters from image data

Parameter estimation is a major challenge in computational modeling of biological processes. This is especially the case in image-based modeling where the inherently quantitative output of the model is measured against image data, which is typically noisy and non-quantitative. In addition, these models can have a high computational cost, limiting the number of feasible simulations, and therefore rendering most traditional parameter estimation methods unsuitable. In this paper, we present a pipeline that uses Gaussian process learning to estimate biological parameters from noisy, non-quantitative image data when the model has a high computational cost. This approach is first successfully tested on a parametric function with the goal of retrieving the original parameters. We then apply it to estimating parameters in a biological setting by fitting artificial in-situ hybridization (ISH) data of the developing murine limb bud. We expect that this method will be of use in a variety of modeling scenarios where quantitative data is missing and the use of standard parameter estimation approaches in biological modeling is prohibited by the computational cost of the model.

A time series driven decomposed evolutionary optimization approach for reconstructing large-scale gene regulatory networks based on fuzzy cognitive maps

BACKGROUND

Reconstructing gene regulatory networks (GRNs) from expression data plays an important role in understanding the fundamental cellular processes and revealing the underlying relations among genes. Although many algorithms have been proposed to reconstruct GRNs, more rapid and efficient methods which can handle large-scale problems still need to be developed. The process of reconstructing GRNs can be formulated as an optimization problem, which is actually reconstructing GRNs from time series data, and the reconstructed GRNs have good ability to simulate the observed time series. This is a typical big optimization problem, since the number of variables needs to be optimized increases quadratically with the scale of GRNs, resulting an exponential increase in the number of candidate solutions. Thus, there is a legitimate need to devise methods capable of automatically reconstructing large-scale GRNs.

RESULTS

In this paper, we use fuzzy cognitive maps (FCMs) to model GRNs, in which each node of FCMs represent a single gene. However, most of the current training algorithms for FCMs are only able to train FCMs with dozens of nodes. Here, a new evolutionary algorithm is proposed to train FCMs, which combines a dynamical multi-agent genetic algorithm (dMAGA) with the decomposition-based model, and termed as dMAGA-FCM_D, which is able to deal with large-scale FCMs with up to 500 nodes. Both large-scale synthetic FCMs and the benchmark DREAM4 for reconstructing biological GRNs are used in the experiments to validate the performance of dMAGA-FCM_D.

CONCLUSIONS

The dMAGA-FCM_D is compared with the other four algorithms which are all state-of-the-art FCM training algorithms, and the results show that the dMAGA-FCM_D performs the best. In addition, the experimental results on FCMs with 500 nodes and DREAM4 project demonstrate that dMAGA-FCM_D is capable of effectively and computationally efficiently training large-scale FCMs and GRNs.

Automated numerical simulation of biological pattern formation based on visual feedback simulation framework

There are various fantastic biological phenomena in biological pattern formation. Mathematical modeling using reaction-diffusion partial differential equation systems is employed to study the mechanism of pattern formation. However, model parameter selection is both difficult and time consuming. In this paper, a visual feedback simulation framework is proposed to calculate the parameters of a mathematical model automatically based on the basic principle of feedback control. In the simulation framework, the simulation results are visualized, and the image features are extracted as the system feedback. Then, the unknown model parameters are obtained by comparing the image features of the simulation image and the target biological pattern. Considering two typical applications, the visual feedback simulation framework is applied to fulfill pattern formation simulations for vascular mesenchymal cells and lung development. In the simulation framework, the spot, stripe, labyrinthine patterns of vascular mesenchymal cells, the normal branching pattern and the branching pattern lacking side branching for lung branching are obtained in a finite number of iterations. The simulation results indicate that it is easy to achieve the simulation targets, especially when the simulation patterns are sensitive to the model parameters. Moreover, this simulation framework can expand to other types of biological pattern formation.

Parameter estimation in large-scale systems biology models: a parallel and self-adaptive cooperative strategy

Background

The development of large-scale kinetic models is one of the current key issues in computational systems biology and bioinformatics. Here we consider the problem of parameter estimation in nonlinear dynamic models. Global optimization methods can be used to solve this type of problems but the associated computational cost is very large. Moreover, many of these methods need the tuning of a number of adjustable search parameters, requiring a number of initial exploratory runs and therefore further increasing the computation times.

Here we present a novel parallel method, self-adaptive cooperative enhanced scatter search (saCeSS), to accelerate the solution of this class of problems. The method is based on the scatter search optimization metaheuristic and incorporates several key new mechanisms: (i) asynchronous cooperation between parallel processes, (ii) coarse and fine-grained parallelism, and (iii) self-tuning strategies.

Results

The performance and robustness of saCeSS is illustrated by solving a set of challenging parameter estimation problems, including medium and large-scale kinetic models of the bacterium E. coli, bakerés yeast S. cerevisiae, the vinegar fly D. melanogaster, Chinese Hamster Ovary cells, and a generic signal transduction network.

The results consistently show that saCeSS is a robust and efficient method, allowing very significant reduction of computation times with respect to several previous state of the art methods (from days to minutes, in several cases) even when only a small number of processors is used.

Conclusions

The new parallel cooperative method presented here allows the solution of medium and large scale parameter estimation problems in reasonable computation times and with small hardware requirements. Further, the method includes self-tuning mechanisms which facilitate its use by non-experts. We believe that this new method can play a key role in the development of large-scale and even whole-cell dynamic models.

Electronic supplementary material

The online version of this article (doi:10.1186/s12859-016-1452-4) contains supplementary material, which is available to authorized users.

MoCha: Molecular Characterization of Unknown Pathways

Automated methods for the reverse-engineering of complex regulatory networks are paving the way for the inference of mechanistic comprehensive models directly from experimental data. These novel methods can infer not only the relations and parameters of the known molecules defined in their input datasets, but also unknown components and pathways identified as necessary by the automated algorithms. Identifying the molecular nature of these unknown components is a crucial step for making testable predictions and experimentally validating the models, yet no specific and efficient tools exist to aid in this process. To this end, we present here MoCha (Molecular Characterization), a tool optimized for the search of unknown proteins and their pathways from a given set of known interacting proteins. MoCha uses the comprehensive dataset of protein-protein interactions provided by the STRING database, which currently includes more than a billion interactions from over 2,000 organisms. MoCha is highly optimized, performing typical searches within seconds. We demonstrate the use of MoCha with the characterization of unknown components from reverse-engineered models from the literature. MoCha is useful for working on network models by hand or as a downstream step of a model inference engine workflow and represents a valuable and efficient tool for the characterization of unknown pathways using known data from thousands of organisms. MoCha and its source code are freely available online under the GPLv3 license.

Inferring regulatory networks from experimental morphological phenotypes: a computational method reverse-engineers planarian regeneration

Transformative applications in biomedicine require the discovery of complex regulatory networks that explain the development and regeneration of anatomical structures, and reveal what external signals will trigger desired changes of large-scale pattern. Despite recent advances in bioinformatics, extracting mechanistic pathway models from experimental morphological data is a key open challenge that has resisted automation. The fundamental difficulty of manually predicting emergent behavior of even simple networks has limited the models invented by human scientists to pathway diagrams that show necessary subunit interactions but do not reveal the dynamics that are sufficient for complex, self-regulating pattern to emerge. To finally bridge the gap between high-resolution genetic data and the ability to understand and control patterning, it is critical to develop computational tools to efficiently extract regulatory pathways from the resultant experimental shape phenotypes. For example, planarian regeneration has been studied for over a century, but despite increasing insight into the pathways that control its stem cells, no constructive, mechanistic model has yet been found by human scientists that explains more than one or two key features of its remarkable ability to regenerate its correct anatomical pattern after drastic perturbations. We present a method to infer the molecular products, topology, and spatial and temporal non-linear dynamics of regulatory networks recapitulating in silico the rich dataset of morphological phenotypes resulting from genetic, surgical, and pharmacological experiments. We demonstrated our approach by inferring complete regulatory networks explaining the outcomes of the main functional regeneration experiments in the planarian literature; By analyzing all the datasets together, our system inferred the first systems-biology comprehensive dynamical model explaining patterning in planarian regeneration. This method provides an automated, highly generalizable framework for identifying the underlying control mechanisms responsible for the dynamic regulation of growth and form.

Developmental and regenerative biology experiments are producing a huge number of morphological phenotypes from functional perturbation experiments. However, existing pathway models do not generally explain the dynamic regulation of anatomical shape due to the difficulty of inferring and testing non-linear regulatory networks responsible for appropriate form, shape, and pattern. We present a method that automates the discovery and testing of regulatory networks explaining morphological outcomes directly from the resultant phenotypes, producing network models as testable hypotheses explaining regeneration data. Our system integrates a formalization of the published results in planarian regeneration, an in silico simulator in which the patterning properties of regulatory networks can be quantitatively tested in a regeneration assay, and a machine learning module that evolves networks whose behavior in this assay optimally matches the database of planarian results. We applied our method to explain the key experiments in planarian regeneration, and discovered the first comprehensive model of anterior-posterior patterning in planaria under surgical, pharmacological, and genetic manipulations. Beyond the planarian data, our approach is readily generalizable to facilitate the discovery of testable regulatory networks in developmental biology and biomedicine, and represents the first developmental model discovered de novo from morphological outcomes by an automated system.

The presence of nuclear cactus in the early Drosophila embryo may extend the dynamic range of the dorsal gradient

In a developing embryo, the spatial distribution of a signaling molecule, or a morphogen gradient, has been hypothesized to carry positional information to pattern tissues. Recent measurements of morphogen distribution have allowed us to subject this hypothesis to rigorous physical testing. In the early Drosophila embryo, measurements of the morphogen Dorsal, which is a transcription factor responsible for initiating the earliest zygotic patterns along the dorsal-ventral axis, have revealed a gradient that is too narrow to pattern the entire axis. In this study, we use a mathematical model of Dorsal dynamics, fit to experimental data, to determine the ability of the Dorsal gradient to regulate gene expression across the entire dorsal-ventral axis. We found that two assumptions are required for the model to match experimental data in both Dorsal distribution and gene expression patterns. First, we assume that Cactus, an inhibitor that binds to Dorsal and prevents it from entering the nuclei, must itself be present in the nuclei. And second, we assume that fluorescence measurements of Dorsal reflect both free Dorsal and Cactus-bound Dorsal. Our model explains the dynamic behavior of the Dorsal gradient at lateral and dorsal positions of the embryo, the ability of Dorsal to regulate gene expression across the entire dorsal-ventral axis, and the robustness of gene expression to stochastic effects. Our results have a general implication for interpreting fluorescence-based measurements of signaling molecules.

Spatial gradients of protein-level time delays set the pace of the traveling segmentation clock waves

The vertebrate segmentation clock is a gene expression oscillator controlling rhythmic segmentation of the vertebral column during embryonic development. The period of oscillations becomes longer as cells are displaced along the posterior to anterior axis, which results in traveling waves of clock gene expression sweeping in the unsegmented tissue. Although various hypotheses necessitating the inclusion of additional regulatory genes into the core clock network at different spatial locations have been proposed, the mechanism underlying traveling waves has remained elusive. Here, we combined molecular-level computational modeling and quantitative experimentation to solve this puzzle. Our model predicts the existence of an increasing gradient of gene expression time delays along the posterior to anterior direction to recapitulate spatiotemporal profiles of the traveling segmentation clock waves in different genetic backgrounds in zebrafish. We validated this prediction by measuring an increased time delay of oscillatory Her1 protein production along the unsegmented tissue. Our results refuted the need for spatial expansion of the core feedback loop to explain the occurrence of traveling waves. Spatial regulation of gene expression time delays is a novel way of creating dynamic patterns; this is the first report demonstrating such a control mechanism in any tissue and future investigations will explore the presence of analogous examples in other biological systems.

Making models match measurements: model optimization for morphogen patterning networks

Mathematical modeling of developmental signaling networks has played an increasingly important role in the identification of regulatory mechanisms by providing a sandbox for hypothesis testing and experiment design. Whether these models consist of an equation with a few parameters or dozens of equations with hundreds of parameters, a prerequisite to model-based discovery is to bring simulated behavior into agreement with observed data via parameter estimation. These parameters provide insight into the system (e.g., enzymatic rate constants describe enzyme properties). Depending on the nature of the model fit desired - from qualitative (relative spatial positions of phosphorylation) to quantitative (exact agreement of spatial position and concentration of gene products) - different measures of data-model mismatch are used to estimate different parameter values, which contain different levels of usable information and/or uncertainty. To facilitate the adoption of modeling as a tool for discovery alongside other tools such as genetics, immunostaining, and biochemistry, careful consideration needs to be given to how well a model fits the available data, what the optimized parameter values mean in a biological context, and how the uncertainty in model parameters and predictions plays into experiment design. The core discussion herein pertains to the quantification of model-to-data agreement, which constitutes the first measure of a model's performance and future utility to the problem at hand. Integration of this experimental data and the appropriate choice of objective measures of data-model agreement will continue to drive modeling forward as a tool that contributes to experimental discovery. The Drosophila melanogaster gap gene system, in which model parameters are optimized against in situ immunofluorescence intensities, demonstrates the importance of error quantification, which is applicable to a wide array of developmental modeling studies.

An extended gene protein/products Boolean network model including post-transcriptional regulation

Background

Networks Biology allows the study of complex interactions between biological systems using formal, well structured, and computationally friendly models. Several different network models can be created, depending on the type of interactions that need to be investigated. Gene Regulatory Networks (GRN) are an effective model commonly used to study the complex regulatory mechanisms of a cell. Unfortunately, given their intrinsic complexity and non discrete nature, the computational study of realistic-sized complex GRNs requires some abstractions. Boolean Networks (BNs), for example, are a reliable model that can be used to represent networks where the possible state of a node is a boolean value (0 or 1). Despite this strong simplification, BNs have been used to study both structural and dynamic properties of real as well as randomly generated GRNs.

Results

In this paper we show how it is possible to include the post-transcriptional regulation mechanism (a key process mediated by small non-coding RNA molecules like the miRNAs) into the BN model of a GRN. The enhanced BN model is implemented in a software toolkit (EBNT) that allows to analyze boolean GRNs from both a structural and a dynamic point of view. The open-source toolkit is compatible with available visualization tools like Cytoscape and allows to run detailed analysis of the network topology as well as of its attractors, trajectories, and state-space. In the paper, a small GRN built around the mTOR gene is used to demonstrate the main capabilities of the toolkit.

Conclusions

The extended model proposed in this paper opens new opportunities in the study of gene regulation. Several of the successful researches done with the support of BN to understand high-level characteristics of regulatory networks, can now be improved to better understand the role of post-transcriptional regulation for example as a network-wide noise-reduction or stabilization mechanisms.

Using evolutionary computations to understand the design and evolution of gene and cell regulatory networks

This paper surveys modeling approaches for studying the evolution of gene regulatory networks (GRNs). Modeling of the design or 'wiring' of GRNs has become increasingly common in developmental and medical biology, as a means of quantifying gene-gene interactions, the response to perturbations, and the overall dynamic motifs of networks. Drawing from developments in GRN 'design' modeling, a number of groups are now using simulations to study how GRNs evolve, both for comparative genomics and to uncover general principles of evolutionary processes. Such work can generally be termed evolution in silico. Complementary to these biologically-focused approaches, a now well-established field of computer science is Evolutionary Computations (ECs), in which highly efficient optimization techniques are inspired from evolutionary principles. In surveying biological simulation approaches, we discuss the considerations that must be taken with respect to: (a) the precision and completeness of the data (e.g. are the simulations for very close matches to anatomical data, or are they for more general exploration of evolutionary principles); (b) the level of detail to model (we proceed from 'coarse-grained' evolution of simple gene-gene interactions to 'fine-grained' evolution at the DNA sequence level); (c) to what degree is it important to include the genome's cellular context; and (d) the efficiency of computation. With respect to the latter, we argue that developments in computer science EC offer the means to perform more complete simulation searches, and will lead to more comprehensive biological predictions.

Systematic parameter estimation in data-rich environments for cell signalling dynamics

Motivation: Computational models of biological signalling networks, based on ordinary differential equations (ODEs), have generated many insights into cellular dynamics, but the model-building process typically requires estimating rate parameters based on experimentally observed concentrations. New proteomic methods can measure concentrations for all molecular species in a pathway; this creates a new opportunity to decompose the optimization of rate parameters.

Results: In contrast with conventional parameter estimation methods that minimize the disagreement between simulated and observed concentrations, the SPEDRE method fits spline curves through observed concentration points, estimates derivatives and then matches the derivatives to the production and consumption of each species. This reformulation of the problem permits an extreme decomposition of the high-dimensional optimization into a product of low-dimensional factors, each factor enforcing the equality of one ODE at one time slice. Coarsely discretized solutions to the factors can be computed systematically. Then the discrete solutions are combined using loopy belief propagation, and refined using local optimization. SPEDRE has unique asymptotic behaviour with runtime polynomial in the number of molecules and timepoints, but exponential in the degree of the biochemical network. SPEDRE performance is comparatively evaluated on a novel model of Akt activation dynamics including redox-mediated inactivation of PTEN (phosphatase and tensin homologue).

Availability and implementation: Web service, software and supplementary information are available at www.LtkLab.org/SPEDRE

Supplementary information:Supplementary data are available at Bioinformatics online.

Contact: LisaTK@nus.edu.sg

Spatial gene expression quantification: a tool for analysis of in situ hybridizations in sea anemone Nematostella vectensis

BACKGROUND

Spatial gene expression quantification is required for modeling gene regulation in developing organisms. The fruit fly Drosophila melanogaster is the model system most widely applied for spatial gene expression analysis due to its unique embryonic properties: the shape does not change significantly during its early cleavage cycles and most genes are differentially expressed along a straight axis. This system of development is quite exceptional in the animal kingdom.In the sea anemone Nematostella vectensis the embryo changes its shape during early development; there are cell divisions and cell movement, like in most other metazoans. Nematostella is an attractive case study for spatial gene expression since its transparent body wall makes it accessible to various imaging techniques.

FINDINGS

Our new quantification method produces standardized gene expression profiles from raw or annotated Nematostella in situ hybridizations by measuring the expression intensity along its cell layer. The procedure is based on digital morphologies derived from high-resolution fluorescence pictures. Additionally, complete descriptions of nonsymmetric expression patterns have been constructed by transforming the gene expression images into a three-dimensional representation.

CONCLUSIONS

We created a standard format for gene expression data, which enables quantitative analysis of in situ hybridizations from embryos with various shapes in different developmental stages. The obtained expression profiles are suitable as input for optimization of gene regulatory network models, and for correlation analysis of genes from dissimilar Nematostella morphologies. This approach is potentially applicable to many other metazoan model organisms and may also be suitable for processing data from three-dimensional imaging techniques.

Modeling of gap gene expression in Drosophila Kruppel mutants

The segmentation gene network in Drosophila embryo solves the fundamental problem of embryonic patterning: how to establish a periodic pattern of gene expression, which determines both the positions and the identities of body segments. The gap gene network constitutes the first zygotic regulatory tier in this process. Here we have applied the systems-level approach to investigate the regulatory effect of gap gene Kruppel (Kr) on segmentation gene expression. We acquired a large dataset on the expression of gap genes in Kr null mutants and demonstrated that the expression levels of these genes are significantly reduced in the second half of cycle 14A. To explain this novel biological result we applied the gene circuit method which extracts regulatory information from spatial gene expression data. Previous attempts to use this formalism to correctly and quantitatively reproduce gap gene expression in mutants for a trunk gap gene failed, therefore here we constructed a revised model and showed that it correctly reproduces the expression patterns of gap genes in Kr null mutants. We found that the remarkable alteration of gap gene expression patterns in Kr mutants can be explained by the dynamic decrease of activating effect of Cad on a target gene and exclusion of Kr gene from the complex network of gap gene interactions, that makes it possible for other interactions, in particular, between hb and gt, to come into effect. The successful modeling of the quantitative aspects of gap gene expression in mutant for the trunk gap gene Kr is a significant achievement of this work. This result also clearly indicates that the oversimplified representation of transcriptional regulation in the previous models is one of the reasons for unsuccessful attempts of mutant simulations.

Systems biology is aimed to develop an understanding of biological function or process as a system of interacting components. Here we apply the systems-level approach to understand how the blueprints for segments in the fruit fly Drosophila embryo arise. We obtain gene expression data and use the gene circuits method which allow us to reconstruct the segment determination process in the computer. To understand the system we need not only to describe it in detail, but also to comprehend what happens when certain stimuli or disruptions occur. Previous attempts to model segmentation gene expression patterns in a mutant for a trunk gap gene were unsuccessful. Here we describe the extension of the model that allows us to solve this problem in the context of Kruppel (Kr) gene. We show that remarkable alteration of gap gene expression patterns in Kr mutants can be explained by dynamic decrease of the activating effect of Cad on a target gene and exclusion of Kr from the complex network of gap gene interactions, that makes it possible for other interactions, in particular between hb and gt, to come into effect.

Efficient reverse-engineering of a developmental gene regulatory network

Understanding the complex regulatory networks underlying development and evolution of multi-cellular organisms is a major problem in biology. Computational models can be used as tools to extract the regulatory structure and dynamics of such networks from gene expression data. This approach is called reverse engineering. It has been successfully applied to many gene networks in various biological systems. However, to reconstitute the structure and non-linear dynamics of a developmental gene network in its spatial context remains a considerable challenge. Here, we address this challenge using a case study: the gap gene network involved in segment determination during early development of Drosophila melanogaster. A major problem for reverse-engineering pattern-forming networks is the significant amount of time and effort required to acquire and quantify spatial gene expression data. We have developed a simplified data processing pipeline that considerably increases the throughput of the method, but results in data of reduced accuracy compared to those previously used for gap gene network inference. We demonstrate that we can infer the correct network structure using our reduced data set, and investigate minimal data requirements for successful reverse engineering. Our results show that timing and position of expression domain boundaries are the crucial features for determining regulatory network structure from data, while it is less important to precisely measure expression levels. Based on this, we define minimal data requirements for gap gene network inference. Our results demonstrate the feasibility of reverse-engineering with much reduced experimental effort. This enables more widespread use of the method in different developmental contexts and organisms. Such systematic application of data-driven models to real-world networks has enormous potential. Only the quantitative investigation of a large number of developmental gene regulatory networks will allow us to discover whether there are rules or regularities governing development and evolution of complex multi-cellular organisms.

Life's attractors : understanding developmental systems through reverse engineering and in silico evolution

We propose an approach to evolutionary systems biology which is based on reverse engineering of gene regulatory networks and in silico evolutionary simulations. We infer regulatory parameters for gene networks by fitting computational models to quantitative expression data. This allows us to characterize the regulatory structure and dynamical repertoire of evolving gene regulatory networks with a reasonable amount of experimental and computational effort. We use the resulting network models to identify those regulatory interactions that are conserved, and those that have diverged between different species. Moreover, we use the models obtained by data fitting as starting points for simulations of evolutionary transitions between species. These simulations enable us to investigate whether such transitions are random, or whether they show stereotypical series of regulatory changes which depend on the structure and dynamical repertoire of an evolving network. Finally, we present a case study-the gap gene network in dipterans (flies, midges, and mosquitoes)-to illustrate the practical application of the proposed methodology, and to highlight the kind of biological insights that can be gained by this approach.

Analysis of gap gene regulation in a 3D organism-scale model of the Drosophila melanogaster embryo

The axial bodyplan of Drosophila melanogaster is determined during a process called morphogenesis. Shortly after fertilization, maternal bicoid mRNA is translated into Bicoid (Bcd). This protein establishes a spatially graded morphogen distribution along the anterior-posterior (AP) axis of the embryo. Bcd initiates AP axis determination by triggering expression of gap genes that subsequently regulate each other's expression to form a precisely controlled spatial distribution of gene products. Reaction-diffusion models of gap gene expression on a 1D domain have previously been used to infer complex genetic regulatory network (GRN) interactions by optimizing model parameters with respect to 1D gap gene expression data. Here we construct a finite element reaction-diffusion model with a realistic 3D geometry fit to full 3D gap gene expression data. Though gap gene products exhibit dorsal-ventral asymmetries, we discover that previously inferred gap GRNs yield qualitatively correct AP distributions on the 3D domain only when DV-symmetric initial conditions are employed. Model patterning loses qualitative agreement with experimental data when we incorporate a realistic DV-asymmetric distribution of Bcd. Further, we find that geometry alone is insufficient to account for DV-asymmetries in the final gap gene distribution. Additional GRN optimization confirms that the 3D model remains sensitive to GRN parameter perturbations. Finally, we find that incorporation of 3D data in simulation and optimization does not constrain the search space or improve optimization results.

Whole-embryo modeling of early segmentation in Drosophila identifies robust and fragile expression domains

Segmentation of the Drosophila melanogaster embryo results from the dynamic establishment of spatial mRNA and protein patterns. Here, we exploit recent temporal mRNA and protein expression measurements on the full surface of the blastoderm to calibrate a dynamical model of the gap gene network on the entire embryo cortex. We model the early mRNA and protein dynamics of the gap genes hunchback, Kruppel, giant, and knirps, taking as regulatory inputs the maternal Bicoid and Caudal gradients, plus the zygotic Tailless and Huckebein proteins. The model captures the expression patterns faithfully, and its predictions are assessed from gap gene mutants. The inferred network shows an architecture based on reciprocal repression between gap genes that can stably pattern the embryo on a realistic geometry but requires complex regulations such as those involving the Hunchback monomer and dimers. Sensitivity analysis identifies the posterior domain of giant as among the most fragile features of an otherwise robust network, and hints at redundant regulations by Bicoid and Hunchback, possibly reflecting recent evolutionary changes in the gap-gene network in insects.

The role of regulated mRNA stability in establishing bicoid morphogen gradient in Drosophila embryonic development

The Bicoid morphogen is amongst the earliest triggers of differential spatial pattern of gene expression and subsequent cell fate determination in the embryonic development of Drosophila. This maternally deposited morphogen is thought to diffuse in the embryo, establishing a concentration gradient which is sensed by downstream genes. In most model based analyses of this process, the translation of the bicoid mRNA is thought to take place at a fixed rate from the anterior pole of the embryo and a supply of the resulting protein at a constant rate is assumed. Is this process of morphogen generation a passive one as assumed in the modelling literature so far, or would available data support an alternate hypothesis that the stability of the mRNA is regulated by active processes? We introduce a model in which the stability of the maternal mRNA is regulated by being held constant for a length of time, followed by rapid degradation. With this more realistic model of the source, we have analysed three computational models of spatial morphogen propagation along the anterior-posterior axis: (a) passive diffusion modelled as a deterministic differential equation, (b) diffusion enhanced by a cytoplasmic flow term; and (c) diffusion modelled by stochastic simulation of the corresponding chemical reactions. Parameter estimation on these models by matching to publicly available data on spatio-temporal Bicoid profiles suggests strong support for regulated stability over either a constant supply rate or one where the maternal mRNA is permitted to degrade in a passive manner.

Mathematical modeling of gene expression: a guide for the perplexed biologist

The detailed analysis of transcriptional networks holds a key for understanding central biological processes, and interest in this field has exploded due to new large-scale data acquisition techniques. Mathematical modeling can provide essential insights, but the diversity of modeling approaches can be a daunting prospect to investigators new to this area. For those interested in beginning a transcriptional mathematical modeling project, we provide here an overview of major types of models and their applications to transcriptional networks. In this discussion of recent literature on thermodynamic, Boolean, and differential equation models, we focus on considerations critical for choosing and validating a modeling approach that will be useful for quantitative understanding of biological systems.

The gap gene network

Gap genes are involved in segment determination during the early development of the fruit fly Drosophila melanogaster as well as in other insects. This review attempts to synthesize the current knowledge of the gap gene network through a comprehensive survey of the experimental literature. I focus on genetic and molecular evidence, which provides us with an almost-complete picture of the regulatory interactions responsible for trunk gap gene expression. I discuss the regulatory mechanisms involved, and highlight the remaining ambiguities and gaps in the evidence. This is followed by a brief discussion of molecular regulatory mechanisms for transcriptional regulation, as well as precision and size-regulation provided by the system. Finally, I discuss evidence on the evolution of gap gene expression from species other than Drosophila. My survey concludes that studies of the gap gene system continue to reveal interesting and important new insights into the role of gene regulatory networks in development and evolution.

A software tool to model genetic regulatory networks. Applications to the modeling of threshold phenomena and of spatial patterning in Drosophila

We present a general methodology in order to build mathematical models of genetic regulatory networks. This approach is based on the mass action law and on the Jacob and Monod operon model. The mathematical models are built symbolically by the Mathematica software package GeneticNetworks. This package accepts as input the interaction graphs of the transcriptional activators and repressors of a biological process and, as output, gives the mathematical model in the form of a system of ordinary differential equations. All the relevant biological parameters are chosen automatically by the software. Within this framework, we show that concentration dependent threshold effects in biology emerge from the catalytic properties of genes and its associated conservation laws. We apply this methodology to the segment patterning in Drosophila early development and we calibrate the genetic transcriptional network responsible for the patterning of the gap gene proteins Hunchback and Knirps, along the antero-posterior axis of the Drosophila embryo. In this approach, the zygotically produced proteins Hunchback and Knirps do not diffuse along the antero-posterior axis of the embryo of Drosophila, developing a spatial pattern due to concentration dependent thresholds. This shows that patterning at the gap genes stage can be explained by the concentration gradients along the embryo of the transcriptional regulators.

Organism-scale modeling of early Drosophila patterning via bone morphogenetic proteins

Advances in image acquisition and informatics technology have led to organism-scale spatiotemporal atlases of gene expression and protein distributions. To maximize the utility of this information for the study of developmental processes, a new generation of mathematical models is needed for discovery and hypothesis testing. Here, we develop a data-driven, geometrically accurate model of early Drosophila embryonic bone morphogenetic protein (BMP)-mediated patterning. We tested nine different mechanisms for signal transduction with feedback, eight combinations of geometry and gene expression prepatterns, and two scale-invariance mechanisms for their ability to reproduce proper BMP signaling output in wild-type and mutant embryos. We found that a model based on positive feedback of a secreted BMP-binding protein, coupled with the experimentally measured embryo geometry, provides the best agreement with population mean image data. Our results demonstrate that using bioimages to build and optimize a three-dimensional model provides significant insights into mechanisms that guide tissue patterning.

Parameter estimation and inference for stochastic reaction-diffusion systems: application to morphogenesis in D. melanogaster

BACKGROUND

Reaction-diffusion systems are frequently used in systems biology to model developmental and signalling processes. In many applications, count numbers of the diffusing molecular species are very low, leading to the need to explicitly model the inherent variability using stochastic methods. Despite their importance and frequent use, parameter estimation for both deterministic and stochastic reaction-diffusion systems is still a challenging problem.

RESULTS

We present a Bayesian inference approach to solve both the parameter and state estimation problem for stochastic reaction-diffusion systems. This allows a determination of the full posterior distribution of the parameters (expected values and uncertainty). We benchmark the method by illustrating it on a simple synthetic experiment. We then test the method on real data about the diffusion of the morphogen Bicoid in Drosophila melanogaster. The results show how the precision with which parameters can be inferred varies dramatically, indicating that the ability to infer full posterior distributions on the parameters can have important experimental design consequences.

CONCLUSIONS

The results obtained demonstrate the feasibility and potential advantages of applying a Bayesian approach to parameter estimation in stochastic reaction-diffusion systems. In particular, the ability to estimate credibility intervals associated with parameter estimates can be precious for experimental design. Further work, however, will be needed to ensure the method can scale up to larger problems.

Reverse engineering a gene network using an asynchronous parallel evolution strategy

BACKGROUND

The use of reverse engineering methods to infer gene regulatory networks by fitting mathematical models to gene expression data is becoming increasingly popular and successful. However, increasing model complexity means that more powerful global optimisation techniques are required for model fitting. The parallel Lam Simulated Annealing (pLSA) algorithm has been used in such approaches, but recent research has shown that island Evolutionary Strategies can produce faster, more reliable results. However, no parallel island Evolutionary Strategy (piES) has yet been demonstrated to be effective for this task.

RESULTS

Here, we present synchronous and asynchronous versions of the piES algorithm, and apply them to a real reverse engineering problem: inferring parameters in the gap gene network. We find that the asynchronous piES exhibits very little communication overhead, and shows significant speed-up for up to 50 nodes: the piES running on 50 nodes is nearly 10 times faster than the best serial algorithm. We compare the asynchronous piES to pLSA on the same test problem, measuring the time required to reach particular levels of residual error, and show that it shows much faster convergence than pLSA across all optimisation conditions tested.

CONCLUSIONS

Our results demonstrate that the piES is consistently faster and more reliable than the pLSA algorithm on this problem, and scales better with increasing numbers of nodes. In addition, the piES is especially well suited to further improvements and adaptations: Firstly, the algorithm's fast initial descent speed and high reliability make it a good candidate for being used as part of a global/local search hybrid algorithm. Secondly, it has the potential to be used as part of a hierarchical evolutionary algorithm, which takes advantage of modern multi-core computing architectures.

Deciphering a transcriptional regulatory code: modeling short-range repression in the Drosophila embryo

A well-defined set of transcriptional regulatory modules was created and analyzed in the Drosophila embryo.

Fractional occupancy-based models were developed to explain the interaction of short range transcriptional repressors with endogenous activators by using quantitative data from these modules.

Our fractional occupancy-based modeling uncovered specific quantitative features of short-range repressors; a complex nonlinear quenching relationship, similar quenching efficiencies for different activators, and modest levels of cooperativity

The extension of the study to endogenous enhancers highlighted several features of enhancer architecture design in Drosophila embryos.

Transcriptional regulatory information, represented by patterns of protein-binding sites on DNA, comprises an important portion of genetic coding. Despite the abundance of genomic sequences now available, identifying and characterizing this information remain a major challenge. Minor changes in protein-binding sites can have profound effects on gene expression, and such changes have been shown to underlie important aspects of disease and evolution. Thus, an important aim in contemporary systems biology is to develop a global understanding of the transcriptional regulatory code, allowing prediction of gene output based on DNA sequence information. Recent studies have focused on endogenous transcriptional regulatory sequences (Janssens et al, 2006; Zinzen et al, 2006; Segal et al, 2008); however, distinct enhancers differ in many features, including transcription factor activity, spacing, and cooperativity, making it difficult to learn the effects of individual features and generalize them to other cis-regulatory elements. We have pursued a bottom up approach to understand the mechanistic processing of regulatory elements by the transcriptional machinery, using a well-defined and characterized set of repressors and activators in Drosophila blastoderm embryos. The study focuses on the Giant, Krüppel, Knirps, and Snail proteins, which have been characterized as short-range repressors, able to act locally to interfere with activator function (quenching) (Gray et al, 1994; Arnosti et al, 1996a). Such repressors have central functions in development.

The aim our study was to enable ab initio predictions of enhancer function, given defined quantities of regulatory proteins and the sequence of the enhancer (Figure 1). We have generated a large quantitative data set using fluorescent confocal laser scanning microscopy to determine the inputs (Giant, Krüppel, and Knirps protein levels) and outputs (lacZ mRNA levels) of the regulatory elements introduced into Drosophila by transgenesis. We analyzed the effect of altering specific features of a set of related gene modules, designed to uncover critical aspects of repression, including quenching distance, cooperativity, and overall factor potency.

We generated specific descriptions for each regulatory element using fractional occupancy-based modeling and identified quantitative values for parameters affecting transcriptional regulation in vivo, and these parameters were used to build and test the model. Through this process, we uncovered earlier unknown features that allow correct predictions of regulation by short-range repressors, including a non-monotonic distance function for quenching, which implicates possible phasing effects, a modest contribution for repressor–repressor cooperativity, and similarity in repression of disparate activators.

By applying these parameters to a model of the endogenous rhomboid enhancer, we uncovered novel insights into the architecture of this enhancer (Figure 8). Our study provides essential quantitative elements of a transcriptional regulatory code that will allow extensive analysis of genomic information in Drosophila melanogaster and related organisms. Extension of these predictive models should facilitate the development of more sophisticated computational algorithms for the identification and functional characterization of novel regulatory elements. The development of such quantitative modeling tools will change our understanding of the genome from essentially a parts list to a dynamically regulated system, and will greatly facilitate studies in disease, population genetics, and evolutionary biology.

Systems biology seeks a genomic-level interpretation of transcriptional regulatory information represented by patterns of protein-binding sites. Obtaining this information without direct experimentation is challenging; minor alterations in binding sites can have profound effects on gene expression, and underlie important aspects of disease and evolution. Quantitative modeling offers an alternative path to develop a global understanding of the transcriptional regulatory code. Recent studies have focused on endogenous regulatory sequences; however, distinct enhancers differ in many features, making it difficult to generalize to other cis-regulatory elements. We applied a systematic approach to simpler elements and present here the first quantitative analysis of short-range transcriptional repressors, which have central functions in metazoan development. Our fractional occupancy-based modeling uncovered unexpected features of these proteins' activity that allow accurate predictions of regulation by the Giant, Knirps, Krüppel, and Snail repressors, including modeling of an endogenous enhancer. This study provides essential elements of a transcriptional regulatory code that will allow extensive analysis of genomic information in Drosophila melanogaster and related organisms.

Inferring Drosophila gap gene regulatory network: pattern analysis of simulated gene expression profiles and stability analysis

Background

Inference of gene regulatory networks (GRNs) requires accurate data, a method to simulate the expression patterns and an efficient optimization algorithm to estimate the unknown parameters. Using this approach it is possible to obtain alternative circuits without making any a priori assumptions about the interactions, which all simulate the observed patterns. It is important to analyze the properties of the circuits.

Findings

We have analyzed the simulated gene expression patterns of previously obtained circuits that describe gap gene dynamics during early Drosophila melanogaster embryogenesis. Using hierarchical clustering we show that amplitude variation and defects observed in the simulated gene expression patterns are linked to similar circuits, which can be grouped. Furthermore, analysis of the long-term dynamics revealed four main dynamical attractors comprising stable patterns and oscillatory patterns. In addition, we also performed a correlation analysis on the parameters showing an intricate correlation pattern.

Conclusions

The analysis demonstrates that the obtained gap gene circuits are not unique showing variable long-term dynamics and highly correlating scattered parameters. Furthermore, although the model can simulate the pattern up to gastrulation and confirms several of the known regulatory interactions, it does not reproduce the transient expression of all gap genes as observed experimentally. We suggest that the shortcomings of the model may be caused by overfitting, incomplete model description and/or missing data.

Modelling the Drosophila embryo

I provide a historical overview on the use of mathematical models to gain insight into pattern formation during early development of the fruit fly Drosophila melanogaster. It is my intention to illustrate how the aims and methodology of modelling have changed from the early beginnings of a theoretical developmental biology in the 1960s to modern-day systems biology. I show that even early modelling attempts addressed interesting and relevant questions, which were not tractable by experimental approaches. Unfortunately, their validation was severely hampered by a lack of specificity and appropriate experimental evidence. There is a simple lesson to be learned from this: we cannot deduce general rules for pattern formation from first principles or spurious reproduction of developmental phenomena. Instead, we must infer such rules (if any) from detailed and accurate studies of specific developmental systems. To achieve this, mathematical modelling must be closely integrated with experimental approaches. I report on progress that has been made in this direction in the past few years and illustrate the kind of novel insights that can be gained from such combined approaches. These insights demonstrate the great potential (and some pitfalls) of an integrative, systems-level investigation of pattern formation.

Gene circuit analysis of the terminal gap gene huckebein

The early embryo of Drosophila melanogaster provides a powerful model system to study the role of genes in pattern formation. The gap gene network constitutes the first zygotic regulatory tier in the hierarchy of the segmentation genes involved in specifying the position of body segments. Here, we use an integrative, systems-level approach to investigate the regulatory effect of the terminal gap gene huckebein (hkb) on gap gene expression. We present quantitative expression data for the Hkb protein, which enable us to include hkb in gap gene circuit models. Gap gene circuits are mathematical models of gene networks used as computational tools to extract regulatory information from spatial expression data. This is achieved by fitting the model to gap gene expression patterns, in order to obtain estimates for regulatory parameters which predict a specific network topology. We show how considering variability in the data combined with analysis of parameter determinability significantly improves the biological relevance and consistency of the approach. Our models are in agreement with earlier results, which they extend in two important respects: First, we show that Hkb is involved in the regulation of the posterior hunchback (hb) domain, but does not have any other essential function. Specifically, Hkb is required for the anterior shift in the posterior border of this domain, which is now reproduced correctly in our models. Second, gap gene circuits presented here are able to reproduce mutants of terminal gap genes, while previously published models were unable to reproduce any null mutants correctly. As a consequence, our models now capture the expression dynamics of all posterior gap genes and some variational properties of the system correctly. This is an important step towards a better, quantitative understanding of the developmental and evolutionary dynamics of the gap gene network.

Inferring Drosophila gap gene regulatory network: a parameter sensitivity and perturbation analysis

BACKGROUND

Inverse modelling of gene regulatory networks (GRNs) capable of simulating continuous spatio-temporal biological processes requires accurate data and a good description of the system. If quantitative relations between genes cannot be extracted from direct measurements, an efficient method to estimate the unknown parameters is mandatory. A model that has been proposed to simulate spatio-temporal gene expression patterns is the connectionist model. This method describes the quantitative dynamics of a regulatory network in space. The model parameters are estimated by means of model-fitting algorithms. The gene interactions are identified without making any prior assumptions concerning the network connectivity. As a result, the inverse modelling might lead to multiple circuits showing the same quantitative behaviour and it is not possible to identify one optimal circuit. Consequently, it is important to address the quality of the circuits in terms of model robustness.

RESULTS

Here we investigate the sensitivity and robustness of circuits obtained from reverse engineering a model capable of simulating measured gene expression patterns. As a case study we use the early gap gene segmentation mechanism in Drosophila melanogaster. We consider the limitations of the connectionist model used to describe GRN Inferred from spatio-temporal gene expression. We address the problem of circuit discrimination, where the selection criterion within the optimization technique is based of the least square minimization on the error between data and simulated results.

CONCLUSION

Parameter sensitivity analysis allows one to discriminate between circuits having significant parameter and qualitative differences but exhibiting the same quantitative pattern. Furthermore, we show that using a stochastic model derived from a deterministic solution, one can introduce fluctuations within the model to analyze the circuits' robustness. Ultimately, we show that there is a close relation between circuit sensitivity and robustness to fluctuation, and that circuit robustness is rather modular than global. The current study shows that reverse engineering of GRNs should not only focus on estimating parameters by minimizing the difference between observation and simulation but also on other model properties. Our study suggests that multi-objective optimization based on robustness and sensitivity analysis has to be considered.

Systems biology: parameter estimation for biochemical models

Mathematical models of biological processes have various applications: to assist in understanding the functioning of a system, to simulate experiments before actually performing them, to study situations that cannot be dealt with experimentally, etc. Some parameters in the model can be directly obtained from experiments or from the literature. Others have to be inferred by comparing model results to experiments. In this minireview, we discuss the identifiability of models, both intrinsic to the model and taking into account the available data. Furthermore, we give an overview of the most frequently used approaches to search the parameter space.

Canalization of gene expression in the Drosophila blastoderm by gap gene cross regulation

Developing embryos exhibit a robust capability to reduce phenotypic variations that occur naturally or as a result of experimental manipulation. This reduction in variation occurs by an epigenetic mechanism called canalization, a phenomenon which has resisted understanding because of a lack of necessary molecular data and of appropriate gene regulation models. In recent years, quantitative gene expression data have become available for the segment determination process in the Drosophila blastoderm, revealing a specific instance of canalization. These data show that the variation of the zygotic segmentation gene expression patterns is markedly reduced compared to earlier levels by the time gastrulation begins, and this variation is significantly lower than the variation of the maternal protein gradient Bicoid. We used a predictive dynamical model of gene regulation to study the effect of Bicoid variation on the downstream gap genes. The model correctly predicts the reduced variation of the gap gene expression patterns and allows the characterization of the canalizing mechanism. We show that the canalization is the result of specific regulatory interactions among the zygotic gap genes. We demonstrate the validity of this explanation by showing that variation is increased in embryos mutant for two gap genes, Krüppel and knirps, disproving competing proposals that canalization is due to an undiscovered morphogen, or that it does not take place at all. In an accompanying article in PLoS Computational Biology (doi:10.1371/journal.pcbi.1000303), we show that cross regulation between the gap genes causes their expression to approach dynamical attractors, reducing initial variation and providing a robust output. These results demonstrate that the Bicoid gradient is not sufficient to produce gap gene borders having the low variance observed, and instead this low variance is generated by gap gene cross regulation. More generally, we show that the complex multigenic phenomenon of canalization can be understood at a quantitative and predictive level by the application of a precise dynamical model.

Animals have an astonishing ability to develop reliably in spite of variable conditions during embryogenesis. More than 60 years ago, it was proposed that this property of development, called canalization, results from genetic interactions that adjust biochemical reactions so as to bring about reliable outcomes. Since then, a great deal of progress has been made in understanding the buffering of genotypic and environmental variation, and individual mutations that reveal variation have been identified. However, the mechanisms by which genetic interactions produce canalization are not yet well understood, because this requires molecular data on multiple developmental determinants and models that correctly predict complex interactions. We make use of gene expression data at both high spatial and temporal resolution for the gap genes involved in the segmentation of Drosophila. We also apply a mathematical model to show that cross regulation among the gap genes is responsible for canalization in this system. Furthermore, the model predicted specific interactions that cause canalization, and the prediction was validated experimentally. Our results show that groups of genes can act on one another to reduce variation and highlights the importance of genetic networks in generating robust development.

DuringDrosophila development, the expression patterns of gap genes are much less variable than the Bicoid morphogen gradient. Modeling and experiments show that this specific instance of canalization or developmental robustness occurs by gap gene cross regulation.

Parameter estimation and determinability analysis applied to Drosophila gap gene circuits

Background

Mathematical modeling of real-life processes often requires the estimation of unknown parameters. Once the parameters are found by means of optimization, it is important to assess the quality of the parameter estimates, especially if parameter values are used to draw biological conclusions from the model.

Results

In this paper we describe how the quality of parameter estimates can be analyzed. We apply our methodology to assess parameter determinability for gene circuit models of the gap gene network in early Drosophila embryos.

Conclusion

Our analysis shows that none of the parameters of the considered model can be determined individually with reasonable accuracy due to correlations between parameters. Therefore, the model cannot be used as a tool to infer quantitative regulatory weights. On the other hand, our results show that it is still possible to draw reliable qualitative conclusions on the regulatory topology of the gene network. Moreover, it improves previous analyses of the same model by allowing us to identify those interactions for which qualitative conclusions are reliable, and those for which they are ambiguous.

Simulating calcium influx and free calcium concentrations in yeast

Yeast can proliferate in environments containing very high Ca(2+) primarily due to the activity of vacuolar Ca(2+) transporters Pmc1 and Vcx1. Yeast mutants lacking these transporters fail to grow in high Ca(2+) environments, but growth can be restored by small increases in environmental Mg(2+). Low extracellular Mg(2+) appeared to competitively inhibit novel Ca(2+) influx pathways and to diminish the concentration of free Ca(2+) in the cytoplasm, as judged from the luminescence of the photoprotein aequorin. These Mg(2+)-sensitive Ca(2+) influx pathways persisted in yvc1 cch1 double mutants. Based on mathematical models of the aequorin luminescence traces, we propose the existence in yeast of at least two Ca(2+) transporters that undergo rapid feedback inhibition in response to elevated cytosolic free Ca(2+) concentration. Finally, we show that Vcx1 helps return cytosolic Ca(2+) toward resting levels after shock with high extracellular Ca(2+) much more effectively than Pmc1 and that calcineurin, a protein phosphatase regulator of Vcx1 and Pmc1, had no detectable effects on these factors within the first few minutes of its activation. Therefore, computational modeling of Ca(2+) transport and signaling in yeast can provide important insights into the dynamics of this complex system.