Enzyme Kinetics : accurate numerical solution compared to the Michaelis-Menten model approximation¶

Two Coupled Reactions: `E + S <-> ES`, and `ES -> E + P` , using real-life kinetic parameters.¶

Scenario with small amount of Enzyme, relative to the initial Substrate concentration.¶

Two other reactions being assumed negligible, and not included, are :

the reverse direction of the 2nd reaction
the non-catalyzed E + S <-> P reaction

The reaction of the enzyme Adenosinedeaminase with the substrate 2,6-Diamino-9-β-D-deoxyribofuranosyl-9-H-purine,
and the initial concentration values choosen below, all satisfy the customary Michaelis-Menten assumptions that
[E] << [S] and that the reaction rate constants satisfy k1_reverse >> k2_forward

Source of kinetic parameters: page 16 of "Analysis of Enzyme Reaction Kinetics, Vol. 1", by F. Xavier Malcata, Wiley, 2023

TAGS : "uniform compartment", "chemistry", "numerical", "enzymes"¶

In [1]:

LAST_REVISED = "Oct. 11, 2024"
LIFE123_VERSION = "1.0.0.beta.39"   # Library version this experiment is based on

In [2]:

#import set_path                    # Using MyBinder?  Uncomment this before running the next cell!

In [3]:

#import sys
#sys.path.append("C:/some_path/my_env_or_install")   # CHANGE to the folder containing your venv or libraries installation!
# NOTE: If any of the imports below can't find a module, uncomment the lines above, or try:  import set_path   

import ipynbname
import pandas as pd

from life123 import check_version, ChemData, UniformCompartment, GraphicLog, PlotlyHelper

In [4]:

check_version(LIFE123_VERSION)

OK

In [5]:

# Initialize the HTML logging (for the graphics)
log_file = ipynbname.name() + ".log.htm"    # Use the notebook base filename for the log file
                                            # IN CASE OF PROBLEMS, set manually to any desired name

# Set up the use of some specified graphic (Vue) components
GraphicLog.config(filename=log_file,
                  components=["vue_cytoscape_2"],
                  extra_js="https://cdnjs.cloudflare.com/ajax/libs/cytoscape/3.21.2/cytoscape.umd.js")

-> Output will be LOGGED into the file 'enzyme_1.log.htm'

In [ ]:

1. Accurate numerical solution¶

In [6]:

chem_data = ChemData(names=["P", "ES*"])

In [7]:

# Our Enzyme
chem_data.add_chemical(name="Adenosinedeaminase", label="E") 

# Our Substrate
chem_data.add_chemical(name="2,6-Diamino-9-β-D-deoxyribofuranosyl-9-H-purine", label="S");

In [8]:

chem_data.all_chemicals()

Out[8]:

	name	label
0	P	P
1	ES*	ES*
2	Adenosinedeaminase	E
3	2,6-Diamino-9-β-D-deoxyribofuranosyl-9-H-purine	S

In [ ]:

Specify the Kinetic Parameters¶

Source: page 16 of "Analysis of Enzyme Reaction Kinetics, Vol. 1", by F. Xavier Malcata, Wiley, 2023

In [9]:

# Reaction E + S <-> ES* , with 1st-order kinetics, 
# and a forward rate that is much faster than its revers one
chem_data.add_reaction(reactants=["E", "S"], products=["ES*"],
                       forward_rate=18., reverse_rate=100.) 

# Reaction ES* <-> E + P , with 1st-order kinetics
chem_data.add_reaction(reactants=["ES*"], products=["E", "P"],
                forward_rate=49.,reverse_rate=0)

chem_data.describe_reactions()

Number of reactions: 2 (at temp. 25 C)
0: E + S <-> ES*  (kF = 18 / kR = 100 / delta_G = 4,250.9 / K = 0.18) | 1st order in all reactants & products
1: ES* <-> E + P  (kF = 49 / kR = 0) | 1st order in all reactants & products
Set of chemicals involved in the above reactions: {'E', 'P', 'ES*', 'S'}

In [10]:

# Send a plot of the network of reactions to the HTML log file
chem_data.plot_reaction_network("vue_cytoscape_2")

[GRAPHIC ELEMENT SENT TO LOG FILE `enzyme_1.log.htm`]

In [ ]:

Set the initial concentrations of all the chemicals¶

in the scenario we're exploring, there's little enzyme relative to initial substrate concentration¶

In [11]:

S0 = 20.
E0 = 1.

In [12]:

# Here we use the "slower" preset for the variable steps, a conservative option prioritizing accuracy over speed
uc = UniformCompartment(chem_data=chem_data, preset="slower")
uc.set_conc(conc={"S": S0, "E": E0})      # Small ampount of enzyme `E`, relative to substrate `S`
uc.describe_state()

SYSTEM STATE at Time t = 0:
4 species:
  Species 0 (P). Conc: 0.0
  Species 1 (ES*). Conc: 0.0
  Species 2 (E). Conc: 1.0
  Species 3 (S). Conc: 20.0
Set of chemicals involved in reactions: {'E', 'P', 'ES*', 'S'}

In [13]:

uc.enable_diagnostics()   # To save diagnostic information about the simulation - in particular, the REACTION RATES

Advance the reactions to equilibrium¶

In [14]:

# Perform the reactions
uc.single_compartment_react(duration=1.5, initial_step=0.05)

Some steps were backtracked and re-done, to prevent negative concentrations or excessively large concentration changes
827 total step(s) taken
Number of step re-do's because of negative concentrations: 2
Number of step re-do's because of elective soft aborts: 1
Norm usage: {'norm_A': 706, 'norm_B': 710, 'norm_C': 706, 'norm_D': 706}

In [15]:

uc.plot_history(colors=['cyan', 'green', 'violet', 'red'], show_intervals=True, title_prefix="Small amout of E relative to S(0)")

plot_pandas() NOTICE: Excessive number of vertical lines (828) - only showing 1 every 6 lines

In [16]:

uc.plot_history(colors=['cyan', 'green', 'violet', 'red'], title_prefix="DETAIL at early times",
                range_x=[0, 0.5], range_y=[0, 2.])

Notice how the bound enzyme ES* quickly builds up at the very beginning, from time 0 to roughly 0.01 ... and that, in the longer term, the enzyme returns to its unbound state E

In [ ]:

What is the initial rate of production of the final reaction product `P`?¶

One could take the numerical derivative (gradient) of the time values of [P] - but no need to! Reaction rates are computed in the course of the simulation, and stored alongside other diagnostic data, if diagnostics were enabled (as we did indeed enable)

In [17]:

rates = uc.get_diagnostics().get_rxn_rates(rxn_index=1)   # We specify the reaction that involves `P`
rates

Reaction:  ES* <-> E + P

Out[17]:

	START_TIME	rate
0	0.000000	0.000000
1	0.000500	8.820000
2	0.000750	12.075109
3	0.000763	12.216716
4	0.000769	12.287060
...	...	...
822	1.468805	0.285397
823	1.475832	0.274566
824	1.482859	0.264142
825	1.489886	0.254109
826	1.496913	0.244454

827 rows × 2 columns

In [18]:

# Let's take a look at how the reaction rate varies with time
PlotlyHelper.plot_pandas(df=rates, 
                         title="Reaction rate, dP/dt, over time",
                         x_var="START_TIME", fields="rate", 
                         x_label="time", y_label="dP/dt")

In [19]:

# A closer peek at it maximum value
PlotlyHelper.plot_pandas(df=rates, 
                         title="Reaction rate, dP/dt, over time (DETAIL at early times)",
                         x_var="START_TIME", fields="rate", 
                         x_label="time", y_label="dP/dt",
                         range_x=[0,0.05], range_y=[33., 34.5])

As we saw earlier, the time it took for ES* to build up was about 0.01
Ignoring the brief initial transient phase, the reaction rate starts at about 34.1

In [ ]:

2. Comparing the results to the Michaelis-Menten model approximation¶

Let's compute some parameters used by the Michaelis-Menten model¶

for background reference, see: https://vallance.chem.ox.ac.uk/pdfs/KineticsLectureNotes.pdf (p. 20)

In [20]:

kM = (chem_data.get_forward_rate(1) + chem_data.get_reverse_rate(0)) / chem_data.get_forward_rate(0)
kM      # (49. + 100.) / 18. in this experiment

Out[20]:

8.277777777777779

In [21]:

kcat = chem_data.get_forward_rate(1)
kcat

Out[21]:

49.0

In [22]:

vmax = kcat * E0
vmax

Out[22]:

49.0

In [23]:

initial_rxn_rate = (vmax * S0) / (kM + S0)
initial_rxn_rate

Out[23]:

34.656188605108056

Not too far from the value of about 34.1 we saw earlier...
To keep in mind that the initial part of the reaction is affected by transients

Now, let's look at rate as a function of [S]; we'll compare what we computed earlier vs. as given by the approximation of the Michaelis-Menten model¶

In [ ]:

In [24]:

# Let's retrieve all the values of `S` at the various simulation time
hist = uc.get_history()[:-1]    # Dropping the last row, because no rate information is known about the next simulation step not taken!
hist

Out[24]:

	SYSTEM TIME	P	ES*	E	S	caption
0	0.000000	0.000000	0.000000	1.000000	20.000000	Initialized state
1	0.000500	0.000000	0.180000	0.820000	19.820000	1st reaction step
2	0.000750	0.002205	0.246431	0.753569	19.751364
3	0.000763	0.002356	0.249321	0.750679	19.748323
4	0.000769	0.002432	0.250756	0.749244	19.746811
...	...	...	...	...	...	...
822	1.468805	19.947437	0.005824	0.994176	0.046738
823	1.475832	19.949443	0.005603	0.994397	0.044954
824	1.482859	19.951372	0.005391	0.994609	0.043237
825	1.489886	19.953228	0.005186	0.994814	0.041586
826	1.496913	19.955014	0.004989	0.995011	0.039997

