1 Introduction

# Enzyme allocation problems in kinetic metabolic networks: Optimal solutions are elementary flux modes

## Abstract

The survival and proliferation of cells and organisms require a highly coordinated allocation of cellular resources to ensure the efficient synthesis of cellular components. In particular, the total enzymatic capacity for cellular metabolism is limited by finite resources that are shared between all enzymes, such as cytosolic space, energy expenditure for amino-acid synthesis, or micro-nutrients. While extensive work has been done to study constrained optimization problems based only on stoichiometric information, mathematical results that characterize the optimal flux in kinetic metabolic networks are still scarce. Here, we study constrained enzyme allocation problems with general kinetics, using the theory of oriented matroids. We give a rigorous proof for the fact that optimal solutions of the non-linear optimization problem are elementary flux modes. This finding has significant consequences for our understanding of optimality in metabolic networks as well as for the identification of metabolic switches and the computation of optimal flux distributions in kinetic metabolic networks.

Enzyme allocation problems in kinetic metabolic networks: Optimal solutions are elementary flux modes

Stefan Müller, Georg Regensburger, Ralf Steuer

Johann Radon Institute for Computational and Applied Mathematics,

Altenbergerstraße 69, 4040 Linz, Austria

CzechGlobe,

Academy of Sciences of the Czech Republic

Bělidla 986/4a, 603 00 Brno, Czech Republic

Institute for Theoretical Biology,

Humboldt University Berlin

Invalidenstraße 43, 10115 Berlin, Germany

stefan.mueller@ricam.oeaw.ac.at

Keywords. metabolic optimization, enzyme kinetics, oriented matroid, elementary vector, conformal sum

## 1 Introduction

Living organisms are under constant evolutionary pressure to survive and reproduce in complex environments. As a direct consequence, cellular pathways are often assumed to be highly adapted to their respective tasks, given the biochemical and biophysical constraints of their environment. Optimality principles have proven to be powerful methods to study and understand the large-scale organization of metabolic pathways [5, 12, 28, 34, 17]. A variety of recent computational techniques, such as flux-balance analysis (FBA), seek to identify metabolic flux distributions that maximize given objective functions, such as ATP regeneration or biomass yield, under a set of linear constraints. As one of their prime merits, FBA and related stoichiometric methods, including the generalization to time-dependent metabolism [16, 1], only require knowledge of the stoichiometry of a metabolic network – data that are available for an increasing number of organisms in the form of large-scale metabolic reconstructions [23, 25].

However, despite their explanatory and predictive success, constraint-based stoichiometric methods also have inherent limits. Specifically, FBA and related methods typically maximize stoichiometric yield. That is, the value of a designated output flux is maximized, given a set of limiting input fluxes. As emphasized in a number of recent studies, the assumption of maximal stoichiometric yield is not necessarily a universal principle of metabolic network function [32, 29, 10, 17]. Quite on the contrary, examples of seemingly suboptimal metabolic behavior, at least from a stoichiometric perspective, are well-known for many decades. Among the most prominent instances are the Warburg and the Crabtree effect [39, 7, 13]. Under certain circumstances, cells utilize a fermentative metabolism rather than aerobic respiration to regenerate ATP, despite the presence of oxygen and despite its significantly lower stoichiometric yield of ATP per amount of glucose consumed.

To account for such seemingly suboptimal behavior, several modifications and extensions of FBA have been developed recently. Conventional FBA is augmented with additional principles concerning limited cytosolic volume [4, 35, 36, 37, 33], membrane occupancy [41], and other, more general capacity constraints [29, 11]. Each of these extensions allows for new insights into suboptimal stoichiometric behavior, and additional constraints often also induce the utilization of pathways with lower stoichiometric yield. However, none of the modifications of FBA addresses a metabolic network as a genuine dynamical system with particular kinetics that depend on a number of parameters. The neglect of the dynamical nature is a direct consequence of the extensive data requirements for parametrizing enzymatic reaction rates. Correspondingly, and despite its importance to understand metabolic optimality, only few mathematically rigorous results are currently available that allow to characterize solutions of constrained non-linear optimization problems arising from kinetic metabolic networks.

