Supplementary information for:

'A quantitative method for reverse engineering gene networks from microarray experiments using regulatory strengths'

Alberto de la Fuente, Paul Brazhnik and Pedro Mendes

alf@vt.edu    brazhnik@vt.edu    mendes@vt.edu

Biochemical Networks Modeling Group

     Virginia Bioinformatics Institute

We propose a method for reverse engineering gene regulatory networks from microarray gene expression data. The method is based on metabolic control analysis and requires data from microarray gene expression experiments where the rate of change of a single gene has been perturbed. Gene regulatory interactions are quantified through "regulatory strengths" which are determined from co-responses of messenger RNA to a common perturbation. The application of the method is illustrated by analyzing data produced with computer simulations of gene regulatory networks. Evidence is shown indicating that the method performs well in some cases even when perturbations are as large as two-fold. A global gene regulatory network can be uncovered by this method if applied systematically to all genes in a genome. Alternatively, it can be applied to a subset of genes in which case it identifies a phenomenological gene network that contains, beside direct interactions, contributions from indirect interactions (via the genes that were not considered or through complex formation of gene products).

To demonstrate the method we here perform the analysis on a model of a network of only three genes (figure 1). The main motivation for employing the method on simulated data is that in this case the structure of the system that produced the data is exactly known, in contrast to data produced by biological systems. We think that all methods used in gene expression analysis should be tested to simulated data before applied with confidence to the much larger biological data sets. The model is described by non-linear ordinary differential equations as given by equation 1.

Figure 1: A three gene regulatory network

Where [A], [B] and [C] are the concentrations of the mRNA species. V_a, V_b and V_c are basal rates of transcription, the K_iB and K_iB' are inhibition constants, K_aC and K_aA are activation constants, k_a, k_b, and k_c are first-order degradation constants. A reference steady state for [A], [B] and [C] was calculated setting all parameter values arbitrarily to unity, exept K_iB , K_iB' , K_aA and K_aC were set to 0.1.

The method requires a first measurement of concentrations of mRNA in a certain reference steady state. Perturbations are made by increasing the rate of transcription of single genes and the effect on the steady state concentrations of all mRNA species are measured in response to that perturbation. Perturbations of where made by increasing the rate constants V by 10% in individual experiments. The data is given in table 1.

Table 1: Steady state responses of mRNA concentrations in the 3-gene network to a 10% increase in transcription rates

where F is the fluorescence measured in the state before the purturbation (reference state) and F' the the fluorescence in the state after the perturbation, which is assumed to be proportional with the concentrations of mRNA in the reference state, [mRNA]_i and in the perturbed state [mRNA]_i' respectively.
While the co-control coefficients are defined for infinitisimal perturbations, the right hand side of the equation represents the finite estimation for the co-control coefficient. This is the form of the equation we used to calculate the co-control coefficients from the data in table 1.

Three different co-control matrices, each with respect to another variable, are calculated.

One with respect to the control coefficients on A. In this matrix the control coefficients are devided by the control coefficients of A:

To infer genetic networks we are only interested in the direct effects of genes on eachother. In the regulatory strength matrices only some rowes represent direct effects, the other elements represent regulatory effects that are indirect. The top row of the regulatory strength matrix in A gives direct effects on A, the middel row of the regulatory strength martix in B gives the direct effects on B and the bottom row in the regulatory strength matrix in C gives the direct effects on C. Mathematically the relevant rows can be joined in a single matrix as follows:

The obtained matrix is one with direct regulatory strengths in its elements. This matrix represent the genetic network and is called the genetic interaction matrix. For the system under investigation we obtained:

Using this matrix we draw the structure of the genetic network. We consider regulatory strengths lower than 0.055 to be zero, because we know that errors are obtained by making the step from infinitisimal to finite perturbations. We obtained the following network:

Figure2: The network obtained from the gene interaction matrix

Table2: The values for the regulatory coefficients determined in three different series of perturbation experiments

The model can be downloaded here as a Gepasi file. The co-control coeffients can be calculated from the steady state concentrations of mRNA using a excel spreadsheet, which can be be downloaded here. The calculated co-control matrices can be inverted to obtain the regulatory strenght matrices by using for example Mathematica or a Java based web program for matrix operations (here).

The model can be downloaded here as a Gepasi file. The co-control coeffients can be calculated from the steady state concentrations of mRNA using a excel spreadsheet, which can be be downloaded here. The calculated co-control matrices can be inverted to obtain the regulatory strenght matrices by using for example Mathematica or a Java based web program for matrix operations (here).

The model can be downloaded here as a Gepasi file. The co-control coeffients can be calculated from the steady state concentrations of mRNA using a excel spreadsheet, which can be be downloaded here. The calculated co-control matrices can be inverted to obtain the regulatory strenght matrices by using for example Mathematica or a Java based web program for matrix operations (here).

The model can be downloaded here as a Gepasi file. The co-control coeffients can be calculated from the steady state concentrations of mRNA using a excel spreadsheet, which can be be downloaded here. The calculated co-control matrices can be inverted to obtain the regulatory strenght matrices by using for example Mathematica or a Java based web program for matrix operations (here).