Optimization and Control of Metabolic Systems

Copyright © IFAC Computer Applications in Biotechnology, Osaka, Japan . 1998

OPTIMIZA TION AND CONTROL OF METABOLIC SYSTEMS

Reinhart Heinrich

Humboldt-University Berlin. Institute of Biology. Theoretical Biophysics Invalidenstrasse 42. 10115 Berlin. Germany

Abstract: Evolutionary optimization principles are applied to explain the catalytic efficiency of single enzymes as well as the structural design of metabolic pathways. Special results concern (I) the optimal values of elementary rate constants of enzymatic mechanisms as functions of the external substrate and product concentrations; (2) the distribution of enzyme concentrations in metabolic pathways for states of maximal fluxes, and (3) the optimal distribution of A TP consuming and A TP producing reactions in glycolysis. Conclusions are drawn for the distribution of flux control coefficients resulting from optimization processes. Copyright © 1998 IFAC Keywords : Bio control, Coefficient, Flux. Optimization, Sequences. Time constant

parameters of metabolic pathways are analyzed . From the methodological point of view, evolutionary optImIzation is related to optImIzation in biotechnology. Also here, relevant objectives concern the increase of metabolic yield, the optimization of stability, and so on . On the other hand, there are some differences in that optimization in biotechnology is aimed at the improvement of one or few specialized functions, whereas biological evolution has mainly acted to achieve a well tuned balance between several functions . In any case, one may expect that the type of models presented in this communication is useful not only for the explanation of the structural design of metabolic pathways but also for the optimization in the context of biotechnology.

I . INTRODUCTION Metabolic systems are characterized by two distinct groups of data. One set is composed of the variables (essentially concentration and fluxes); the other set comprises the system parameters. Traditional simulation models serve to compute the system variables on the basis of given values for parameters. The question arises whether the latter quantities are also amenable to theoretical explanation. To answer this question, one should consider times scales on which the kinetic properties and stoichiometry of enzymatic systems have changed, that is, the dimension of evolution. A certain degree of understanding of the evolution of metabolic pathways may be gained by considering it as an optimization process. This view implies that metabolic systems found in living cells show some fitness properties which may be described by extremum principles. Investigation of optimization principles is meaningful at the level of individual reaction steps as well as on the level of multienzyme systems. In this communication it is studied whether the kinetic parameters of enzyme kinetic mechanisms may be explained on the assumption of maximal catalytic power. Furthermore. the implications of flux maximization concerning the distribution of kinetic

It will be shown that optimization is closely related to problems which are tackled in the framework of metabolic control analysis (MeA). From the distribution of flux control coefficients one may derive. for example, which enzymes should be amplified by genetic manipulation to give the highest effect in increasing the synthesis rate of a target biosynthetic product. On the other hand, evolutionary optimization may lead to special distributions of control coefficients.

331

L. : k_l = ~2P/q(1 + S),

2. OPTIMIZATION OF THE CATALYTIC PROPERTIES OF SINGLE ENZYMES

k_ 1 =~(1 + S);2Pq

The hypothesis that evolutionary pressure on the enzyme function was mainly directed toward an increase of the catalytic activity (Fersht, 1974; Crowley, 1975; Albery and Knowles, 1976; CornishBowden, 1976; Brocklehurst, 1977; Mavrovouniotis et al., 1990; Heinrich and Hoffmann, 1991 ; Pettersson, 1992), is strongly supported by the fact that the rates of enzi:matically catalyzed reactions are typically \06 - \0 2-fold higher than those of the corresponding uncatalyzed reactions. Obviously, such high reaction rates may only be achieved if the kinetic properties of the enzymes meet certain requirements. Some authors have stated that enzymes with optimal catalytic activity should have Michaelis constants close to the concentrations of their substrates in vivo. An alternative view is that the Km

(c) three solutions with submaximal values of one backward rate constant and one forward rate constant:

L,: kl =~2q(I+P)/S, k_l

constants k"

,

k_1 ), and EP

H

k_, = ~2(1 + S)/ Pq k, = ~2q(S + p)/q

and (d) one solution with all backward rate constants being submaximal : L III