In this work, we formulate and study constrained enzyme allocation problems in metabolic networks with general kinetics. In particular, we are interested in enzyme distributions that maximize a designated output flux, given a limited total enzymatic capacity. We show that the optimal distributions of metabolic fluxes differ from solutions obtained by FBA and related stoichiometric methods. Most importantly, we give a rigorous proof for the fact that optimal flux distributions are elementary flux modes. Therein, we make use of results from the theory of oriented matroids that were hitherto only scarcely applied to metabolic networks [3], but offer great potential to unify and advance metabolic network analysis, as mentioned in [9, 20]. Our finding has significant consequences for the understanding of metabolic optimality as well as for the efficient computation of optimal fluxes in kinetic metabolic networks.

The paper is organized as follows. In Section 2, we introduce kinetic metabolic networks and state the enzyme allocation problem of interest. In Section 3, we illustrate our mathematical results and the ideas underlying our proofs by a conceptual example of a minimal metabolic network. In Section 4, we address the connection between metabolic network analysis and the theory of oriented matroids. In particular, we reformulate the optimization problem and show that, if the enzyme allocation problem has an optimal solution, then it has an optimal solution which is an elementary flux mode. Finally, we provide a discussion of our results in the context of metabolic optimization problems.

## 2 Problem statement

After introducing the necessary mathematical notation, we define kinetic metabolic networks and elementary flux modes, and state the metabolic optimization problem that we investigate in the following.

#### Mathematical notation

We denote the positive real numbers by and the non-negative real numbers by . For a finite index set , we write for the real vector space of vectors with , and and for the corresponding subsets. Given , we write if and if . We denote the support of a vector by . For , we denote the component-wise (or Hadamard) product by , that is, .

#### Kinetic metabolic networks

A metabolic network consists of a set of internal metabolites, a set of reactions, and the stoichiometric matrix , which contains the net stoichiometric coefficients for each metabolite in each reaction . The set of reactions is the disjoint union of the sets of reversible and irreversible reactions, and , respectively.

In the following, we assume that each reaction can be catalyzed by an enzyme. Let denote the vector of metabolite concentrations, the vector of enzyme concentrations, and a vector of parameters such as turnover numbers, equilibrium constants, and Michaelis-Menten constants. We write the vector of rate functions as

 v(x;c,p)=c∘κ(x,p)

with a function . In other words, each reaction rate is the product of the corresponding enzyme concentration with a particular kinetics .

A kinetic metabolic network is a metabolic network together with rate functions as defined above. The dynamics of is determined by the ODEs

 dxdt=Nv(x;c,p).

 0=N¯v.

#### Elementary flux modes

A flux mode is a non-zero steady-state flux with non-negative components for all irreversible reactions. In other words, a flux mode is a non-zero element of the flux cone

 C={f∈RR∣Nf=0 and fr≥0 for % all r∈R→}.

An elementary flux mode (EFM) is a flux mode with minimal support:

 f∈C with f≠0 and supp(f)⊆supp(e)⇒supp(f)=supp(e). (efm)

In fact, further implies with . Otherwise, one can construct another flux mode from and with smaller support. As a consequence, there can be only finitely many EFMs (up to multiplication with a positive scalar). For references on elementary flux modes and related computational issues, see [30, 31, 14, 9, 38, 15].

#### Enzyme allocation problem

We are now in a position to state the optimization problem that we study in this work.

Let be a kinetic metabolic network with flux cone . Fix a reaction , positive weights , a subset , and parameters . Maximize the component of the steady-state flux by varying the steady-state metabolite concentrations and the enzyme concentrations . Thereby, fix the weighted sum of enzyme concentrations and require the steady-state flux to be a flux mode:

 max¯x∈X,c∈RR≥¯vr∗ (1a) subject to ∑r∈Rwrcr=ctot, (1b) ¯v∈C,¯vr∗>0. (1c)

Note that the constraint implies which is non-linear, in general.

Further, note that the enzyme allocation problem may be unfeasible, in particular, the constraint (1c) may be unsatisfiable. In this case, the flux cone and the kinetics are “incompatible”. Moreover, even if the problem is feasible, the maximum may not be attained at finite metabolite concentrations. On the other hand, if the problem is feasible, the kinetics is continuous in , and is compact (bounded and closed), then the maximum is attained.

