Calculation of Equilibrium Compositions of Large Systems of

A more global view of the thermodynamics of a system of biochemical reactions can be obtained by .... systems of biochemical reactions in a living cel...
2 downloads 0 Views 287KB Size
J. Phys. Chem. B 2000, 104, 4807-4814

4807

Calculation of Equilibrium Compositions of Large Systems of Biochemical Reactions Robert A. Alberty Department of Chemistry, Massachusetts Institute of Technology, Cambridge, Massachusetts 02139 ReceiVed: January 13, 2000; In Final Form: February 28, 2000

A more global view of the thermodynamics of a system of biochemical reactions can be obtained by specifying concentrations of coenzymes such as adenosine triphosphate and adenosine diphosphate. These calculations are facilitated by use of conservation matrices that can be arranged in such a way that the coenzymes are components. In calculations at specified pH in dilute aqueous solutions, neither hydrogen atoms or oxygen atoms are conserved, and so the applicable conservation matrices are calculated from stoichiometric number matrices that omit the stoichiometric numbers for H2O. For a system of biochemical reactions, the number of components is generally greater than the number of elements other than hydrogen and oxygen because of constraints that arise in the mechanisms of the enzyme-catalyzed reactions. When concentrations of coenzymes are specified, the system can be described as being made up of pseudoisomer groups, each made up of one or several reactants. The values of apparent equilibrium constants K′′ for reactions between the pseudoisomer groups can be calculated if the standard transformed Gibbs energies of formation of the reactants are known at the desired pH or the apparent equilibrium constants K′ are all known. These concepts and equations are applied to glycolysis plus the citric acid cycle, and it is shown that quantitative thermodynamic calculations can be made by working with 7 reactants, rather than 32, and 5 reactions, rather than 20.

Introduction The calculation of the equilibrium composition of a large system of biochemical reactions in dilute aqueous solutions is complicated by the large number of reactants. When the pH is specified, the reaction system can be described in terms of sums of species, like the sum of species that make up citrate, but there is a need for a still more global way to make thermodynamic calculations on large reaction systems. An additional problem in treating the thermodynamics of such systems is dealt with here, and that is the fact that oxygen atoms are not conserved in reactions involving H2O because the solvent provides a reservoir that is assumed to be infinite. When the pH is held constant, the Gibbs energy G no longer provides the criterion for spontaneous change and equilibrium, and it is necessary to define a transformed Gibbs energy G′ using a Legendre transform to obtain a thermodynamic potential that does provide the criterion for spontaneous change and equilibrium.1,2 At a specified pH, a reactant, such as citrate, has a standard transformed Gibbs energy of formation ∆fG′o that can be calculated from the standard Gibbs energies of the species it includes. At a specified pH there is no conservation equation for hydrogen atoms in the reaction system, and no conservation equation for oxygen atoms. But there are conservation equations for other elements. (This statement excludes calcium ions, for example, that are bound reversibly and can be treated in the same way as hydrogen ions.) A feature that distinguishes many enzyme-catalyzed reactions from chemical reactions is that they may involve constraints in addition to element balances. For example, when two chemical reactions without a common reactant are coupled together by an enzymatic mechanism, the number of reactions in the system is decreased by one, but the number of reactants is unchanged. A system of

biochemical reactions may involve a large number of independent constraints of this typesthat is, constraints that are represented by linear conservation equations. Both element balances and constraints of this second type can be represented by conservation equations written in terms of amounts of components that are selected from the reactants in the system. The term component is used in a special way when reactions occur in a system. Components are the things that all of the species in a system can be considered to be made up of,3,4 and they are conserved. In a chemical reaction system the components can be chosen to be atoms of elements and electric charges. But, components can also be taken to be arbitrary collections of atoms, provided the conservation equations for the particular collections of atoms are independent. The number C of components is fixed, but many choices of components can be made. The operational method for testing a set of components for independence, or obtaining a suitable set of components, is to write the conservation matrix (see below) for the reaction system and carry out a row reduction.5,6 When the columns in the conservation matrix are arranged in such an order that an identity matrix is obtained on the left side of the conservation matrix, the species listed in the first group of columns are components, and the columns for noncomponents show the compositions of the noncomponents in terms of components. In considering biochemical reactions at specified pH, the row for hydrogen atoms is deleted. When this is done, the columns for some species become identical, which means that they are pseudoisomers. A pseudoisomer group is a group of species in which the mole fractions within the group are independent of the other reactions that occur in the system. Pseudoisomer groups can be treated as reactants; that is, reactions in a system can be written in terms of pseudoisomer groups, as is done in

10.1021/jp000185w CCC: $19.00 © 2000 American Chemical Society Published on Web 04/18/2000

4808 J. Phys. Chem. B, Vol. 104, No. 19, 2000

Alberty