E + P (rate

k1 , k_l :::; km '

(S, P)-affinity solution) . When the substrate is present at a high concentration, it is weakly bound to the enzyme in the optimal state (solution L ,: low S-

k_I' k, :::; k,) the mathematical analysis shows that maximization of the reaction rate V yields, depending on the concentrations Sand P, ten different optimal solutions L J :

affinity solution) . An analogous statement applies to the product (solution L. : low P-affinity solution) .

(a) three solutions with a submaximal value of one backward rate constant:

L1

:

In the special case that both the normalized reactant concentration Sand P equal unity, always solution Lw is obtained. The fourth-order equation given in

k_I=I/q , k_1 =I/q ,

L, : k_,=I/q .

eq. (Id) can then be solved analytically. One obtains the two real solutions, k_1 = ffq , and k _1 = -I , and

( la)

the (b) three solutions with submaximal values of two backward rate constants :

two

complex

solutions

k _1 =(- I/2±i.J3j2)ffq . As only the positive real solution is relevant, one obtains for the optimal rate constants

Ls : k_1= ~(S + P);q(1 + p), k_l

(Id)

The following conclusions may be derived . At very low substrate and product concentrations, optimal enzymic activity is achieved by improving the binding of Sand P to the enzyme (solution L,, : high

k_i and (2) that there are upper bounds for the

L1 :

k41+k ~I-Pk _ Jq-SP/q=O.

(all kinetic constants which do not enter these relations assume their maximal values for the indicated solutions; furthermore, rate constants are normalized with respect to their maximal values , substrate and product concentrations are normalized as khS/k, ~S, khP/k, ~P).

(evolutionary) changes of the rate constants k., and

,

:

k _2 = P/(qk :l ) , k_, = I/(qk _1k _1)

equilibrium constant remains unaffected by

individual rate constants (k l , k_, :::; k h

(Ic)

L, : k_1 =~2q(S + p)/q ,

k_,). Taking into account (I) that the

thermodynamic q =klk lk,/k _lk _lk_ ,

=~2(1 + p)/Sq

L. : kl =~2q(I+S)/P,

values tend to be large relative to the respective substrate concentrations. We have extended previous studies by considering a reaction mechanism which involves reversible binding processes of the substrate S and product P to the enzyme E and one or more reversible transformations of enzyme - intermediate complexes. Explicit expressions for the optimal values of rate constants can be derived for the three step enzymatic mechanism (Haldane mechanism): E + S H ES (rate constants k l , k_I' ES H EP (rate constants k1

(I b)

k - I =k - 2 ;;:k - ,

=~(1 + p)/q(S + p) 332

;;:~~q .

,

(2)

enzymes where the steady state flux J may be expressed analytically in the following way

Solution (2) shows some correspondence to the result of the descending staircase model proposed by Stackhouse et al. (1985) .

(6)

The Michaelis-Menten constants for a reversible three step mechanism may be rewritten in terms of elementary rate constants

where Sand P denote the concentrations of the pathway substrate and of the end product, respectively, and f ) is the characteristic time of step

=

K mS

K mP

kl , + k _,k, + k _,k_2 k,(k 2+k, +k _2)

= k 2k, + k _,k, + k_,k_2 k _, (k 2 + k _, + k _2)

j in a reference state with enzyme concentration E)

(3)

=I,

KmS

P -=1

Kmp

1987). The enzyme

may be determined by the method of Lagrange mUltipliers. From

these expressions. For the special case S = P = lone derives by using solution (2)

S

