from __future__ import annotations
import matplotlib.pyplot as plt
from example_models import (
get_upper_glycolysis,
)
from modelbase2 import plot
Metabolic control analysis¶
Metabolic control analysis answers the question: what happens to the concentrations and fluxes if I slightly perturb the system?
It is thus a local measurement about which reactions hold the most control.
The most common measurements are elasticities and response coefficients.
The main difference between them is that response coefficients are calculated at steady-state, while elasticities can be calculated at every arbitrary state.
All the required functionality can be found in the modelbase2.mca module.
As an example, we will use the upper glycolysis model from the Systems Biology textbook by Klipp et al. (2005).
from modelbase2 import mca
Variable elasticities¶
Variable elasticities are the sensitivity of reactions to a small change in the concentration of a variable.
They are not a steady-state measurement and can be calculated for any arbitrary state.
Both the concs and variables arguments are optional.
If concs is not supplied, the routine will use the initial conditions from the model.
If variables is not supplied, the elasticities will be calculated for all variables.
elas = mca.variable_elasticities(
get_upper_glycolysis(),
concs={
"GLC": 0.3,
"G6P": 0.4,
"F6P": 0.5,
"FBP": 0.6,
"ATP": 0.4,
"ADP": 0.6,
},
variables=["GLC", "F6P"],
)
_ = plot.heatmap(elas)
plt.show()
Parameter elasticities¶
Similarly, parameter elasticities are the sensitivity of reactions to a small change in the concentration of a variable.
They are not a steady-state measurement and can be calculated for any arbitrary state.
Both the concs and parameters arguments are optional.
If concs is not supplied, the routine will use the initial conditions from the model.
If parameters is not supplied, the elasticities will be calculated for all parameters.
elas = mca.parameter_elasticities(
get_upper_glycolysis(),
concs={
"GLC": 0.3,
"G6P": 0.4,
"F6P": 0.5,
"FBP": 0.6,
"ATP": 0.4,
"ADP": 0.6,
},
parameters=["k1", "k2"],
)
_ = plot.heatmap(elas)
plt.show()
Response coefficients¶
Response coefficients show the sensitivity of variables and reactions at steady-state to a small change in a parameter.
If the parameter is proportional to the rate, they are also called control coefficients.
Calculation of these is run in parallel by default.
crcs, frcs = mca.response_coefficients(
get_upper_glycolysis(),
parameters=["k1", "k2", "k3", "k4", "k5", "k6", "k7"],
)
fig, (ax1, ax2) = plt.subplots(1, 2, figsize=(12, 4))
_ = plot.heatmap(crcs, ax=ax1)
_ = plot.heatmap(frcs, ax=ax2)
plt.show()
0%| | 0/7 [00:00<?, ?it/s]
14%|█▍ | 1/7 [00:00<00:03, 1.70it/s]
71%|███████▏ | 5/7 [00:01<00:00, 5.22it/s]
100%|██████████| 7/7 [00:01<00:00, 6.12it/s]
First finish line
With that you now know most of what you will need from a day-to-day basis about metabolic control analysis in modelbase2.Congratulations!
Advanced concepts / customisation¶
Normalisation¶
By default the elasticities and response coefficients are normalised, so e.g. for the response coefficients the following is calculated:
$$C^J_{v_i} = \left( \frac{dJ}{dp} \frac{p}{J} \right) \bigg/ \left( \frac{\partial v_i}{\partial p}\frac{p}{v_i} \right)$$
You can also obtain the non-normalised coefficients, by setting normalized=False, which amounts to the following calculation:
$$C^J_{v_i} = \left( \frac{dJ}{dp} \frac{p}{J} \right)$$
_ = mca.response_coefficients(
get_upper_glycolysis(),
parameters=["k1", "k2", "k3", "k4", "k5", "k6", "k7"],
normalized=False,
)
0%| | 0/7 [00:00<?, ?it/s]
14%|█▍ | 1/7 [00:00<00:02, 2.41it/s]
71%|███████▏ | 5/7 [00:00<00:00, 8.78it/s]
100%|██████████| 7/7 [00:00<00:00, 8.25it/s]
Displacement¶
We wrote above that we slightly change the value, but by how much exactly?
By default the relative change is set to 1e-4.
modelbase2 uses a central finite difference approximation, which means that in this case change the value to
$$\textrm{value} \cdot \left(1 \pm 10^{-4} \right)$$
which amounts to a change of 0.01 %.
You can set that using the displacement argument.
_ = mca.response_coefficients(
get_upper_glycolysis(),
parameters=["k1", "k2", "k3", "k4", "k5", "k6", "k7"],
displacement=1e-2,
)
0%| | 0/7 [00:00<?, ?it/s]
14%|█▍ | 1/7 [00:00<00:03, 1.63it/s]
71%|███████▏ | 5/7 [00:00<00:00, 6.25it/s]
86%|████████▌ | 6/7 [00:01<00:00, 6.35it/s]
100%|██████████| 7/7 [00:01<00:00, 6.03it/s]