considering biochemical reactions at specified pH. This is referred to as a level 2 treatment in terms of reactants to distinguish it from the level 1 treatment in terms of species. This process can be taken a step further by holding the concentrations of some reactants constant in a level 3 treatment of a system of biochemical reactions in terms of pseudoisomer groups of reactants.7,8 This article deals with the advantages of selecting coenzymes, like adenosine triphosphate (ATP) and adenosine diphosphate (ADP), as components in considering the thermodynamics of systems of biochemical reactions in a living cell where the concentrations of coenzymes are in steady states. The term coenzyme is not very specific, but it is used here for reactants that are not the principal carbon compounds being oxidized in metabolism. When the concentrations of coenzymes are constant, their conservation equations are eliminated, and as a consequence the thermodynamics of a system of biochemical reactions can be discussed in terms of pseudoisomer groups of reactants that are involved in fewer reactions. Pseudoisomer groups have standard further transformed Gibbs energies of formation ∆fG′′o that can be calculated from the standard transformed Gibbs energies of formation ∆fG′o of the reactants they include.7,8 Thus, a more global view of the thermodynamics of the system of biochemical reactions can be obtained at specified concentrations of coenzymes. A choice can be made as to how many and which coenzymes are held constant. The equilibrium composition of the system can be calculated in terms of pseudoisomer groups. The equilibrium concentrations of these pseudoisomer groups can then be used to calculate the equilibrium concentrations of the reactants that are not in steady states. Thus, no thermodynamic information is lost in taking this more global view. For a system of chemical reactions, the conservation matrix A has a row for each element and a column for each species. Electric charge may be a component, but electric charge and some elements that are present may not be components because the conservation equations for components must be independent. Various choices of components are possible, and this choice is important because amounts of components are conserved and amounts of noncomponents are not. The corresponding stoichiometric number matrix ν has a row for each species and a column for each independent reaction. The number of components is represented by C and the number of independent reactions is represented by R, and so the number N of species is given by N ) C + R. Matrix A is C × N, and matrix ν is N × R. These two types of matrices are related by

Aν ) 0

(1)

so that a possible stoichiometric number matrix ν can be obtained by calculating6 a basis for the null space of A. Inversely, a basis for the transpose AT of a conservation matrix A can be calculated from the transpose νT of the stoichiometric number matrix because

νTAT ) 0

(2)

Neither A nor ν is unique. When the pH is held constant, the conservation matrix for the reaction system is altered because hydrogen atoms are no longer conserved. When the row for the hydrogen component in A is deleted, groups of species become pseudoisomers because they have the same numbers of atoms other than hydrogen atoms. Redundant columns are deleted and the columns in the apparent conservation matrix A′ are now labeled

with names of reactants (pseudoisomer groups). The apparent conservation matrix A′ is C′ × N′, where C′ is the number of apparent components and N′ is the number of reactants. When the pH is specified, C′ ) C - 1. The corresponding stoichiometric number matrix ν′ calculated from A′ by use of an equation like eq 1 is N′ × R′, where R′ is the number of independent biochemical reactions written in terms of reactants. The number N′ of reactants is given by N′ ) C′ + R′. Conversely, the apparent conservation matrix A′ can be obtained by calculating the null space of (ν′)T by use of an equation like eq 2. As pointed out in the third paragraph of the Introduction, the remarkable feature of biochemical reaction systems is that they usually have more components than the number of elements involved because the mechanisms of enzyme-catalyzed reactions may introduce additional conservation equations. The approach to the thermodynamics of biochemical reactions that is described here is based on the fact that conservation equations can be expressed in terms of reactants, rather than elements and constraints arising from enzymatic mechanisms. When the concentrations of reactants that have been chosen as components are held constant, the conservation equations for these components are eliminated. When there are fewer conservation equations, a reaction system can be described in a simpler way. When the concentrations of coenzymes such as ATP and ADP are held constant in addition to the pH, the conservation matrix A′ needs to be arranged so that the coenzymes are components. This can be done by arranging the columns in a particular order and making a row reduction (Gaussian elimination). When the concentrations of coenzymes are held constant, their rows and columns are deleted from A′, and consequently, groups of reactants become pseudoisomers because they have the same compositions in terms of the other components. Redundant columns for pseudoisomers are deleted and the columns in the new apparent conservation matrix A′′ now apply to pseudoisomer groups. Thus, the A′′ matrix is C′′ × N′′. The corresponding stoichiometric number matrix ν′′ calculated from A′′ by use of an equation like eq 1 is N′′ × R′′, where R′′ is the number of independent reactions written in terms of pseudoisomer groups. The number N′′ of pseudoisomer groups is given by N′′ ) C′′ + R′′. It might appear that decreasing the number of components would increase the number of reactions required to describe the system, but the number N′′ of pseudoisomer groups is smaller than the number of reactants. As we have seen before, the apparent conservation matrix A′′ can be obtained by calculating the null space of the (ν′′)T by use of an equation like eq 2. When the concentrations of coenzymes that have been chosen as components are held constant, the conservation equations for these components are eliminated. When there are fewer conservation equations, a reaction system can be described in a simpler way. The maximum number of components that can be held constant at specified pH is equal to C′′ - 1 because at least one component must remain. The complete Legendre transform for a system involves all of the components and yields the Gibbs-Duhem equation that gives the relation between the intensive variables for the system. Since the concentration of H2O is not included in the equilibrium expressions for reactions involving H2O in dilute aqueous solutions, eqs 1 and 2 for chemical reactions and their counterparts at constant pH and constant [ATP] and [ADP], for example, are not applicable if the stoichiometric numbers are taken from balanced biochemical equations involving H2O. The solution to this problem that is used here is to always start with

Equilibrium Compositions of Biochemical Reactions