(Heinrich et aI.,

concentrations maximizing the steady-state flux J under the constraint that the total enzyme concentration E." for a metabolic pathway is limited

Optimal values for these Michaelis-constants are obtained by introducing k" from eqs. (la-d) into

for q = I : -

=E

~(J-?'(~ E -E )J=~-?=o dE L '" "" dE

(4)

)

111= 1

(7)

j

( ? : Lagrange multiplier) one obtains S _ P _ _'I' lor q » 1. - - = 2, - - = q . K mS Kmp

•

.

(5)

(8)

These relations bear the interesting fact that the optimal Michaelis constant K,nS of the substrate is of

This equation expresses the fact that in states of maximal steady state activity poor catalysts should be present in high concentrations. Moreover, it demonstrates that the optimal distribution of enzyme concentrations depends crucially on the equilibrium constants; e.g. for q, > 1 formula (8) predicts that

the same order of magnitude as the substrate concentration S, irrespective of the equilibrium constant q. This gives a strong support of the hypothesis that there is a mutual evolutionary adaptation of substrate concentrations and corresponding Michaelis constants (Heinrich and Hoffmann, 1991 ; Wilhelm et aI., 1994). This optimization procedure of single enzymes has been mechanisms of inorganic triosephosphate isomerase 1994).

states of maximal steady state activity are characterized by a decrease of enzyme concentrations toward the end of the chain.

for the catalytic activity applied to the reaction pyrophosphatase and (Klipp and Heinrich,

It is worth mentioning that optimization of the steadystate flux under the constraint of fixed total enzyme concentration is mathematically equivalent to the problem of minImIzing the total enzyme concentration at fixed steady-state flux (Heinrich et aI. , 1987; Brown, 1991).

3. OPTIMIZATION OF MUL TIENZYME SYSTEMS

According to coefficients

the

definition of flux control condition (7) may be

C; of reactions j

Maximization of Steady State Fluxes: The maximization of catalytic efficiencies as studied for single enzymes remains relevant also in the context of enzymic networks. However, due to the non linearity of most rate equations, the mathematical treatment is hampered by the fact that there are generally no explicit expressions for the parameter dependence of the performance functions . Let us consider the simple case of an unbranched pathway of r nonsaturated

rewritten as

Cl = E j dJ ) J dE )

=).,!i , C l = E J

)

j

EM,

(9 a,b)

Expression (9b) follows from eq . (9a) by taking into account the summation theorem of metabolic control analysis . Relation (9b) means that in states of

333

location of ATP-consuming and A TP-producing steps. If, for simplicity's sake, glycolysis is considered as an unbranched chain of reactions, the analysis may be based on eq. (6) for the steady state flux J (the glycolytic flux in this case). This formula can be applied also to chains with bimolecular reactions involving cofactors, if they are considered as external reactants .

maximal steady state activity the normalized control coefficients and the optimal enzyme concentrations in unbranched pathways show the same distribution . Recently, it has been shown that relation (9b) may be reformulated in such a way that it is applicable for metabolic systems of any complexity (Klipp and Heinrich, 1997, unpublished results). In states of minimal total enzyme concentration and at given steady state fluxes the control coefficients and the enzymes concentration are related as follows

Denoting by a and b the number of ATP consuming reactions and ATP producing reactions, respectively, the ATP-production rate is related to the glycolytic flux by J ATP= (b - a)J. To identify the optimal

(10)

where C denotes a controL matrix, that is, its elements C,~ are the flux control coefficients of the enzyme

structural design acccording to the principle J ATP = max, the kinetic properties of chains with

E J for the steady state flux J , ; the superscript T

different numbers and different locations of coupled reactions are compared. Taking into account that a coupling of a reaction to ATP consumption or ATP production will change its thermodynamic properties one may derive the following conclusions:

stands for the transpose of the matrix and E denotes the vector of enzyme concentrations E = (El, .... E,); cf. Heinrich and Schuster, 1996. Optimal Stoichiometries: It may be argued that evolutionary optimization of the kinetic properties was nothing else than a fine-tuning which guaranteed the efficient interplay of enzymes within the pathways whose basic structure had evolved in a much earlier stage of evolution. The question arises of whether the topoLogy of enzymatic systems may also be described as a result of evolutionary optimization.

[ I a] The replacement of an uncoupLed reaction by an ATP consuming reaction (i.e.. a ~ a + I) at any site increases the glycolytic rate J. [ I b] The replacement of an uncoupled reaction by an A TP producing reaction (i.e .. b ~ b + I) decreases J .

[2] J as well as J ATP are increased first by an exchange of an ATP producing reaction i for an uncoupled reaction at site m with i < m and second by an exchange of an ATP consuming reaction at site j for an uncoupled reactions at site m with m < j .

