ME/CFS, a mathematical model

Robert Phair and the trap

Many of my readers are probably aware of the attempts that are currently being made to mathematically simulate the energy metabolism of ME/CFS patients, integrating metabolic data with genetic data. In particular, dr. Robert Phair has developed a mathematical model of the main metabolic pathways involved in energy conversion, from energy stored in the chemical bonds of big molecules like glucose, fatty acids, and amino acids, to energy stored in adenosine triphosphate (ATP), ready to use. Phair, who is an engineer, determined the differential equations that rule this huge amount of chemical reactions and adapted them to the genetic profile found in ME/CFS patients. Genetic data are the result of the analysis of the exomes of more than 40 patients (I am among them). He found, in particular, that the enzyme ATP synthase (one step of the mitochondrial electron transport chain, called complex V) presents a damaging variant that is found in less than 30% of healthy controls while being detected in more than 60% of ME/CFS patients. He found other enzymes (I don’t know their exact number) of energy metabolism meeting these same criteria. Setting a reduced activity for these enzymes in his mathematical model, he found that, after particularly stressing events (such as infections), ME/CFS patients metabolism can fall into a low-energy state without being able to escape from it (this theory has been called the “metabolic trap hypothesis”) (R, R, R). Phair’s hypothesis is currently being experimentally tested and has not been published yet, but some years ago two physicists published an interesting mathematical model of energy metabolism during and after exercise, in ME/CFS patients compared to healthy controls (Lengert N. et Drossel B. 2015). In what follows I will describe this model and its predictions and we will also see how these differential equations look like. 

Metabolic pathways that have been analyzed

The model by Lengert and Drossel extends two previously published systems of differential equations that describe the behaviour of glycolysis, Krebs cycle (hugely simplified as a single reaction!), mitochondrial electron transport chain (described in detail), creatine kinase system, and of conversion of adenosine diphosphate (ADP) in ATP, in skeletal muscles (Korzeniewski B. et Zoladz JA. 2001), (Korzeniewski B. et Liguzinski P. 2004). They added equations for lactate accumulation and efflux out of the cell, for de novo synthesis of inosine monophosphate (IMP) during recovery, for degradation of adenosine monophosphate (AMP) into IMP, for degradation of IMP to inosine and hypoxanthine. All the pathways involved are collected in figure 1. These reactions are described by 15 differential equations and the solution is a set of 15 functions of time that represent the concentration of the main metabolites involved (such as lactate, pyruvate, ATP etc). Let’s now take a closer look at one of these equations and at the general structure of the whole system of equations.

Figure 1. This is a schematic representation of the metabolic pathways described by the mathematical model developed by Lengert and Drossel. In detail: cytosolic synthesis and degradation of ADP, AMP and IMP (left), protein kinase pathway and glycolysis (centre), electron transport chain and TCA cycle (right). From Lengert N. et Drossel B. 2015.

Figure 2. Lactate dehydrogenase is the enzyme involved in the catalysis of the conversion of lactate to pyruvate. This reaction goes in both directions.

Differential equations for chemical reactions

Let’s consider the equation used by the Authors for the reaction catalyzed by lactate dehydrogenase (the transformation of pyruvate into lactate, figure 2) where they also took into account the efflux of lactate from the cytosol. The differential equation is the following one:

where the three parameters are experimentally determined and their values are

The first one describes the activity of the enzyme lactate dehydrogenase: the higher this parameter is, the more active the enzyme is. The second one describes the inverse reaction (from lactate to pyruvate). The third one is a measure of how much lactate the cell is able to transport outside its membrane. You have probably recognized that the equation of lactate is a first-order ordinary differential one. We say “first-order” because in the equation there is only the first-order derivative of the function that we have to determine (lactate, in this case); “ordinary” refers to the fact that lactate is a function only of one variable (time, in this case). You can easily realize that an equation like that can be written as follows:

Suppose now that we had other two differential equations of this type, one for pyruvate  and one for protons (the other two functions of time that are present in the equation):

We would have a system of three ordinary differential equations like this one

The initial values of the functions that we have to determine are collected in the last row: these are the values that they have at the beginning of the simulation (t=0). In this case, these values are the concentrations of lactate, pyruvate and protons in the cytosol, at rest. The three functions of time are called the solution to the system. This kind of system of equations is an example of a Cauchy’s problem, and we know from mathematical theory that not only it has a solution, but that this solution is unique. Moreover, whereas this solution can’t be always easily found with rigorous methods, it is quite easy to solve the problem with computational methods, like the Runge-Kutta method or the Heun’s method. All that said, the system of ordinary differential equations proposed by Lengert and Drossel for energy metabolism is just like this one, with the exception that it comprises 15 equations instead of three. So, the main difficulty in this kind of simulation is not the computational aspect but the determination of the parameters (like the enzymatic ones) and of the initial values, that have to be collected from the medical literature or have to be determined experimentally, if not already available. The other problem is how to design the equations: there are several ways to build a model for a chemical reaction or for any other biological process.

The mathematical model of ME/CFS