J. Phys. Chem. B, Vol. 104, No. 19, 2000 4809

the stoichiometric number matrix omitting the stoichiometric numbers of H2O and to calculate the corresponding conservation matrix using eq 2. When this is done, oxygen atoms are not conserved. However, it is still necessary to use ∆fG′o(H2O) in calculating the apparent equilibrium constant K′ for a biochemical reaction at a specified pH. Thermodynamics of Pseudoisomer Groups When the pH is held constant, a Legendre transform is used to define a transformed Gibbs energy G′ that provides the criterion for spontaneous change and equilibrium. This Legendre transform is1,2

G′ ) G - nc(H)µ(H+)

(3)

where nc(H) is the amount of the hydrogen component in the system and µ(H+) is the specified chemical potential of hydrogen ions. A Legendre transform always involves subtraction of a product of conjugate variables.9 There are corresponding transformed entropies S′ and transformed enthalpies H′. It can be shown that the standard transformed Gibbs energy of formation ∆fGi′o of a species at a specified pH can be calculated from its standard Gibbs energy of formation ∆fGio by using

∆fGi′° ) ∆fGi° - NH(i)(∆fG°(H+) + RT ln 10-pH) (4) where NH(i) is the number of hydrogen atoms in species i and pH ) -log[H+]. Since the equilibria in biochemical reaction systems are generally affected by the ionic strength, thermodynamic properties such as ∆fGio and ∆fGi′o are taken to be functions of the ionic strength.10 Therefore, the expressions for equilibrium constants K′ in terms of species and apparent equilibrium constants K′ in terms of reactants (sums of species) are written in terms of concentrations. K is a function of T, P, and ionic strength, and K′ is a function of T, P, pH, and ionic strength. When a biochemical reactant is made up of several species that are pseudoisomers at a specified pH, the standard transformed Gibbs energy of formation ∆fG′o(iso) of the pseudoisomer group at equilibrium can be calculated from the standard transformed Gibbs energies of formation (see eq 4) of the several species by use of Niso

∆fG′°(iso) ) -RT ln

exp(-∆fGi′°/RT) ∑ i)1

(5)

where Niso is the number of pseudoisomers in the group. The equilibrium mole fractions ri within the pseudoisomer group at a specified pH are given by

ri ) exp((∆fG′°(iso) - ∆fGi′°)/RT)

N′

νi′∆fGi′° ) -RT ln K′ ∑ i)1

Cs

G′′ ) G′ -

(7)

where the prime is used on the stoichiometric numbers to distinguish them from the stoichiometric numbers of the underlying reactions in terms of species. Here N′ is the number

nc(i)µ′(i) ∑ i)1

(8)

where nc(i) is the amount of component i in the system and Cs is the number of components with specified concentrations. There are corresponding further transformed entropies S′′ and further transformed enthalpies H′′. It can be shown that the standard further transformed Gibbs energy of formation ∆fGi′′o of a reactant at specified pH and concentrations of certain components (coenzymes) can be calculated from its standard transformed Gibbs energy of formation ∆fGi′o at specified pH by using Cs

∆fGi′′° ) ∆fGi′° -

Ncj(i)(∆fG′°(j) + RT ln[j]) ∑ i)1

(9)

where Ncj(i) is the number of component j involved in reactant i. When these reactants are pseudoisomers, the standard further transformed Gibbs energy of formation ∆fG′′o(iso) of the pseudoisomer group can be calculated by using Niso

∆fG′′o(iso) ) -RT ln

exp(-∆fGi′′o/RT) ∑ i)1

(10)

where Niso is the number of pseudoisomers in the group. The equilibrium mole fractions ri within the pseudoisomer group at specified concentrations of coenzymes can be calculated by using

ri ) exp((∆fG′′°(iso) - ∆fG′′°(i))/RT)

(11)

If the standard further transformed Gibbs energies of formation of all the pseudoisomer groups in a reaction are known, the apparent equilibrium constant K′′ at the specified concentrations of coenzymes such as ATP and ADP can be calculated using

(6)

If the standard transformed Gibbs energies of formation of all the reactants in a biochemical reaction are known, the apparent equilibrium constant K′ for the reaction at a specified pH can be calculated using

∆rG′° )

of reactants in a reaction, but more generally N′ is used for the number of reactants (pseudoisomer groups) in a system of reactions. The expression for the apparent equilibrium constant for a reaction involving H2O in dilute aqueous solutions does not have a term in the concentration of H2O, but ∆fG′o(H2O) is involved in the summation in eq 7. When concentrations of coenzymes such as ATP and ADP are held constant, a Legendre transform is used to define a further transformed Gibbs energy G′′ that provides the criterion for spontaneous change and equilibrium.7,8 The concentrations of a number of reactants can be specified, but these reactants must be components, rather than noncomponents. The Legendre transform that is used to define the further transformed Gibbs energy G′′ is

N′′

∆rG′′° )

νi′′∆fGi′′° ) -RT ln K′′ ∑ i)1

(12)

where the double prime is used on the stoichiometric numbers to distinguish them from the stoichiometric numbers of the underlying reactions in terms of reactants. Here N′′ is the number of reactants in a reaction, but more generally N′′ is used for the number of reactants (pseudoisomer groups) in a system of reactions. When H2O is a reactant, its ∆fG′′o value has to be included in eq 12, but ∆fG′′o(H2O) ) ∆fG′o(H2O) because H2O is not made up of other components. Equilibrium compositions can be calculated at level 3 using apparent equilibrium constants