The problem of predicting optimal stoichiometric properties has been stressed in the pioneering work of Melendez-Hevia and Torres (1988) . Analyzing the structure of the nonoxidative phase of the pentose phosphate pathway, they came to the conclusion that the reduction of the number of reaction steps in the transformation of an initial substrate S into an end product P may be considered as a general principle of evolutionary optimization of metabolic pathways.

From [2] it follows that J ATP becomes maximum when all ATP producing reactions are located at the lower end of the chain and all ATP consuming reactions are located at the upper end of the chain (cf. Scheme below)

Recently, we have studied whether some structural properties of glycolysis may be understood on the basis of optimization principles (Heinrich et af.. 1997). A remarkable feature of the stoichiometry of glycolysis is that it involves ATP-consuming reactions, despite the fact that its main biological function consists in the production of ATP. It is, furthermore, striking that the two ATP-consuming reactions, hexokinase (HK) and phosphofructokinase (PFK), are located in the upper part, whereas the two ATP-producing reactions, phosphoglycerate kinase (PGK) and pyruvate kinase (PK), belong to the lower part of this pathway. Despite the fact that there are certainly various constraints concerning the chemical possibilities of converting glucose into lactate, which are in favor of this special design (cf. MelendezHevia et al., 1997), it seems worthwhile to consider also the possible kinetic advantages of such a

u ATP

"ADP

'----/

hA DP

h ATP

~

An optimum for the ATP-production rate J ATP is not only obtained by proper localization of the coupled reactions at the two ends of the chain but also by variation of their numbers a and b. Using realistic thermodynamic parameters it may be shown that an optimum for the ATP-production rate is achieved for a = I and b = 3 or a = 2 and b = 4. The latter optimum is in accord with the stoichiometric structure of glycolysis (Heinrich et af., 1997).

334

Recently an extended theory of ATP producing pathways has been presented (Heinrich and Stephani, 1997) by eliminating several oversimplifications of the previous analysis. This concerns mainly the following points: (I) the interaction of the ATP producing pathway with an external ATP consuming process is taken into account. In this way, the concentrations of adenine nucleotides may be treated as variable quantities in contrast to the previous model; (2) in the optimization studies a multitude of different pathways is taken into account which follow from strict rules concerning all possibilities of mutual arrangements of coupled and uncoupled reactions; (3) various optimal stoichometric designs are studied depending on the charactristic times of the participating processes.