Problem (1) defines a very general metabolic optimization problem, in which the set of enzyme concentrations is adjusted in order to maximize a specific metabolic flux within the steady-state flux vector. In this respect, the weighted sum (1b) may encode different enzymatic constraints, such as limited cellular or membrane surface space, limited nitrogen or transition metal availability, as well as other constraints for the abundance of certain enzymes. In each case, the weight factors denote the fraction of the resource used per unit enzyme. Likewise, the flux component may stand for diverse metabolic processes, ranging from the synthesis rate of a particular product within a specific pathway to the rate of overall cellular growth. The optimization problem seeks to identify the maximum value of , the associated enzyme and steady-state metabolite concentrations, and , as well as the corresponding flux .

## 3 A conceptual example

In order to illustrate our mathematical results and the ideas in its proofs, we first study a minimal metabolic network for the production of a precursor molecule from glucose. In particular, we consider two alternative pathways: fermentation (low yield) and respiration (high yield), cf. Figure 1. The actual metabolic network is further simplified and involves the internal metabolites , the set of reactions consisting of

 1 : Glc\text{ex} ⇌Glc 2 : Glc ⇌2ATP 3 : O2,\text{ex} ⇌O2 4 : Glc+O2 ⇌10ATP 5 : Glc+2ATP →Pre\text{ex},

and the resulting stoichiometric matrix

 Misplaced &

The external substrates/products \ceGlc_ex, \ceO_2,ex, \cePre_ex do not appear in , but their constant concentrations can enter the rate functions as parameters. In a kinetic metabolic network , the rate of reaction is given as

 vr=crκr,

that is, as a product of the corresponding enzyme concentration and a particular kinetics . In this example, we use the kinetics

 κ1 =k1([Glc\text{ex}]−K1[Glc]) κ2 =k2([Glc]−K2[ATP]) κ3 =k3([O2,\text{ex}]−K3[O2]) κ4 =k4([Glc][O2]−K4[ATP]) κ5 =k5[Glc][ATP],

however, our mathematical results do not depend on the kinetics. In vector notation, we write

 v(x;c,p)=c∘κ(x,p),

thereby introducing the concentrations of the internal metabolites and the parameters :

 x =([Glc],[O2],[ATP])T, p =([Glc\text{ex}],[O2,\text{ex}];k1,k2,k3,k4,k5;K1,K2,K3,K4)T.

The dynamics of the network is governed by the ODEs

 dxdt=Nv(x;c,p).

 0=N¯v.

The goal is to maximize the production rate of the precursor \cePre_ex, that is, the component of the steady-state flux , by varying the steady-state metabolite concentrations and the enzyme concentrations . Thereby, the (weighted) sum of enzyme concentrations is fixed and the steady-state flux must be a flux mode:

 max¯x,c¯v5 (2a) subject to 5∑r=1cr=ctot, (2b) N¯v=0,¯v5>0. (2c)

For convenience, we use equal weights in the sum constraint.

To simplify the problem, we consider a restriction on the steady-state metabolite concentrations , in particular, we require . In chemical terms, we assume the thermodynamic feasibility of a situation where all reactions can proceed from left to right. Since , this implies . That is, every component of the steady-state flux is non-negative; in fact, it is zero if and only if the corresponding enzyme concentration is zero.

Every feasible steady-state flux is a flux mode. In particular, where is a flux mode with and . From , we further obtain

 cr=¯v5frκr.

Now, we can rewrite the constraint on the enzyme concentrations as

 ctot=5∑r=1cr=¯v55∑r=1frκr.

Instead of maximizing , we can minimize by varying the steady-state metabolite concentrations and the flux mode . Hence, the enzyme allocation problem (2) with the restriction is equivalent to:

 min¯x,f5∑r=1frκr(¯x,p) (3a) subject to κ(¯x,p)>0,f≥0, (3b) Nf=0,f5=1. (3c)

As shown in Subsection 4.2, every flux mode is a non-negative linear combination of elementary flux modes (EFMs) . In fact, there are two such EFMs,

 e1 =(2,1,0,0,1)T e2 =(65,0,15,15,1)T,

representing fermentation and respiration, and hence

 f=α1e1+α2e2with α1,α2≥0.

Note that we have scaled the EFMs such that . The condition implies . As a result, we obtain another equivalent formulation of the restricted enzyme allocation problem:

 min¯x,α1,α25∑r=1α1e1r+α2e2rκr(¯x,p) (4a) subject to κ(¯x,p)>0, (4b) α1+α2=1. (4c)