4810 J. Phys. Chem. B, Vol. 104, No. 19, 2000

Alberty

K′′, and further thermodynamic calculations, such as the effect of temperature, can be made. To show how these concepts and equations apply to actual systems, glycolysis is treated first, the citric acid cycle is treated second, and then glycolysis and the citric acid cycle are treated together. Glycolysis Glycolysis provides a pathway between glucose and pyruvate. It involves 10 biochemical reactions and 16 reactants. Water is not counted as a reactant in writing the stoichiometric number matrix or the conservation matrix for reasons described above. Thus, there are six components because C′ ) N′ - R′ ) 16 10 ) 6. From a chemical standpoint this is a surprise because the reactants involve only C, H, O, N, and P. Since H and O are not conserved at specified pH in dilute aqueous solution, there are only three conservation equations based on elements. Thus, three additional conservation relations arise from the mechanisms of the enzyme-catalyzed reactions in glycolysis. Some of these are discussed in ref 1. At specified pH in dilute aqueous solutions the reactions in glycolysis are represented by

Glc + ATP ) G6P + ADP

(13)

G6P ) F6P

(14)

F6P + ATP ) FBP + ADP

(15)

FBP ) DHAP + GAP

(16)

DHAP ) GAP

(17)

GAP + Pi + NADox ) 13BPG + NADred

(18)

13BPG + ADP ) 3PG + ATP

(19)

3PG ) 2PG

(20)

2PG ) PEP (+ H2O)

(21)

PEP + ADP ) Pyr + ATP

(22)

where the following abbreviations are used: glucose (Glc), adenosine triphosphate (ATP), adenosine diphosphate (ADP), glucose 6-phosphate (G6P), fructose 6-phosphate (F6P), fructose 1,6-biphosphate (FBP), glyceraldehyde 3-phosphate (GAP), dihydroxyacetone phosphate (DHAP), 1,3-bisphosphoglycerate (13BPG), nicotinamide adenine dinucleotide-oxidized (NADox), nicotinamide adenine dinucleotide-reduced (NADred), 3-phosphoglycerate (3PG), 2-phosphoglycerate (2PG), phosphoenolpyruvate (PEP), and pyruvate (Pyr). If reactions 18-22 are each multiplied by 2 and the reactions are added, the net reaction is

Glc + 2Pi + 2ADP + 2NADox ) 2Pyr + 2ATP + 2NADred (+ 2H2O) (23) When a computer is used, a net reaction is obtained more conveniently by use of a matrix multiplication (see Appendices A and B). H2O is put in parentheses because its stoichiometric number is not used in the stoichiometric number matrix, but it is involved in the calculation of K′ using eq 7. In writing the stoichiometric number matrix for glycolysis, there is a choice as to the order of the reactants. To make Glc,

Figure 1. Apparent stoichiometric number matrix ν′ for the 10 reactions of glycolysis at specified pH in dilute aqueous solutions.

ATP, ADP, NADox, NADred, and Pi components, they are put first in the columns for the transposed stoichiometric matrix, followed by the rest of the reactants ending with Pyr. Other choices of components can be made, but not any arbitrary choice can be made because the conservation equations for components have to be independent. The stoichiometric number matrix for glycolysis is shown as (ν′)T in Figure 1. To check that these 10 reactions are indeed independent, a row reduction can be used. Mathematica is very useful11 for performing matrix operations, and RowReduce is used for that purpose. Since it is easy to make errors in setting up a stoichiometric number matrix, it is useful to check the matrix by using it to print out the reactions (see Appendix B). Conservation matrix A′ that corresponds to this stoichiometric matrix is obtained by calculating the null space of (ν′)T, as indicated by eq 2. To obtain a conservation matrix with identifiable rows, RowReduce is used again and the result is shown in Figure 2. This shows that Glu, ATP, ADP, NADox, NADred, and Pi can be taken as the six components for glycolysis. This is not the only possible choice for components, but it is the one to use to obtain the most global view of the thermodynamics of the reaction system. The conservation matrix in Figure 2 illustrates the statement that conservation equations can be expressed in terms of reactants, rather than elements and constraints arising from enzyme mechanisms. When the concentration of a component is held constant in an equilibrium calculation, its row and column in the conservation matrix A′ are deleted. When the rows and columns for ATP, ADP, NADox, NADred, and Pi are deleted, the remaining stoichiometric number matrix is {{1,1,1,1,1/2,1/2,1/2,1/2,1/2,1/2,1/ 2}}, which applies to the conservation of the component Glu. Under these conditions, the reactants in glycolysis consist of two pseudoisomer groups. The first contains Glc, G6P, F6P, and FBP and will be referred to as C6, where the C refers to the element carbon. The second contains DHAP, 13BPG, 3PG, 2PG, PEP, GAP, and Pyr and will be referred to as C3. Deleting the redundant columns in the conservation matrix for glycolysis yields the apparent conservation matrix A′′ ) {{1,1/2}}, where the 1 is for C6 and the 1/2 is for C3. This is the conservation matrix that applies when [ATP], [ADP], [NADox], [NADred], and [Pi] are held constant. Calculating the null space for this conservation matrix yields (ν′′)T ) {{-1/2,1}}, which indicates that there is a single reaction (1/2) C6 ) C3 or C6 ) 2C3 at level 3. Thus, specifying the concentrations of 5 coenzymes has reduced the system of 10 reactions and 16 reactants to one

Equilibrium Compositions of Biochemical Reactions

J. Phys. Chem. B, Vol. 104, No. 19, 2000 4811

Figure 2. Row-reduced form of the apparent conservation matrix A′ for glycolysis at a specified pH in dilute aqueous solutions.

reaction and two reactants. The expression for the apparent equilibrium constant KGLY′′ under these conditions can be written

KGLY′′ ) [C3] /[C6] 2

(24)

where KGLY′′ is a function of T, P, pH, [ATP], [ADP], [NADox], [NADred], [Pi], and ionic strength. Equilibrium constants are dimensionless, but the reference concentration c° ) 1 M is omitted in the denominator of eq 24 as a simplification. This equilibrium expression provides the most global view of the thermodynamics of glycolysis. The value of the apparent equilibrium constant KGLY′′ can be calculated when ∆fGi′o values are known for the 16 reactants at the desired T, P, pH, and ionic strength. These values can be used in eq 10 to calculate the ∆fGi′′o values of the reactants at the desired concentrations of the coenzymes. Figure 2 shows the content of each of these reactants in terms of components; in other words, the columns give the values of Ncj(i) for each of the components held constant. Equation 10 can then be used to calculate the ∆fG′′o(iso) values of the C6 and C3 pseudoisomer groups. To calculate the numerical value for KGLY′′ using eq 12, the reaction between pseudoisomer groups should be written C6 ) 2C3 (+ 2H2O) because ∆fG′′o(H2O) is involved in calculating KGLY′′. Note that ∆fG′′o(H2O) ) ∆fG′o(H2O) because H2O is not made up of components. The value of KGLY′′ is not calculated here because some of the needed ∆fG′o values are currently unknown. The equilibrium concentrations of C6 and C3 can be calculated using eq 24, and then the concentrations of C6 and C3 can be divided into the equilibrium concentrations of each of the reactants (noncomponents) by use of eq 11. Thus, KGLC′′ provides the means to calculate the equilibrium concentrations of the 11 reactants for which concentrations have not been specified. This more global view of the thermodynamics of a system of biochemical reactions provides different information than the net reaction (eq 23) for the system because it deals with the pseudoisomer groups C6 and C3. The expression for the apparent equilibrium constant for the net reaction for glycolysis given in eq 23 can be used to calculate the equilibrium value of [Pyr]2/ [Glu] by setting the concentrations of the coenzymes equal to their steady-state values. The apparent equilibrium constant for the net reaction at 298.15 K, pH 7, and 0.25 M ionic strength calculated from ∆fG′o values12,13 is 1.41 × 1014. When the steady-state concentrations of ATP and NADred are 0.01 M, of ADP and NADox are 10-5 M, and of Pi is 0.001 M, the equilibrium concentration of glucose will exceed that of pyruvate. The advantage of the equilibrium expression in eq 24 is that it yields the equilibrium concentrations of pseudoisomer groups C6 and C3; that is, it accounts for all of the reactants in glycolysis. Citric Acid Cycle At a specified pH in dilute aqueous solution, the biochemical reactions in the citric acid cycle can be written as

Pyr + CoA + NADox (+ H2O) ) totCO2 + NADred + AcetCoA (25) AcetCoA + Oxalac (+ H2O) ) CoA + Cit

(26)

Cit ) cisAcon (+ H2O)

(27)

cisAcon (+ H2O) ) isoCit

(28)

isoCit + NADox (+ H2O) ) NADred + Ketoglu + totCO2 (29) Ketoglut + CoA + NADox (+ H2O) ) SuccCoA + totCO2 + NADred (30) SuccCoA + GDP + Pi ) GTP + CoA + Succ (31) Succ + FADox ) Fum + FADred

(32)

Fum (+ H2O) ) Mal

(33)

Mal + NADox ) Oxaloac + NADred

(34)

The abbreviations in addition to those used earlier are coenzyme A (CoA), total carbon dioxide in solution (totCO2), acetyl-CoA (AcetCoA), oxaloacetate (Oxalac), citrate (Cit), cis-aconitate (cisAcon), isocitrate (isoCit), R-ketoglutarate (Ketoglu), succinyl-CoA (SuccCoA), guanosine triphosphate (GTP), guanosine diphosphate (GDP), succinate (Succ), flavin adenine dinucleotide-oxidized (FADox), flavin adenine dinucleotide-reduced (FADred), fumarate (Fum), and L-malate (Mal). The net reaction obtained by adding these 10 reactions is

Pyr + 4NADox + GDP + Pi + FADox (+ 5H2O) ) 3totCO2 + 4NADred + FADred + GTP (35) Rather than writing these reactions in terms of gaseous CO2, they are written in terms of the total concentration of carbon dioxide in aqueous solution at specified, T, P, pH, and ionic strength.14,15 The four species making up totCO2 are CO2(aq), H2CO3, HCO3-, and CO32-. When a biochemical reaction is balanced with CO2(g) and it is desired to balance it with totCO2, a H2O has to be inserted on the other side of the biochemical equation. The advantage of writing biochemical reactions in terms of totCO2, rather than CO2(g), is that the expression for the apparent equilibrium constant K′ is then written in terms of concentrations of all of the reactants in the aqueous solution, except H2O. In writing the stoichiometric number matrix, Pyr is put first, followed by coenzymes and then noncomponents, ending with totCO2. The stoichiometric number matrix ν′ for the citric acid cycle is shown in Figure 3. This matrix can be used to calculate the apparent conservation matrix A′ for the citric acid cycle at specified pH in dilute aqueous solutions.

4812 J. Phys. Chem. B, Vol. 104, No. 19, 2000

Alberty yields (ν′)T ) {{1,0,0,0,0,-3},{0,1,0,0,-1,1},{0,0,1,0,0,-2}, {0,0,0,1,-1,-1}}. This stoichiometric number matrix has the same row-reduced form as that obtained from the following four reactions:

Pyr (+ H2O) ) AcetCoA + totCO2

(36)

AcetCoA + C4 (+ 2H2O) ) C6′

(37)

C6′ (+ H2O) ) totCO2 + Ketoglu

(38)

Ketoglu (+ H2O) ) totCO2 + C4

(39)

Specifying the concentrations of 8 coenzymes has reduced the system of 10 reactions to 5 and 20 reactants to 7 but still allows quantitative calculations. The net reaction for the citric acid cycle under these conditions is given by the sum of reactions 36-39: Figure 3. Stoichiometric number matrix ν′ for the 10 reactions of the citric acid cycle at specified pH in dilute aqueous solutions.

Pyr (+ 5H2O) ) 3 totCO2

(40)

The apparent equilibrium constant for this reaction The conservation matrix A′ obtained from the null space of (ν′)T is shown in Figure 4 after it has been row reduced. This shows that there are 10 components, as expected from the fact that there are 20 reactants and 10 independent reactions: C′ ) N′ - R′ ) 20 - 10 ) 10. The 10 components are Pyr, GTP, GDP, NADox, NADred, FADox, FADred, Pi, CoA, and SuccCoA. This analysis reveals an interesting difference between glycolysis and the citric acid cycle. In glycolysis all of the coenzymes are represented by terms in the expression for the apparent equilibrium constant for the net reaction. In the citric acid cycle the components CoA and SuccCoA do not appear in the expression for the apparent equilibrium constant for the net reaction. Therefore, the equilibrium value of [totCO2]3/[Pyr] at specified concentrations of GTP, GDP, NADox, NADred, FADox, FADred, and Pi is independent of the equilibrium concentrations of [CoA] and [SuccCoA]. However, if [CoA] or [SuccCoA] is specified, the concentrations of reactants other than Pyr and totCO2 will be affected. As mentioned at the end of the section on glycolysis, the equilibrium concentrations of all the noncomponents can be calculated when various sets of components have specified concentrations, provided that the number of specified components is less than C′ - 1. In considering the citric acid cycle, it may be of interest to calculate the equilibrium compositions for various concentrations of CoA and to leave the equilibrium concentration of SuccCoA to be calculated. This is done by specifying the concentrations of GTP, GDP, NADox, NADred, FADox, FADred, Pi, and CoA. When this is done, the conservation matrix becomes A′′ ) {{1,0,0,2/3,2/3,2/3,2/3,1/3,0,0,0,1/3}, {0,1,1,0,1,1,1,1,1,1,1,0}}, where the columns are for Pyr, SuccCoA, Oxaloac, AcetCoA, Cit, cisAcon, isoCit, Ketoglu, Succ, Mal, Fum, and totCO2. Several of these columns are redundant and can be deleted. This reduces the conservation matrix to {{1,0,2/3,2/3,1/3,1/3},{0,1,0,1,1,0}}, where the columns are identified as Pyr, C4 (SuccCoA, Oxaloac, Succ, Mal, Fum), AcetCoA, C6′ (Cit, cisAcon, isoCit), Ketoglu, and totCO2. In other words, the conservation matrix A′′ has six columns, four of which represent single reactants (Pyr, AcetCoA, Ketoglu, and totCO2), and two of which are pseudoisomer groups C4 and C6′. The prime is put on the C6′ to distinguish it from the C6 pseudoisomer group in glycolysis. The corresponding stoichiometric number matrix can be obtained by calculating the null space of this conservation matrix and row reducing it. This

K′′ ) [totCO2]3/[Pyr]

(41)

is a function of T, P, pH, [GTP], [GDP], [NADox], [NADred], [FADox], [FADred], [Pi], [CoA], and ionic strength. It can only be calculated when the ∆fG′o values are known for all of the reactants in the citric acid cycle or all the apparent equilibrium constants K′ are known. The apparent equilibrium constants K′′ for reactions 36, 37, 38, and 39 can also be calculated. Then the equilibrium composition can be calculated in terms of Pyr, C4, AcetCoA, C6, Ketoglu, and totCO2. The distribution of reactants within the two pseudoisomer groups can be calculated using eq 11. The effects of varying the specified concentrations of GTP, GDP, NADox, NADred, FADox, FADred, Pi, and CoA can be studied. Glycolysis Plus the Citric Acid Cycle These two reaction systems provide a pathway between glucose and carbon dioxide. The first 10 reactions are (13)(22) and the last 10 reactions are (25)-(34). This reaction system involves 32 different reactants since glycolysis involves 16 reactants, the citric acid cycle involves 20, and there is an overlap of 4 reactants, namely Pyr, NADox, NADred, and Pi. Again H2O is not counted as a reactant. The net reaction for glycolysis plus the citric acid cycle is

Glc + 10NADox + 2ADP + 2GDP + 2FADox + 4Pi (+ 8H2O) ) 6totCO2 + 10NADred + 2ATP + 2GTP + 2FADred (42) The stoichiometric number matrix is set up with the reactants previously selected as components in glycolysis being taken first, the reactants previously selected as components in the citric acid cycle being taken second, then the noncomponents of glycolysis, and the noncomponents of the citric acid cycle, with totCO2 last. It does not seem necessary to give this 32 × 20 stoichiometric matrix ν′ since it is made up of the matrices in Figures 1 and 3 put together properly. The corresponding conservation matrix A′ is obtained by calculating the null space of (ν′)T. Since this conservation matrix is 12 × 32, there are 12 components, as expected from C′ ) N′ - R′ ) 32 - 20 ) 12. Another way to look at these 12 components is that glycolysis

Equilibrium Compositions of Biochemical Reactions

J. Phys. Chem. B, Vol. 104, No. 19, 2000 4813

Figure 4. Row-reduced conservation matrix A′ for the citric acid cycle at specified pH in dilute aqueous solutions.

Figure 5. Row-reduced conservation matrix A′ to the right of the 12 × 12 identify matrix for glycosis plus the citric acid cycle.

Figure 6. Stoichiometric number matrix (ν)T for glycolysis plus the citric acid cycle at specified T, P, pH, and specified concentrations of 10 components.

has 6 components, the citric acid cycle has 10, and there is an overlap of 4 components (NADox, NADred, Pi, and carbon). This can be contrasted with the fact that only C, N, P and S atoms are conserved. The other 8 conservation relations arise in the mechanisms of the enzyme-catalyzed reactions. The 12 components Glu, ATP, ADP, NADox, NADred, Pi, GTP, GDP, FADox, FADred, CoA, and SuccCoA are independent; that is, none of them can be made from the other 11. The row-reduced conservation matrix A′ has a 12 × 12 identity matrix in the left side, and since that is known, only the remaining 12 × 20 part of the reduced conservation matrix is shown in Figure 5. The columns for the noncomponents in glycolysis give the same information about how to calculate ∆fG′′o of these reactants as Figures 2 and 4. The columns for the noncomponents in the citric acid cycle are different from the corresponding columns in the section on the citric acid cycle because they can be made up from reactants in glycolysis like ATP, which were not available when the citric acid cycle was treated by itself. When the concentrations of the 10 components other than Glu and SuccCoA are held constant, the rows and columns for these 10 components in A′ are deleted. When the redundant columns are also deleted, the apparent conservation matrix A′′ is {{1,0,1/2,1/2,1/3,1/6,1/6},{0,1,0,0,1,1,0}}. Three of these columns refer to single reactants (AcetCoA, Ketoglu, and totCO2) and four of the columns refer to pseudoisomer groups, which are referred to as C6, C4, C3, and C6′. The pseudoisomer group Glc, G6P, F6P, and FBP is referred to as C6. The pseudoisomer group SuccCoA, Oxaloac, Succ, Mal, and Fum is referred to as C4.

The pseudoisomer group DHAP, 13BPG, 3PG, 2PG, PEP, and Pyr is referred to as C3. The pseudoisomer group Cit, cisAcon, and isoCit is referred to as C6′. The corresponding stoichiometric number matrix ν′′ can be obtained by calculating the null space of A′′, and its row-reduced form is given in Figure 6. This stoichiometric number matrix is consistent with the following five reactions:

C6 ) 2C3 (+ H2O)

(43)

C3 (+ H2O) ) AcetCoA + totCO2

(44)

AcetCoA + C4 (+ H2O) ) C6′

(45)

C6′ (+ H2O) ) totCO2 + Ketoglu

(46)

Ketoglu (+ H2O) ) totCO2 + C4

(47)

When the last four reactions are multiplied by 2 and the reactions are added, the following net reaction is obtained for glycolysis plus the citric acid cycle when the concentrations of ATP, ADP, NADox, NADred, GTP, GDP, FADox, FADred, Pi, and CoA are specified:

C6 (+ 8H2O) ) 6totCO2

(48)

The standard further transformed Gibbs energies of formation of the seven reactants in reactions 43-47 can be calculated by using eqs 9 and 10. These standard further transformed Gibbs energies of formation can then be used to calculate the equilibrium composition in terms of these seven reactants for any desired steady-state concentrations of ATP, ADP, NADox, NADred, GTP, GDP, FADox, FADred, Pi, and CoA. The equilibrium concentrations of pseudoisomer groups can then be divided into the concentrations of the reactants in each pseudoisomer group by use of eq 11. Thus, quantitative calculations can be made on glycolysis plus the citric acid cycle by working with 7 reactants, rather than 32, and 5 reactions, rather than 20.

4814 J. Phys. Chem. B, Vol. 104, No. 19, 2000

Alberty

A still more global view of glycolysis plus the citric acid cycle can be obtained by specifying the concentration of SuccCoA in addition.

number matrix for the system and the path vector:16

Discussion

The path vector gives the number of times each reaction in the system has to go to yield the net reaction.

(ν′)T

When the transpose of the stoichiometric number matrix for a system of biochemical reactions is used to calculate the corresponding conservation matrix A′, the conservation matrix is not unique because calculating the null space yields a basis for the stoichiometric number matrix. One of the simplest ways to change the form of the conservation matrix without changing the information in it is to make a row reduction (Gaussian elimination). If the conservation equations for the reactants that are listed first are independent, there will be an identity matrix in the left side of the conservation matrix. In this case, the rows of the conservation matrix give the coefficients in the conservation equations for components that are identical with the reactants listed first. The key idea of this approach is that when the concentrations of one or several of these components are held constant, their rows and columns in the conservation matrix are deleted. The maximum number of components that can be held constant is C′ - 1, where C′ is the apparent number of components when the pH is specified. It is important to understand that systems of enzyme-catalyzed reactions generally involve constraints in addition to element balances. These additional constraints arise in the enzymatic mechanisms. The smaller conservation matrix that is obtained by deleting rows and columns for some components can be further reduced in size by deleting redundant columns to obtain the apparent conservation matrix A′′. Reactants that have the same entries at specified concentrations of coenzymes are pseudoisomers, and the standard further transformed Gibbs energies of formation ∆fGi′′o of these pseudoisomer groups can be calculated using eqs 9 and 10. Thus, when the concentrations of several coenzymes are held constant, the reactants and reactions can be reconceptualized to obtain a simpler description of the thermodynamics of the system. The equilibrium composition of this simpler system can be calculated using equilibrium expressions for apparent equilibrium constants K′′ if the standard Gibbs energies of formation ∆fG′o are known for all of the reactants. Equilibrium concentrations have not been calculated for glycolysis and the citric acid cycle because ∆fG′o values for only 22 of the 32 reactants are currently known. This approach provides complete information on the thermodynamics of the system of enzymatic reactions because the equilibrium concentrations of all the reactants other than coenzymes can be calculated. The treatments given here for glycolysis, the citric acid cycle, and glycolysis and the citric acid cycle combined show that larger and more complicated systems of biochemical reactions can be treated in this way. Appendix A The stoichiometric number matrix for the net reaction for a reaction system is given by the dot product of the stoichiometric

ν‚path ) νnet

Appendix B The following Mathematica function will take a vector of stoichiometric numbers (c_List) for a reaction and a vector of names of reactants in quotes (s_List) and print out the reaction: 16,17

mkeqm[c_List, s_List] :) Map[Max[#, 0] &, -c].s f Map[Max[#, 0] &, c].s The following function (in the same notebook as mkeqm) will take the transpose of the stoichiometric number matrix for a reaction sytem (m_List) and a vector of names of reactants in quotes (s_List) and print out all of the reactions:

nameMatrix[m_List, s_List]: ) Map[mkeqm[#, s] &, m] For example, if nu2T is the transposed stoichiometric number matrix and names3 is the corresponding list of names in quotes, the following function will print out a list of all of the reactions:

nameMatrix[nu2T, names3] These functions are very useful for checking a stoichiometric number matrix to be sure that there are no errors. Acknowledgment. I am indebted to Dr. Robert N. Goldberg and Prof. Irwin Oppenheim for many helpful discussions. I am indebted to NIH for support of this research by award number 2-R01-GM48358-04A2. References and Notes (1) Alberty, R. A. Biophys. Chem. 1992, 42, 117. (2) Alberty, R. A. Biophys. Chem. 1992, 43, 239. (3) Beattie, J. A.; Oppenheim, I. Principles of Thermodynamics; Elsevier: Amsterdam, 1979. (4) Alberty, R. A. Biophys. Chem. 1994, 49, 251. (5) Smith, W. R.; Missen, R. W. Chemical Reaction Equilibrium Analysis: Theory and Algorithms; Wiley-Interscience: New York, 1982. (6) Strang, G. Linear Algebra and Its Applications, 3rd ed.; Harcourt Brace Jovanovich: San Diego, 1988. (7) Alberty, R. A. Biophys. J. 1993, 65, 1243. (8) Alberty, R. A. J. Phys. Chem. B 2000, 104, 650. (9) Callen, H. B. Thermodynamics and an Introduction to Thermostatistics, 2nd ed.; Wiley: New York, 1985; p 465. (10) Goldberg, R. N.; Tewari, Y. Biophys. Chem. 1991, 40, 241. (11) Wolfram Research, Inc., 100 World Trade Center, Champaign, IL 61820-7237. (12) Alberty, R. A. Arch. Biochem. Biophys. 1998, 353, 116. (13) Alberty, R. A. Arch. Biochem. Biophys. 1998, 358, 25. (14) Alberty, R. A. J. Phys. Chem. 1995, 99, 1028. (15) Alberty, R. A. Arch. Biochem. Biophys. 1997, 348, 116. (16) Alberty, R. A. Biophys. J. 1996, 71, 507. (17) Alberty, R. A. Biophys. J. 1997, 72, 2349.