Recently, genetic algorithms have been applied for the determination of optimal reaction sequences (Stephani, A., Nu n 0, J.C., Heinrich, R., unpublished results). These algorithms proved to be very efficient also for cases where the number of possible sequences is extremely high (>

With the characteristic times being equal for all types of reactions the result of the optimization procedure is the reaction sequence

REFERENCES

(AAUUPAUPPAUPAA)

10\

Obviously, biological systems have to meet not only one, but several optimization criteria. In the modelling of evolutionary optimization, it is sensible to begin with analysing particular extremum principles separately, as outlined above. The next step may then be to combine these principles and search for compromise solutions of multi-criteria problems (cf. Schuster and Heinrich, 1991).

Albery, WJ . and Knowles, J.R. (1976). Evolution of enzyme function and the development of catalytic efficiency, Biochemistry, 15, 56315640. Brocklehurst, K. (1977). Evolution of enzyme catalytic power, Biochem. 1., 163, 111-116. Brown, G.c. (1991). Total cell protein concentration as an evolutionary constraint on the metabolic control distribution in cells, J. theor. Bioi. , 153, 195-203. Cornish-Bowden, A. (1976). The effect of natural selection on enzymic catalysis, J. Mol. BioI., 101, 1-9. Crowley, P.H. (1975). Natural selection and the Michaelis constant, J. Theor. BioI., 50, 461-475. Fersht, A.R. (1974). Catalysis, binding and enzymesubstrate complementarity, Proc. R. Soc. Lond. B, 187, 379-407. Heinrich, R., Holzhtitter, H.-G. and Schuster, S. (1987). A theoretical approach to the evolution and structural design of enzymatic networks; linear enzymatic chains, branched pathways and glycolysis of erythrocytes, Bull. Math. Bioi., 49, 539-595 . Heinrich, R. and Hoffmann, E. (1991 ). Kinetic parameters of enzymatic reactions in states of maximal activity . An evolutionary approach, J. Theor. BioI., 151, 249-283. Heinrich, R., Montero, F., Klipp, E., Waddell, T. G., and Melendez-Hevia, E. (1997). Theoretical approaches to the evolutionary optimization of glycolysis; Thermodynamic and Kinetic constraints, Eur. 1. Biochem. 243, 191-20 I. Heinrich, R. and Schuster, S. (1996). The Regulation of Cellular Systems, Chapman & Hall, New York. Heinrich, R. and Stephani, A. (1997). Kinetic and thermodynamic principles determining the structural design of ATP producing systems, in: Theoretical Biophysics and Biomathematics,

(11)

where the following notations are used : U for uncoupled reactions, P and P for phosphorylation and dephosphorylation reactions, respectively, involving inorganic phosphate, and A and A for ATP consuming reactions and ATP producing reactions, respectively. This optimal structure shows that the variable concentrations of cofactors do not affect the main conclusions derived in the previous analysis. In particular, the sequence shown above has much in common with the design of real glycolysis. This concerns the location of the ATP consuming reactions in the upper part and the location of the ATP producing reactions in the lower part of the pathway, as well as the numbers of ATP molecules consumed and produced per molecule glucose (see also Stephani and Heinrich, 1998). The stoichiometric properties of the optimized sequences become more apparent if they are compared with those of nonoptimized chains, in particular with the worst one. For the same set of parameters as used above for the calculation of the optimal sequence the sequence with the lowest ATP production rate reads (PPAUAPPUAAAUAPUAPP) .

(12)

Compared with the optimal stoichiometry this reaction sequence is mainly characterized by the properties that (I) all A TP producing reactions precede the ATP consuming reactions and (2) that the phosphorylation reactions taking place by an uptake of inorganic phosphate which are thermodynamically extremely unfavourable, are located at the very beginning of the pathway.

335

Schuster, S. and Heinrich, R. (1991). Minimization of intermediate concentrations as a suggested optimality principle for biochemical networks. I. Theoretical analysis, J. Math. Bioi., 29, 425442. Stackhouse, 1., Nambiar, K.P., Burbaum, 1.1., Stauffer, D.M. and Benner, S.A. (1985). Dynamic transduction of energy and internal equilibria in enzymes: A reexamination of pyruvate kinase, J. Am. Chem. Soc., 107, 27572763 . Stephani, A. and Heinrich, R. (1998). Kinetic and thermodynamic principles determining the structural design of ATP producing systems, Bull. Math. Bioi., in press. Wilhelm, T. , Hoffmann-Klipp, E. and Heinrich, R. (1994). An evolutionary approach to enzyme kinetics: optimization of ordered mechanisms, Bull. Math. BioI., 56, 65-106.

Luo, L. , Li, Q. and Lee, W. (eds.), Hohhot, Inner Mongolia University Press, pp.133-140. Klipp, E. and Heinrich, R. (1994). Evolutionary optimization of enzyme kinetic parameters; The effect of constraints, J. theor. BioI. 171, 309323. Mavrovouniotis, M.L., Stephanopoulos, G. and Stephanopoulos, G. (1990). Estimation of upper bounds for the rates of enzymatic reactions, Chem. Eng. Commun. 93, 211-236. Melendez-Hevia, E. and Torres, N.V . (1988). Economy of design in metabolic pathways: further remarks on the game of the pentose phosphate cycle, J. theor. Bioi., 132,97-111. Melendez - Hevia, E., Waddell , T.G., Heinrich, R., and Montero, F. (1997). Theoretical approaches to the evolutionary optimization of glycolysis; Chemical analysis, Eur. J. Biochem. , 244, 527543 . Pettersson, G. (1992). Evolutionary optimization of the catalytic efficiency of enzymes, Eur. 1. Biochem., 206, 289-295 .

336