We observe that the objective function is linear in and :

 α15∑r=1e1rκr(¯x,p)g1+α25∑r=1e2rκr(¯x,p)g2

with

 g1(¯x,p) =2κ1(¯x,p)+1κ2(¯x,p)+1κ5(¯x,p) g2(¯x,p) =65κ1(¯x,p)+15κ3(¯x,p)+15κ4(¯x,p)+1κ5(¯x,p).

Clearly, and , since . Assume that the minima of and are attained at and , respectively. That is, and . If , then the objective function attains its minimum at , and , that is, for . To see this, assume ; then, for all ,

 α1g1(¯x,p)+α2g2(¯x,p)≥α1^g1+α2^g2>α1^g1+α2^g1=(α1+α2)^g1=^g1.

Conversely, if , the minimum is attained for , and finally, if , both and are optimal. In the degenerate case where at the same minimum point , any (with and ) is optimal.

We can summarize our result as follows: generically, the steady-state flux related to an optimal solution of the restricted enzyme allocation problem (3) is an EFM. The same holds for all appropriate restrictions and hence for the full enzyme allocation problem (2).

For variable external substrate concentrations and , we are interested in which EFM is optimal and when a switch between EFMs occurs. To this end, we determine the optimal solution for each EFM. The optimization problem restricted to EFM is equivalent to

 min¯xg1(¯x,p)

subject to

 κ1(¯x,p)>0,κ2(¯x,p)>0,κ5(¯x,p)>0.

In EFM , reactions 3 and 4 do not carry any flux, that is, . Hence, the corresponding enzyme concentrations are zero, that is, , and there are no constraints involving and . From the optimal metabolite concentrations , we determine the optimal enzyme concentrations , , and as

 cr=¯v5e1rκr=ctote1rκr(¯x,p)∑r=1,2,5e1rκr(¯x,p).

Explicitly, the optimization problem for EFM amounts to

 min[Glc],[ATP]⎛⎜ ⎜⎝2k1([Glc\text{ex}]−K1[Glc])+1k2([Glc]−K2[ATP])+1k5[Glc][ATP]⎞⎟ ⎟⎠

subject to

 [Glc\text{ex}]−K1[Glc]>0,[Glc]−K2[ATP]>0,[Glc][ATP]>0,

where we omit the bar over the steady-state metabolite concentrations. From the optimal metabolite concentrations , we determine the optimal enzyme concentrations , , and . For example,

 c1=ctot2k1([Glc\text{ex}]−K1[Glc])2k1([Glc\text{ex}]−K1[Glc])+1k2([Glc]−K2[ATP])+1k5[Glc][ATP].

The optimization problem restricted to EFM is treated analogously.

Finally, we vary and , solve the restricted optimization problems for EFMs and , and compare the resulting maximum values of , cf. Figure 2. Clearly, the optimal solution of the enzyme allocation problem switches between EFMs and which involves a discontinuous change of enzyme and metabolite concentrations, cf. Figure 3, where we fix and vary .

## 4 Mathematical results

We reformulate the enzyme allocation problem (1) and characterize its solutions. To this end, we employ concepts from the theory of oriented matroids like elementary vectors, sign vectors, and conformal sums.

Realizable oriented matroids arise from vector subspaces. Essentially, a realizable oriented matroid is the set of sign vectors of a subspace or, equivalently, all sign vectors with minimal support. Abstract oriented matroids can be characterized by axiom systems for (co-)vectors (satisfied by the sign vectors of a subspace), (co-)circuits (satisfied by the sign vectors with minimal support), or, equivalently, chirotopes. For an introduction to oriented matroids, we refer to the survey [26], the textbooks [2] and [42, Chapters 6 and 7], and the encyclopedic treatment [6].

In applications to metabolic network analysis, the involved oriented matroids are realizable. For example, the sign vector of a thermodynamically feasible steady-state flux must be orthogonal to all (internal) circuits [3]. In the original proof, the circuit axioms for oriented matroids are used explicitly; however, the result also follows from basic facts about the orthogonality of sign vectors of subspaces [42, Chapter 6]; alternatively, it can be proved using linear programming duality [19, 22]. We note that oriented matroids also appear in the study of directed hypergraph and Petri net models of biochemical reactions [24] and in the theory of chemical reaction networks with generalized mass action kinetics [21].

### 4.1 Elementary vectors

An elementary vector (EV) of a vector subspace is a non-zero vector with minimal support [27]:

 f∈S with f≠0 and supp(f)⊆supp(e)⇒supp(f)=supp(e). (ev)

It is easy to see that EFMs are exactly those EVs of that are flux modes. To our knowledge, this fact has not been clarified before.

###### Lemma 1.

Let be a metabolic network and . The following statements are equivalent:

• is an EFM.

• is an EV of and a flux mode.

###### Proof.

(i) (ii): We have to show that EFM  is an EV of . Suppose ((ev)) is violated, that is, there exists with and . If for all , then in contradiction to ((efm)). Otherwise, consider with the largest scalar such that for all . Then, with and in contradiction to ((efm)). (ii) (i): Let be a flux mode and an EV of . Clearly, ((ev)) implies ((efm)), since implies . ∎

### 4.2 Sign vectors and conformal sums

We define the sign vector of a vector by applying the sign function component-wise. The relations and induce a partial order on : we write for , if the inequality holds component-wise. For , we say that conforms to , if . Analogously, for and , we say that conforms to , if .

The following fundamental result about vectors and EVs will be rephrased for flux modes and EFMs. For a proof, see [27, Theorem 1][2, Proposition 5.35] or [42, Lemma 6.7].

###### Theorem 2.

Let be a subspace. Then every vector is a conformal sum of EVs. That is, there exists a finite set of EVs conforming to such that

 f=∑e∈Ee.

The set can be chosen such that every has a component which is non-zero in , but zero in all other elements of . Hence, and .

It is easy to see that every flux mode is the conformal sum of EFMs. For later use, we present a slightly rephrased version of this result. We note that, if is an EFM, then any element of the ray

 {λe|λ>0}

is an EFM. Hence, we may refer to one representative EFM on each ray.

###### Corollary 3.

Let be a metabolic network, be a sign vector and be a set of representative EFMs conforming to . Then, every flux mode conforming to is a non-negative linear combination of elements of :

 f=∑e∈Eταee with αe≥0.
###### Proof.

Clearly, implies . By Theorem 2, is the conformal sum of EVs of . However, for an EV to conform to , it is required that . Hence, by Lemma 1, is an EFM, which can be written as a positive scalar multiple of a representative EFM. ∎

### 4.3 Problem reformulation

We start with the formal statement of an intuitive argument. Consider a feasible solution of the enzyme allocation problem (1) and the corresponding steady-state flux: if a reaction does not carry any flux, then the corresponding optimal enzyme concentration is zero. We add an appropriate constraint to the enzyme allocation problem and obtain an equivalent optimization problem.

###### Lemma 4.

Let be a kinetic metabolic network. The enzyme allocation problem (1) is equivalent to the following optimization problem:

 max¯x∈X,c∈RR≥¯vr∗ (5a) subject to cr=0if ¯vr=0, (5b) ∑r∈Rwrcr=ctot, (5c) ¯v∈C,¯vr∗>0. (5d)
###### Proof.

For every feasible solution of (1) with objective function , we construct a feasible solution of (5) with objective function : Let . Using , we set

 c′r={λcrif r∈S¯v,0if r∉S¯v.

Clearly,

 ∑r∈Rwrc′r=∑r∈S¯vwrc′r=λ∑r∈S¯vwrcr=ctot.

Further, implies and , that is, . Hence, fulfills constraints (5bcd) and . ∎

As a consequence, variation over enzyme concentrations can be replaced by variation over flux modes.

###### Lemma 5.

Let be a kinetic metabolic network. The enzyme allocation problem (5) is equivalent to the following optimization problem over and :

 min¯x∈X,f∈C∑r∈supp(κ)wrfrκr(¯x,p) (6a) subject to σ(f)≤σ(κ),fr∗=1. (6b)

Let and be corresponding feasible solutions of (5) and (6), respectively. The product of the related objective functions amounts to

 ¯vr∗∑r∈supp(κ)wrfrκr=ctot. (7)
###### Proof.

We show that, for every feasible solution of (5), there exists a feasible solution of (6), and vice versa. Moreover, that the related objective functions fulfill Equation (7).

Assume that is a feasible solution of (5). Define . Clearly, , , and . Hence, is a feasible solution of (6). Let . Using and hence for , we rewrite the sum constraint and obtain the desired Equation (7):

 ctot=∑r∈Rwrcr=∑r∈Sκwrcr=¯vr∗∑r∈Sκwrfrκr.

Conversely, assume that is a feasible solution of (6). Since , we can define by Equation (7), and we set for and for . Clearly, if . Further,

 ∑r∈Rwrcr=∑r∈Sκwrcr=¯vr∗∑r∈Sκwrfrκr=ctot.

By definition, for , and, since , and for . That is, . Hence, is a feasible solution of (5). ∎

We note that the inequality constraints involving the kinetics may be unfeasible. For given flux mode , the existence of steady-state metabolite concentrations such that is equivalent to the existence of chemical potentials such that . Whereas conventional FBA has to be augmented with thermodynamic constraints [3, 20], they are incorporated in the definition of a metabolic network with known kinetics.

### 4.4 Main results

The next statement characterizes optimal solutions of the enzyme allocation problem for fixed metabolite concentrations. Its proof involves the result on conformal sums obtained in Subsection 4.2.

###### Proposition 6.

Let be a kinetic metabolic network. Consider the enzyme allocation problem (1) for fixed . If this restricted optimization problem is feasible, then it has an optimal solution for which the corresponding steady-state flux is an EFM.

###### Proof.

By Lemmas 4 and 5, the enzyme allocation problem (1) is equivalent to optimization problem (6). We consider (6) for fixed and assume that this restricted problem is feasible.

We write short for and introduce and . In (6), we vary over such that and . By Corollary 3, every flux mode conforming to is a non-negative linear combination of elements of , which is a set of representative EFMs conforming to . We assume the EFMs to be scaled by component and divide the set into two subsets, , such that implies and implies . We have:

 f=∑e∈E1αee+∑e∈E0βee with% αe,βe≥0.

From , we obtain the constraint

 1=fr∗=∑e∈E1αeer∗+∑e∈E0βeer∗=∑e∈E1αe.

Using the conformal sum for in (6), we obtain an equivalent formulation of the restricted problem:

 minαe,βe∑r∈Sκwr(∑e∈E1αeer+∑e∈E0βeer)κr (8a) subject to ∑e∈E1αe=1. (8b)

We observe that the objective function is linear in and :

 ∑e∈E1αe∑r∈Sκwrerκrge+∑e∈E0βe∑r∈Sκwrerκrge

with

 ge=∑r∈Sκwrerκrfor e∈Eτ.

Since all conform to , that is, , we have for all and . Moreover, for all , there is such that and hence . Consequently, for all .

Since there is no further restriction on , the minimum of the objective function is attained at for all . In other words, EFMs do not contribute to the optimal solution.

Let be an EFM such that for all . Since and , we have

 ∑e∈E1αege≥∑e∈E1αege′=ge′,

and the minimum of the objective function is attained at and for all other , that is, for . To conclude, we consider a degenerate case: If there are several for which is minimal, then any (with and ) is optimal. ∎

The following statement is the main result of this work.

###### Theorem 7.

Let be a kinetic metabolic network. If the enzyme allocation problem (1) has an optimal solution, then it has an optimal solution for which the corresponding steady-state flux is an EFM.

###### Proof.

Let an optimal solution of the enzyme allocation problem (1) be attained at . Clearly, optimization problem (1) restricted to this particular is feasible. By Proposition 6, this restricted problem has an optimal solution for which the corresponding steady-state flux is an EFM. ∎

In applications, we use Theorem 7 to study the switching behavior of kinetic metabolic networks. Depending on external parameters, the optimal solution of the enzyme allocation problem may switch from one EFM to another, involving a discontinuous change of enzyme and metabolite concentrations. In a first approach, one may vary the external parameters and determine the optimal solution for each EFM in order to find the optimal solution of the full problem. To this end, we transform the optimization problem restricted to an EFM.

###### Corollary 8.

Let be a kinetic metabolic network. In the enzyme allocation problem (1), let the steady-state flux be restricted to , where is an EFM with . Then, this restricted optimization problem is equivalent to the following optimization problem over :

 min¯x∈X∑r∈supp(e)w