In [ ]:
################################ NOTES ##############################ex
# Lines of code that are to be excluded from the documentation are #ex
# marked with `#ex` at the end of the line. #ex
# #ex
# To ensure that figures are displayed correctly together with widgets #ex
# in the sphinx documentation we will include screenshots of some of #ex
# the produced figures. #ex
# Do not run cells with the `display(Image('path_to_image'))` code to #ex
# avoid duplication of results in the notebook. #ex
# #ex
# Some reStructuredText 2 (ReST) syntax is included to aid in #ex
# conversion to ReST for the sphinx documentation. #ex
#########################################################################ex
notebook_dir = %pwd #ex
import pysces #ex
import psctb #ex
import numpy #ex
from os import path #ex
from IPython.display import display, Image #ex
from sys import platform #ex
%matplotlib inline
Thermokin is used to assess the kinetic and thermodynamic aspects of enzyme catalysed reactions in metabolic pathways [7,8]. It provides the functionality to automatically separate the rate equations of reversible reactions into a mass-action ($v_{ma}$) term and a combined binding ($v_{\Theta}$) and rate capacity ($v_{cap}$) term, however rate equations may be manually split into any arbitrary terms if more granularity is required. Additionally $\Gamma/K_{eq}$ is calculated automatically for reversible reactions. Subsequently, elasticity coefficients for the different rate equation terms are automatically calculated. Similar to symbolic control coefficient and control pattern expressions of Symca, the term and elasticity expressions generated by Thermokin can be inspected and manipulated with standard SymPy functionality and their values are automatically recalculated upon a steady-state recalculation.
.. note:: Here we use the word "term" to refer to the terms of the logarithmic form of a rate equation as well as to the corresponding factors of its linear (conventional) form. While not technically correct, this terminology is used in accordance to the original publication [8].
ScanFig.Assessing the kinetic and thermodynamic aspects of enzyme catalysed reactions using Thermokin requires the following steps:
Thermokin object using a PySCeS model object and (optionally) a file in which the rate equations of the model has been split into separate terms. reaction_results and the corresponding reaction name, reaction term name, or $\Gamma/K_{eq}$ name.ec_results and the corresponding elasticity coefficient name.As previously mentioned, Thermokin will attempt to automatically split the rate equations of reversible reactions into separate terms. While this feature should work for most common rate equations and does not require any user intervention or knowledge of the parameter names used in the model file, it is limited in two significant ways:
Additionally, the algorithm can fail in some instances.
For these reasons the separate rate equation terms can be manually defined in a .reqn file using a relatively simple syntax. Below follows such a file as automatically generated for the model lin4_fb.psc:
# Automatically parsed and split rate equations for model: lin4_fb.psc
# generated on: 13:49:07 12-01-2017
# Note that this is a best effort attempt that is highly dependent
# on the form of the rate equations as defined in the model file.
# Check correctness before use.
# R1 :successful separation of rate equation terms
!T{R1}{ma} X0 - S1/Keq_1
!T{R1}{bind_vc} 1.0*Vf_1*(S1/S1_05_1 + X0/X0_05_1)**(h_1 - 1.0)*(a_1*(S3/S3_05_1)**h_1 + 1)/(X0_05_1*(a_1*(S3/S3_05_1)**h_1*(S1/S1_05_1 + X0/X0_05_1)**h_1 + (S3/S3_05_1)**h_1 + (S1/S1_05_1 + X0/X0_05_1)**h_1 + 1))
!G{R1}{gamma_keq} S1/(Keq_1*X0)
# R2 :successful separation of rate equation terms
!T{R2}{ma} S1 - S2/Keq_2
!T{R2}{bind_vc} 1.0*S2_05_2*Vf_2/(S1*S2_05_2 + S1_05_2*S2 + S1_05_2*S2_05_2)
!G{R2}{gamma_keq} S2/(Keq_2*S1)
# R3 :successful separation of rate equation terms
!T{R3}{ma} S2 - S3/Keq_3
!T{R3}{bind_vc} 1.0*S3_05_3*Vf_3/(S2*S3_05_3 + S2_05_3*S3 + S2_05_3*S3_05_3)
!G{R3}{gamma_keq} S3/(Keq_3*S2)
# R4 :rate equation not included - irreversible or unknown formTwo types of "terms" can be defined in a .reqn file. The first type denoted by !T, is factor of the rate equation. When the !T terms for a reaction are multiplied together, they should result in the original rate equation.
Secondly !G terms are any arbitrary terms that could contain some useful information. Unlike the !T terms, the !G are not subject to any restrictions in terms of the value of their product or otherwise. For instance, the !G terms are used for define $\Gamma/K_{eq}$ for reversible reactions.
The syntax for !T and !G terms are as follows:
!T{%reaction_name}{%term_name} %term_expression
!G{%reaction_name}{%term_name} %term_expression%reaction_name - The name of the reaction to which the term belongs as defined in the .psc file (see the PySCeS MDL documentation).%term_name - The name of the term. While this name is arbitrary, there can be no duplication for any single reaction.%term_expression - The expression of the term. Thus using the example provided above for reaction 3 the line !T{R3}{ma} S2 - S3/Keq_3 specifies a !T term belonging to reaction 3 with the name ma and the expression S2 - S3/Keq_3.
Instantiation of a Thermokin analysis object requires PySCeS model object (PysMod) as an argument. Optionally a .reqn file can be provided that includes specifically slit rate equations. If path is provided, Thermokin will attempt to automatically split the reversible rate equations as described above and save a .reqn file at ~/Pysces/psc/%model_name.reqn. If this file already exists, ThermiKin will load it instead. Using the included lin4_fb.psc model a Thermokin session is instantiated as follows:
In [3]:
mod = pysces.model('lin4_fb')
tk = psctb.ThermoKin(mod)
Now that ThermoKin has automatically generated a .reqn file for lin4_fb.psc, we can load that file manually during instantiation as follows:
In [4]:
# This path leads to the provided rate equation file file
path_to_reqn = '~/Pysces/psc/lin4_fb.reqn'
# Correct path depending on platform - necessary for platform independent scripts
if platform == 'win32' and pysces.version.current_version_tuple() < (0,9,8):
path_to_reqn = psctb.utils.misc.unix_to_windows_path(path_to_reqn)
else:
path_to_reqn = path.expanduser(path_to_reqn)
tk = psctb.ThermoKin(mod,path_to_reqn)
If the path specified does not exist, a new .reqn file will be generated there instead.
Finally, ThermoKin can also be forced to regenerate a the .reqn file by setting the overwrite argument to True:
In [5]:
tk = psctb.ThermoKin(mod,overwrite=True)
Unlike RateChar and Symca, ThermoKin generates results immediately after instantiation. Results are organised similar to the other two modules, however, and can be found in the reaction_results and ec_results objects:
In [6]:
tk.reaction_results
Out[6]:
In [7]:
tk.ec_results
Out[7]:
Each results object contains a variety of fields containing data related to a specific term or expression and may be accessed in a similar way to the results of Symca:
In [8]:
# The binding*v_cap term of reaction 1
tk.reaction_results.J_R1_bind_vc
Out[8]:
SymPy expressions can be accessed via the expression field
In [9]:
tk.reaction_results.J_R1_bind_vc.expression
Out[9]:
In [10]:
tk.reaction_results.J_R1_bind_vc.value
Out[10]:
Additionally the latex_name, latex_expression, and parent model mod can also be accessed
In order to promote a logical and exploratory approach to investigating data generated by ThermoKin, the results are also arranged in a manner in which terms and elasticity coefficients associated with a certain reaction can be found nested within the results for that reaction. Using reaction 1 (called J_R1 to signify the fact that its rate is at steady state) as an example, results can also be accessed in the following manner:
In [11]:
# The reaction can also be accessed at the root level of the ThermoKin object
# and the binding*v_cap term is nested under it.
tk.J_R1.bind_vc
Out[11]:
In [12]:
# A reaction or term specific ec_results object is also available
tk.J_R1.bind_vc.ec_results.pecR1_X0_bind_vc
Out[12]:
In [13]:
# All the terms of a specific reaction can be accessed via `terms`
tk.J_R1.terms
Out[13]:
While each reaction/term/elasticity coefficient may be accessed in multiple ways, these fields are all references to the same result object. Modifying a term accessed in one way, therefore affects all references to the object.
In [14]:
# Original value of J_R3
tk.J_R3
Out[14]:
In [15]:
mod.reLoad()
# mod.Vf_3 has a default value of 1000
mod.Vf_3 = 0.1
# calculating new steady state
mod.doState()
In [16]:
# New value (original was 44.618)
tk.J_R3
Out[16]:
In [17]:
# resetting to default Vf_3 value and recalculating
mod.reLoad()
mod.doState()
Parameter scans can be performed in order to determine the effect of a parameter change on a
reaction rate and its individual terms or on the elasticity coefficients relating to a particular reaction and its related term elasticity coefficients (denoted as pec%reaction_%modifier_%term see Basic Usage - Syntax) . The procedures for both the "value" and "elasticity" scans are very much the same and rely on the same principles as described under Basic Usage - Plotting and Displaying Results.
To perform a parameter scan the do_par_scan method is called. This method has the following arguments:
parameter: A String representing the parameter which should be varied.scan_range: Any iterable representing the range of values over which to vary the parameter (typically a NumPy ndarray generated by numpy.linspace or numpy.logspace).scan_type: Either "elasticity" or "value" as described above (default: "value").init_return: If True the parameter value will be reset to its initial value after performing the parameter scan (default: True).par_scan: If True, the parameter scan will be performed by multiple parallel processes rather than a single process, thus speeding performance (default: False).par_engine: Specifies the engine to be used for the parallel scanning processes. Can either be "multiproc" or "ipcluster". A discussion of the differences between these methods are beyond the scope of this document, see here for a brief overview of Multiprocessing in Python. (default: "multiproc").Below we will perform a value scan of the effect of $V_{f^3}$ on the terms of reaction 1 for 200 points between 0.01 and 100000 in log space:
In [18]:
valscan = tk.J_R1.do_par_scan('Vf_3',scan_range=numpy.logspace(-2,5,200),scan_type='value')
In [19]:
valplot = valscan.plot()
# Equivalent to clicking the corresponding buttons
valplot.toggle_category('J_R1', True)
valplot.toggle_category('J_R1_bind_vc', True)
valplot.toggle_category('J_R1_gamma_keq', True)
valplot.toggle_category('J_R1_ma', True)
valplot.interact()
#remove_next
In [20]:
# To avoid duplication - do not run #ex
display(Image(path.join(notebook_dir,'images','tk_valscan.png'))) #ex
Similarly, we can perform an elasticity scan using the same parameters:
In [21]:
ecscan = tk.J_R1.do_par_scan('Vf_3',scan_range=numpy.logspace(-2,5,200),scan_type='elasticity')
.. note:: Elasticity coefficients with expression equal to zero (which will by definition have zero values regardless of any parameter values) are ommitted from the parameter scan results even though they are included in the ec_results objects.
In [22]:
ecplot = ecscan.plot()
# All term elasticity coefficients are enabled
# by default, thus only the "full" elasticity
# coefficients need to be enabled. Here we
# switch on the elasticity coefficients
# representing the sensitivity of R1 with
# respect to the substrate S1 and the inhibitor
# S3.
ecplot.toggle_category('ecR1_S1', True)
ecplot.toggle_category('ecR1_S3', True)
# The y limits are adjusted below as the elasticity
# values of this parameter scan have extremely
# large magnitudes at low Vf_3 values
ecplot.ax.set_ylim((-20,20))
ecplot.interact()
#remove_next
In [23]:
# To avoid duplication - do not run #ex
display(Image(path.join(notebook_dir,'images','tk_ecscan.png'))) #ex
In addition to being able to save parameter scan results (as previously described in Basic Usage - ScanFig), a summary of the results found in reaction_results and ec_results can be saved using the save_results method. This saves a csv file (by default) to disk to any specified location. If no location is specified, a file named tk_summary_N is saved to the ~/Pysces/$modelname/thermokin/ directory, where N is a number starting at 0:
In [24]:
tk.save_results()
save_results has the following optional arguments:
file_name: Specifies a path to save the results to. If None, the path defaults as described above.separator: The separator between fields (default: ",")The contents of the saved data file is as follows:
In [25]:
# the following code requires `pandas` to run
import pandas as pd
# load csv file at default path
results_path = '~/Pysces/lin4_fb/thermokin/tk_summary_0.csv'
# Correct path depending on platform - necessary for platform independent scripts
if platform == 'win32' and pysces.version.current_version_tuple() < (0,9,8):
results_path = psctb.utils.misc.unix_to_windows_path(results_path)
else:
results_path = path.expanduser(results_path)
saved_results = pd.read_csv(results_path)
# show first 20 lines
saved_results.head(n=20)
Out[25]: