September 25, 2017, Monday
University of Colorado at Boulder Search A to Z Campus Map University of Colorado at Boulder CU 
Search Links


Substrate-Enzyme Kinetics

From MathBio

Jump to: navigation, search

by Ben Noe

University of Colorado at Boulder

Updated April 29, 2012


In collaboration with Dr. David Bortz and other members of the Mathematical Biology Undergraduate Research Team, my
Change in Enzyme Concentration
research has focused on accurately modeling the kinetics of any simple substrate – enzyme reaction (in the example below, I used the hydrolysis of benzoyl-L-arginine ethyl ester by

trypsin). More specifically, I created a computer program with Matlab that demonstrates how the concentrations of each chemical species changes over time, shows the reaction rate over time, and determines when the reaction is in the standard quasi-steady-state. For the sake of completeness, the documented code for my program will follow a description of the principles and ideas involved.

We’ll start with a basic discussion of simple substrate-enzyme kinetics. In 1913, Leonor Michaelis and Maud Menten proposed the following chemical reaction to relate the substrate, enzyme, substrate-enzyme complex, and product.

                                           S+E\leftrightarrow SE\rightarrow E+P

k_{{1}},k_{{-1}}, and k_{{2}} correspond to the arrows above, which represent the rates of each step in the reaction.

Using the Law of Mass Action, where the rate of a reaction is directly proportional to the product of the concentrations of
Change in Substrate Concentration
the reactants, we can derive four differential equations to model how the concentration of each chemical species changes.
                                               {\frac  {d[S]}{dt}}=-k_{{1}}[E][S]+k_{{-1}}[ES]
                                        {\frac  {d[E]}{dt}}=-k_{{1}}[E][S]+k_{{-1}}[ES]+k_{{2}}[ES]
                                        {\frac  {d[ES]}{dt}}=k_{{1}}[E][S]-k_{{-1}}[ES]-k_{{2}}[ES]
                                                     {\frac  {d[P]}{dt}}=k_{{2}}[ES]

In 1925, G.E. Briggs and J.B.S. Haldane defined the standard quasi-steady-state (sQSSA) as the part of the reaction where the concentration of the intermediate substrate-enzyme complex does not change appreciably. In mathematical terms, Failed to parse (lexing error): \frac{d[S]}{dt} ≈ 0 . From this approximation, we can use the enzyme conservation law,

                                                         [E]+[ES]=[E]_{{o}}

and substitute to derive the velocity, or rate of product formation, of the substrate-enzyme reaction.

                                 V={\frac  {d[P]}{dt}}={\frac  {V_{{max}}[S]}{K_{{m}}+[S]}}
               Where  K_{{m}}={\frac  {k_{{-1}}+k_{{2}}}{k_{{1}}}} andV_{{max}}=k_{{2}}*[E]_{{o}}
Armed with the above knowledge, I proceeded to create a computer program in Matlab to easily study the progression of the reaction over time. By utilizing a numerical solver such as ode15s, one can easily obtain solutions for the system of four differential equations described above by inputting the appropriate reaction rates and initial concentrations.
Change in Substrate-Enzyme Complex Concentration

Unfortunately, values for k_{{1}} and k_{{-1}} are not readily available in many documented substrate-enzyme reactions. However, values for the Michaelis_Menten constant, k_{{m}}, and the rate-limiting turnover number,k_{{cat}}, have been widely studied. Therefore, the next step in creating a useful simulation was to calculate k_{{1}} and k_{{-1}} in terms of k_{{cat}} and k_{{m}} so that the user could input readily available values.

When calculating the values of k_{{1}} and k_{{-1}}, I derived the following solutions:

                                 k_{{-1}}={\frac  {-k_{{cat}}-k_{{cat}}*k_{{m}}*[SE]}{[SE]-[SE]*[S]}}
                                        k_{{1}}={\frac  {-k_{{cat}}*[SE]}{[SE]-[E]*[S]}}
However, regarding the differential equation solver ode15s in Matlab, the solutions simply had too many variables with the concentrations of the chemical species. Therefore, it extended beyond the scope of the function. Thus, we decided to approximate k_{{1}}={\frac  {k_{{cat}}}{k_{{m}}}} and take k_{{-1}} as zero. We successfully tested many reactions and found that this
Change in Product Concentration
approximation was valid with nearly identical solutions since k_{{1}} is usually negligible in comparison to the other rates.

In addition, I also created another similar function that can use k_{{1}},k_{{1}}, andk_{{2}} values if those values are available or if the above assumption does not hold for a given reaction, thus providing a slightly more accurate solution to the given problem.

The final step in creating a useful program was to develop a test to see when the reaction can be considered in the quasi-steady-state. In previous work on Michaelis-Menten kinetics, mathematicians and scientists have focused on under what conditions the steady state assumption would be valid, but have not to our knowledge developed a test with a specific
Change in Reaction Rate
epsilon value to see when the reaction was indeed in the steady state. But first we’ll look at a review of what others have worked on.

In 1989, Lee and Slemrod 1 proposed the conditions that {\frac  {E_{{o}}}{S_{{o}}}}+Km<<1 with S(0)=S_{{o}}. They also proposed that the substrate-enzyme complex could be modeled by the equation C={\frac  {[E]_{{o}}*[S]}{K_{{m}}+[S]}} during the QSSA nullcline. Furthermore, the stated the condition that t_{{c}}<<t_{{s}}, where t_{{c}} is the duration of the initial fast transient and t_{{s}} is the time for significant change in [S] in the post transient period. Lastly, they calculated under what conditions the errors in the quasi-steady state approximation would be small (for example, {\frac  {S_{{o}}}{K_{{m}}}}>>1).

In 1996, Borghans et al. 2 proposed the condition E_{{t}}+S_{{t}}>>{\frac  {k_{{2}}}{k_{{1}}}} and also introduced the following three steady states:

                                                        sQSSA      Failed to parse (lexing error):  \frac{d[C]}{dt} ≈ 0 


                                                         rQSSA      Failed to parse (lexing error):  \frac{d[S]}{dt} ≈ 0 


                               tQSSA      Failed to parse (lexing error):  \frac{d[S_{t}]}{dt} ≈ 0 
  where ,math> S_{t} = S + C + P </math>

In 2003, A.R. Tzafriri 2 proposed that the value {\frac  {K}{2}}*S_{{t}}(f(r(S_{{t}}))<<1, where K={\frac  {k_{{2}}}{k_{{1}}}} and f(r), or the Van Slyke Cullen constant, equals (1-r)-{\frac  {1}{2}}-1.

Lastly, Banghe Li, Bo Li, and Yuefeng Shen 3 proposed two quasi-steady state laws. First, they stated, in much more rigorous mathematical language with appropriate proofs than stated here, that S_{{o}}>>E_{{o}} and C=constant. Second, they said the {\frac  {d[C]}{dt}} also has to be relatively equal to zero. However, when it comes to a particular epsilon value to determine the steady state, they say, “Epsilon can be any positive number that depends on the requirement of the experiments”. Thus, former work has focused largely on under what reaction conditions will the quasi-steady state hold, but not particularly when such reactions are in the QSS. Thus, we attempt to do so here. Keeping in mind the work of former authors, our method of determining when the reaction is in the QSS is as follows:

We take the maximum value of the substrate-enzyme complex, which occurs within the steady state, and put the epsilon value at 1% of that value. Thus, anything within 1% of the maximum concentration is an indication of the solution remaining relatively constant and may be considered within the steady state.

Furthermore, we find the value of the derivative of the complex at a point very near the maximum concentration (for our example, we chose a time step of 0.1 seconds away). The reason for choosing this value is that the derivative at the maximum is zero. Then using that value, we say that if the derivative is within 1% of that value, thus indicating that the derivative is very near zero, then the reaction may be in the steady state.

In addition, if one takes the limit as the substrate approaches infinity, the concentration of the complex approaches the initial enzyme concentration E_{{o}},C={\frac  {E_{{o}}*[S]}{k_{{m}}+[S]}}. Thus, our third QSS criteria is that the concentration of the substrate enzyme complex is within 0.75 % of the initial enzyme concentration.

When all three conditions are met, the Matlab program returns the time values for when the reaction in question is in the quasi-steady state.

(In our example reaction, the quasi-steady-state was determined to last from 0.1 - 5.7 seconds)

After having established the criteria for determining when the reaction is in the QSS, we decided to explore adding forcing to the substrate differential equation. To start, we added the constant forcing alpha:

{\frac  {d[S]}{dt}}=\alpha -k_{{1}}[E][S]+k_{{-1}}[ES]

Cognizant of the fact that {\frac  {d[C]}{dt}}=-{\frac  {d[S]}{dt}}-k_{{2}}[C], I found an expression for alpha:

\alpha =k_{{1}}[E][S]-k_{{2}}[C]

Note: This result is based on previous approximations that k_{{-1}}=0.

Further analysis revealed that as one increases the value of alpha above the calculated value provided by the preceding equation, the amount of time that the reaction is in the QSS grows exponentially.

Armed with this tool, my current work is focused on sinusoidal forcing in place of alpha, thus modeling a heartbeat or periodic eating.

The code for both programs, one with a user input of k_{{m}} and k_{{cat}} and the other with a user input of k_{{1}},k_{{-1}}, and k_{{2}}, follows:


% This code is written for the hydrolysis of benzoyl-L-arginine ethyl ester by

% trypsin. However, it can be modified to fit any single substrate-enzyme reaction.


% s = substrate

% e = enzyme

% se = substrate enzyme complex

% p = product

% conc = concentration

% k = reaction rates


global k_m k_cat;  % Allowing variables to be used in other functions


% Reaction rates

k_m = 10e-6;  % Michaelis-Menten Constant

k_cat = 15;  % Rate-limiting turnover number (equal to k2)


% Initial concentrations (mol)

conce_init = 10e-5;  % Usually nanomolar

concs_init = 10e-3;  % Usually 5-6 orders of magnitude larger than conce_init

concse_init = 0;  % Standard reaction analysis assumes the condition that

concp_init = 0;  % the complex and product concentrations are initially zero


t_period = .1;  % Time period between function evaluation points

t_end = 10;  % Total time evaluated


% Calculating the solutions

[t,y] = ode15s ('system_template1',0:t_period:t_end,[concs_init,conce_init,concse_init,concp_init]);

concs = y(:,1);

conce = y(:,2);

concse = y(:,3);  % Assigning values to a column of the solution array

concp = y(:,4);


% Setting up time for graphs (in seconds) - must be same as time evaluation % in ode15s above

T = 0:t_period:t_end;


% Quasi-steady-state Approximation (Velocity of the Reaction)

V_qss =(k_cat*conce_init*concs)/(k_m + concs);


%Graphing

figure(1);

plot(T, conce);

ylabel('Enzyme Concentration (mol)');

xlabel('Time (s)');

figure(2);

plot(T, concs);

ylabel('Substrate Concentration (mol)');

xlabel('Time (s)');

figure(3);

plot(T, concse);

ylabel('Complex Concentration (mol)');

xlabel('Time (s)');

figure(4);

plot(T, concp);

ylabel('Product Concentration (mol)');

xlabel('Time (s)');

figure(6);

plot(T, V_qss);

ylabel('Reaction rate');

xlabel('Time');

Numerical = conce_init - concse;

Percentage = (Numerical/conce_init)*100;


% Test for when reaction is in quasi-steady-state

k1 = k_cat/k_m;

k_neg1 = k1*k_m - k_cat; % Equals 0

for i = 1:t_end/t_period + 1

concseprime(i) = k1*conce(i)*concs(i)-k_neg1*concse(i)-k_cat*concse(i);  % Calculating dC/dt

end

max_concse = max(concse);  % Obtaining the maximum value of the complex

interval = max_concse * 0.01;  % Setting epsilon to 1% of max

for i = 1:t_end/t_period + 1

difference(i) = max_concse - concse(i);  % Calculating the difference between the maximum value of the complex and %each value of the function

if difference(i) < interval  % Sets condition for QSS

[C, I] = max(concse);

value = concseprime(I-1) * 0.01;  % Going one element to the left since the derivative at the max should be 0

differenceprime = concseprime(I) - concseprime(i);

if differenceprime < value

if Percentage(i) < .75

disp(i*t_period - t_period)  % Returns exact time values for when reaction is in QSS

end

end

end


% This function is used in MichaelisMenten.m. It assigns differential

% equations for each chemical species to a column in a solution array,

% which is filled in MichaelisMenten by an ode solver.


function yprime = system_template1(t,y)


% Import global parameters

global k_m k_cat;  % The Michaelis-Menten constant and rate-limiting turnover number


% Assigning the concentrations of each chemical species to a column in the array

concs = y(1);  % Substrate

conce = y(2);  % Enzyme

concse = y(3);  % Substrate-Enzyme complex

concp = y(4);  % Product


k1 = k_cat/k_m;

k_neg1 = k1*k_m - k_cat; % Equals 0


% Michaelis-Menten Equations

yprime(1,1) = -k1*conce*concs+k_neg1*concse;

yprime(2,1) = -k1*conce*concs+k_neg1*concse+k_cat*concse;

yprime(3,1) = k1*conce*concs-k_neg1*concse-k_cat*concse;

yprime(4,1) = k_cat*concse;

end


% Note: k_cat = k2 for a simple substrate-enzyme reaction as explained by % Michaelis and Menten (1913)


% This code is written for the hydrolysis of benzoyl-L-arginine ethyl ester by

% trypsin. However, it can be modified to fit any single substrate-enzyme reaction.


% s = substrate

% e = enzyme

% se = substrate enzyme complex

% p = product

% conc = concentration

% k = reaction rates


global k_1 k_neg1 k_2;  % Allowing variables to be used in other functions


% Reaction rates

k_1 = 4e6;

k_neg1 = 25;

k_2 = 15;


% Initial concentrations (mol)

conce_init = 10e-5;  % Usually nanomolar

concs_init = 10e-3;  % Usually 5-6 orders of magnitude larger than conce_init

concse_init = 0;  % Standard reaction analysis assumes the condition that

concp_init = 0;  % the complex and product concentrations are initially zero


t_period = .1;  % Time period between function evaluation points

t_end = 10;  % Total time evaluated


% Calculating the solutions

[t,y] = ode15s ('system_template2',0:t_period:t_end,[concs_init,conce_init,concse_init,concp_init]);

concs = y(:,1);

conce = y(:,2);

concse = y(:,3);  % Assigning values to a column of the solution array

concp = y(:,4);


% Setting up time for graphs (in seconds) - must be same as time evaluation % in ode15s above

T = 0:t_period:t_end;

k_m = (k_neg1 + k_2)/k_1;  % Michaelis-Menten constant


% Quasi-steady-state Approximation (Velocity of the Reaction)

V_qss =(k_2*conce_init*concs)/(k_m + concs);


%Graphing

figure(1);

plot(T, conce);

ylabel('Enzyme Concentration (mol)');

xlabel('Time (s)');

figure(2);

plot(T, concs);

ylabel('Substrate Concentration (mol)');

xlabel('Time (s)');

figure(3);

plot(T, concse);

ylabel('Complex Concentration (mol)');

xlabel('Time (s)');

figure(4);

plot(T, concp);

ylabel('Product Concentration (mol)');

xlabel('Time (s)');

figure(6);

plot(T, V_qss);

ylabel('Reaction rate');

xlabel('Time');

Numerical = conce_init - concse;

Percentage = (Numerical/conce_init)*100;


% Test for when reaction is in quasi-steady-state

for i = 1:t_end/t_period + 1

concseprime(i) = k_1*conce(i)*concs(i)-k_neg1*concse(i)-k_2*concse(i);  % Calculating dC/dt

end

max_concse = max(concse);  % Obtaining the maximum value of the complex

interval = max_concse * 0.01;  % Setting epsilon to 1% of max

for i = 1:t_end/t_period + 1

difference(i) = max_concse - concse(i);  % Calculating the difference between the maximum value of the complex and %each value of the function

if difference(i) < interval  % Sets condition for QSS

[C, I] = max(concse);

value = concseprime(I-1) * 0.01;  % Going one element to the left since the derivative at the max should be 0

differenceprime = concseprime(I) - concseprime(i);

if differenceprime < value

if Percentage(i) < .75

disp(i*t_period - t_period)  % Returns exact time values for when reaction is in QSS

end

end

end


% This function is used in MichaelisMenten1.m. It assigns differential

% equations for each chemical species to a column in a solution array,

% which is filled in MichaelisMenten1 by an ode solver.


function yprime = system_template2(t,y)


% Import global parameters

global k_1 k_neg1 k_2;  % The reaction rates


% Assigning the concentrations of each chemical species to a column in the array

concs = y(1);  % Substrate

conce = y(2);  % Enzyme

concse = y(3);  % Substrate-Enzyme complex

concp = y(4);  % Product


% Michaelis-Menten Equations

yprime(1,1) = -k_1*conce*concs+k_neg1*concse;

yprime(2,1) = -k_1*conce*concs+k_neg1*concse+k_2*concse;

yprime(3,1) = k_1*conce*concs-k_neg1*concse-k_2*concse;

yprime(4,1) = k_2*concse;

end