How do we adapt to ME/CFS patients a model of energy metabolism that has been set with parameters taken from experiments performed on healthy subjects? This is a very good question, and we have seen that Robert Phair had to use genetic data from ME/CFS patients on key enzymes of energy metabolism, in order to set his model. But this data was not available when Lengert and Drossel designed their equations. So what? They looked for studies about the capacity of oxidative phosphorylation in ME/CFS patients in comparison with healthy subjects, and they found that it had been measured with different experimental settings by various groups and that the common denominator was a reduction ranging from about 35% (Myhill S et al. 2009), (Booth, N et al 2012), (Argov Z. et al. 1997), (Lane RJ. et al. 1998) to about 20% (McCully KK. et al. 1996), (McCully KK. et al. 1999). So the idea of the Authors was to multiply the enzymatic parameter of each reaction belonging to the oxidative phosphorylation by a number ranging from 0.6 (severe ME/CFS) to 1.0 (healthy person). In particular, they used a value of 0.7 for ME/CFS, in their in silico experiments.

Predictions of the mathematical model

The mathematical model was used to perform in silico exercise tests with various length and intensities. What they found was that the time of recovery in the average ME/CFS patient was always greater if compared to a healthy person. The time of recovery is defined as the time that a cell needs to replenish its content of ATP (97% of the level in resting state) after exertion. In Figure 3 you see the results of the simulation for a very short (30 seconds) and very intense exercise. As you can see, in the case of a healthy cell (on the left) the recovery time is of about 600 minutes (10 hours) whereas a cell from a person with ME/CFS (on the right) requires more than 1500 minutes (25 hours) to recover.

Figure 3. Results of the simulation for an exercise with a duration of 30 seconds and a high intensity (initial ATP consumption 300 times the resting value). On the left, the case of a healthy skeletal muscle cell, on the right the case of a cell from a person with ME/CFS whose mitochondrial reactions have a velocity reduced to 70% of the velocity of healthy control. The plot has been obtained by using the online version of the software, available here.

Another interesting result of the simulation is an increase in AMP in patients vs control (figure 3, orange line). This is due to the compensatory use of the two metabolic pathways in figure 4: the reaction catalyzed by adenylate kinase, in which two molecules of ADP are used to produce a molecule of ATP and a molecule of AMP; and the reaction catalyzed by AMP deaminase, that degrades AMP to IMP (that is then converted to inosine and hypoxanthine). These two reactions are used by ME/CFS patients more than in healthy control, in order to increase the production of ATP outside mitochondria.

Figure 4. The metabolic pathway on the left is used by ME/CFS patients more than in control to increase the production of ATP outside mitochondria, according to the present model. The pathway on the right then degrades AMP to IMP.

If we give a closer look at the concentrations of AMP and IMP in the 4 hours following the exertion (figure 5), we actually see an increased production of IMP (green line) and AMP (orange line) in skeletal muscles of ME/CFS patients (on the right) vs controls (left).

Figure 5. The same as figure 3, but magnified in order to have a closer look at concentrations during the 4 hours following the exertion. The healthy cell is on the left, whereas the cells from a person with ME/CFS is on the right.

A further compensatory pathway used by patients (according to this model) is the production of ATP from ADP by the enzyme creatine kinase (figure 6). This is another way that we have on our cells to produce ATP in the cytosol without the help of mitochondria. In this model of ME/CFS, there is an increase in the use of this pathway, which leads to a decrease in the cellular concentration of phosphocreatine and an increase in the cellular concentration of creatine (figure 7).

Figure 6. The reaction catalyzed by creatine kinase: a molecule of ADP has converted to ATP thanks to the phosphate group carried by phosphocreatine.

Figure 7. The concentration of phosphocreatine in the cytosol of skeletal muscle cells is lower in ME/CFS (right) versus control (left) during and after exercise. This is due to the greater use of this molecule to produce ATP anaerobically in ME/CFS metabolism vs control. Parameters for this simulation are the same as described in figure 3.

Comparison with available metabolic data

I am curious to see if data from the various metabolomic studies done after the publication of the model by Lengert and Drossel are in agreement with it. I will discuss this topic in another article because I still have to study this aspect. I would just point out that if we assumed true the high rate of IMP degradation proposed in this model, we would probably find a high level of hypoxanthine in the blood of patients, compared to controls, whereas this metabolite is decreased in patients, according to one study (Armstrong CW et al. 2015).

Comparison with the model by Phair

The metabolic model by Robert Phair will probably give a more accurate simulation of the energy metabolism of patients if compared with the system of ordinary differential equations that we have discussed in this article. And there are two reasons for that. The first one is that Phair has included equations also for fatty acid beta-oxidation, pentose phosphate pathway, and NAD synthesis from vitamin niacin. The other one is that, whereas the two German physicists reduced the velocities of all the enzymatic reactions that happen in mitochondria, Phair has genetic data for every enzyme involved in these reactions for a group of ME/CFS patients and thus he can determine and set the actual degree of activity for each enzyme. But there is a further level required in order to bring the mathematical simulation closer to the reality: gene expression. We know, for instance, that in ME/CFS patients there is a higher than normal expression of aconitase (an enzyme belonging to the TCA cycle) and of ATP synthase (Ciregia F et al 2016) and this should be taken into account in a simulation of ME/CFS patients energy metabolism. Note that ATP synthase is exactly the enzyme that Phair has found to be genetically damaged in patients, and this makes perfect sense: if an enzyme has reduced activity, its reaction can be speeded up by expressing more copies of the enzyme itself.

One could expect that, in a near future, genetic data and gene expression data from each of us will be used to set mathematical models for metabolic pathways, in order to build a personalized model of metabolism that might be used to define, study and correct human diseases in a personalized fashion. But we would need a writer of sci-fiction in order to tell this chapter of the future of medicine.