827 rows × 6 columns

In [25]:

# Now, let's put together in the same table (Pandas dataframe) the `S` values from above, and the `rates` we extracted earlier

assert len(hist) == len(rates)   # They'd better match up!

rate_and_substrate = pd.DataFrame({
    "SYSTEM TIME": hist["SYSTEM TIME"],   # We don't actually need the time, but just for clarity
    "S": hist["S"],
    "computed_rate": rates["rate"]
})

rate_and_substrate

Out[25]:

	SYSTEM TIME	S	computed_rate
0	0.000000	20.000000	0.000000
1	0.000500	19.820000	8.820000
2	0.000750	19.751364	12.075109
3	0.000763	19.748323	12.216716
4	0.000769	19.746811	12.287060
...	...	...	...
822	1.468805	0.046738	0.285397
823	1.475832	0.044954	0.274566
824	1.482859	0.043237	0.264142
825	1.489886	0.041586	0.254109
826	1.496913	0.039997	0.244454

827 rows × 3 columns

In [26]:

# Let's a column with the rate estimated by the Michaelis-Menten model
rate_and_substrate["Michaelis_rate"] = vmax * rate_and_substrate["S"] / (kM + rate_and_substrate["S"])
rate_and_substrate

Out[26]:

	SYSTEM TIME	S	computed_rate	Michaelis_rate
0	0.000000	20.000000	0.000000	34.656189
1	0.000500	19.820000	8.820000	34.564299
2	0.000750	19.751364	12.075109	34.528950
3	0.000763	19.748323	12.216716	34.527380
4	0.000769	19.746811	12.287060	34.526599
...	...	...	...	...
822	1.468805	0.046738	0.285397	0.275111
823	1.475832	0.044954	0.274566	0.264664
824	1.482859	0.043237	0.264142	0.254610
825	1.489886	0.041586	0.254109	0.244934
826	1.496913	0.039997	0.244454	0.235622

827 rows × 4 columns

In [27]:

# Let's see how our computed rate compares with the approximations from the Michaelis-Menten model
PlotlyHelper.plot_pandas(df=rate_and_substrate, x_var="S", fields=["computed_rate", "Michaelis_rate"],
                         title="Reaction rate, dP/dt, as a function of Substrate concentration",
                         y_label="dP/dt")

Virtually overlapped plots, except at the very right!

In this scenario, the Michaelis-Menten model is remarkable accurate EXCEPT at the very early times of the reaction (the large values of `S`), which is the brief initial transient phase when `ES*` builds up from zero¶

This is as expected, because the kinetic parameters of this reaction satisfy the customary Michaelis-Menten assumptions that `[E] << [S]` and that the rates satisfy `k1_reverse >> k2_forward`¶

In [ ]:

Enzyme Kinetics : accurate numerical solution compared to the Michaelis-Menten model approximation¶

Two Coupled Reactions: E + S <-> ES*, and ES* -> E + P , using real-life kinetic parameters.¶

Scenario with small amount of Enzyme, relative to the initial Substrate concentration.¶

TAGS : "uniform compartment", "chemistry", "numerical", "enzymes"¶

1. Accurate numerical solution¶

Specify the Kinetic Parameters¶

Set the initial concentrations of all the chemicals¶

in the scenario we're exploring, there's little enzyme relative to initial substrate concentration¶

Advance the reactions to equilibrium¶

What is the initial rate of production of the final reaction product P?¶

2. Comparing the results to the Michaelis-Menten model approximation¶

Let's compute some parameters used by the Michaelis-Menten model¶

Now, let's look at rate as a function of [S]; we'll compare what we computed earlier vs. as given by the approximation of the Michaelis-Menten model¶

In this scenario, the Michaelis-Menten model is remarkable accurate EXCEPT at the very early times of the reaction (the large values of S), which is the brief initial transient phase when ES* builds up from zero¶

This is as expected, because the kinetic parameters of this reaction satisfy the customary Michaelis-Menten assumptions that [E] << [S] and that the rates satisfy k1_reverse >> k2_forward¶

Two Coupled Reactions: `E + S <-> ES`, and `ES -> E + P` , using real-life kinetic parameters.¶

What is the initial rate of production of the final reaction product `P`?¶

In this scenario, the Michaelis-Menten model is remarkable accurate EXCEPT at the very early times of the reaction (the large values of `S`), which is the brief initial transient phase when `ES*` builds up from zero¶

This is as expected, because the kinetic parameters of this reaction satisfy the customary Michaelis-Menten assumptions that `[E] << [S]` and that the rates satisfy `k1_reverse >> k2_forward`¶