How to get the Stoichiometry Matrix using Tellurium

 Here is a simple need, given a reaction model how do we get hold of the stoichiometry matrix?

Consider the following simple model:

import tellurium as te
import roadrunner

r = te.loada("""
     $Xo -> S1; k1*Xo;
      S1 -> S2; k2*S1;
      S2 -> S3; k3*S2;
      S3 -> $X1; k4*S3;
      
      k1 = 0.1; k2 = 0.4;
      k3 = 0.5; k4 = 0.6;
      Xo = 1;
""")

print r.getFullStoichiometryMatrix()

Running this script by clicking on the green arrow in the toolbar will yield:

      _J0, _J1, _J2, _J3
S1 [[   1,  -1,   0,   0],
S2  [   0,   1,  -1,   0],
S3  [   0,   0,   1,  -1]]

The nice thing about this output is that the columns and rows are labeled so you know exact what is what.

What about much larger models? For example the iAF1260.xml model from the Bigg database (http://bigg.ucsd.edu:8888/models/iAF1260). This is a model of E. Coli that includes 1668 metabolites and 2382 reactions. We can download the iAF1260.xml file and load it into libRoadRunner using:

 
r = roadrunner.RoadRunner ('iAF1260.xml')

This might take up to a minute to load depending on how fast your computer is. We are assuming here that the file is located in the current directory (os.getcwd()). If not, move the file, change the current directory (using os.chdir), or use the appropriate path in the call.

Rather than print out the stoichiometry matrix (don't even try) to the screen, we'll save it to a file. Because the stoichiometry matrix is so large we will use numpy to write the matrix out as a text file:

import numpy as np
r = roadrunner.RoadRunner ('iAF1260.xml')
st = r.getFullStoichiometryMatrix()
print "Number of metabolites = ", r.getNumFloatingSpecies()
print "Number of reactions = ", r.getNumReactions()
np.savetxt ('stoich.txt', st)

Number of metabolites = 1668
Number of reactions = 2382

One can change the formatting of the output using savetxt, for example, the following will output the individual stoichiometry coefficient using 3 decimal places, 5 characters minimum, and separated by a comma.

np.savetxt ('c:\\tmp\\st.txt', st, delimiter=',', fmt='%5.3f',)

You can get the labels for the rows and columns by calling r.getFloatingSpeciesIds() and r.getReactionIds() respectively.