- Email: [email protected]
Calculation of chemical equilibrium concentrations of complexing ligands and metals
Taianta, Vol. 20, pp. 1009-1016.
Pergamon
Press, 1913. Printed in Great Britain.
CALCULATION CONCENTRATIONS
OF CHEMICAL EQUILIBRIUM OF COMPLEXING LIGANDS AND METALS
A FLEXIBLE COMPUTER PROGRAMME TAKING UNCERTAINTY IN FORMATION CONSTANT
ACCOUNT VALUES*
OF
G. A. CUMME Institut fiir PhysiologischeChemie, Friedrich-Schiller-Universitat,Jena, G.D.R. (Received 6 September 1912. Accepted 14 February 1973)
Summary-Metal complexesoften act as substratesin enzymicreactions. Becausetheir concentrations cannot be measured directly, they must be calculated from the formation constants and the total concentrations of the different metal and ligand species. An iteration algorithm devised by Sayce has been used to formulate a FORTRAN (ICT 1900) programme, which meets the following requirements. (I) The programme may be made to calculate the total concentrations for the parent species if the free concentration of any metal or ligand species is known, or to calculate the free concentrations of all species for which the total concentrations are known. This feature of the programme is useful for calculation of the total amounts which must be put into a solution in order to achieve prescribed free concentrations. (2) To allow for the uncertainty in the formation constants of the complex species the programme permits the input values of the formation constants to be varied automatically and the resulting variation of the calculated concentrations to be evaluated. Consider a solution which contains given total molar concentrations of metals 1, 2, . . . M, and ligands 1, 2, . . . L, and assume the pH to have been measured. The metals, ligands, Hf and OH- may form complex species, each of which may be described uniquely by the number of metal, ligand, H+ and OH- ions contained in it. At chemical equilibrium and given ionic strength, CK(k) = BK(k) . CHNHck)*4
CM(m)NM@“3 k, *,B CL(I)NL(f*k,
(1)
where CK(k) BK(k) CH NH(k)
= = = =
CM(m) = NM(m, k) = CL(Z) = NL(1, k) =
concentration of the kth complex species cumulative formation constant concentration of H+ number of H+ ions contained in the kth complex species [NH(k) < 0 means that the kth species contains -NH(k) OH- ions] free concentration of mth metal ion species number of ions of the mth metal contained in the kth complex species free concentration of the Ith ligand species number of particles of the Ith ligand contained in the kth complex species
* Part of this paper (without/the Information nationales 1970.
Berliner Symposium
iiber Struktur
on programme structure) was reported at the “VI. Interund Funktion der Erythroryten,” Berlin, G.D.R., August
1010
G. A.
tiMME
Here, the term “complex species ” is used in a generalized manner, so as to include protonated ligands and metal hydroxo-complexes. The concentrations of the complex species, CK, can be calculated easily according to equation (1) from BK, CH, CM and CL. By summing the products NM(m, k) * CK(k) for all k values, and adding CM(m), the total concentration CTM(m) of the mth metal ion is obtained. The total ligand concentration CTL(1) may be calculated similarly. If the free concentrations are to be calculated from given total concentrations, as is usually the case, a system of (M + L) non-linear equations giving CTM(l), . . . , CTM(M), CTL(l), . . . , CTL(L) as functions of the free concentrations has to be solved for the (M + L) unknown free concentrations. Computer programmes capable of solving this problem have been published by various authors, e.g., Perrin and Sayce,’ and Sill& et al.’ The algorithm devised by Sayce, has proved satisfactory in our programme for the largest system so far treated, with 3 metals, 22 ligands, and 93 complex species. Sillen’s more complicated programme’ permits the treatment of solid phases and the use of the mass-balance equation for the proton (Sayce’s very simple algorithm needs the pH to be given, and his programme does not consider solid phases), but these features are not necessary for the biochemical application in question. The algorithm of Sayce has been incorporated into a main programme which had to meet the following requirements. (1) In theoretical studies the question often arises as to how the concentration of one free or complex species depends on the free rather than on the total concentration of another given species. Thus, the programme should permit the free concentrations to be given (and the total concentrations to be calculated) for any subgroup of metals and ligands. (2) Often it is desired to repeat the calculation of free and complex species with only a few parameters altered (e.g., a single total concentration increased). In such cases, it should not be necessary to read in all input parameters anew. (3) The cumulative formation constants available are usually subject to varying degrees of uncertainty. Therefore, the programme should be able to vary the constants within their given limits and to calculate the deviations of free and complex concentrations produced thereby. Because for given total concentrations the free and complex concentrations are non-linear functions of the complex formation constants, their partial derivatives may change their sign. If they do not change sign, the largest and smallest possible free and complex concentrations occur at the BK limits. If, however, the partial derivatives do change in sign, maxima and minima of the free and complex concentrations may occur between the BK limits as well as at them, and it is more complicated to find them. The task of deciding of whether or not this is the case and of determining the absolute extrema has not been solved by our programme. Instead, different estimates of the variation of free and complex concentrations are made, as described later. In the second part of the paper, a description of the programme is given. An application is presented in the third section. The block scheme of the programme, giving detailed information about programme structure, is presented in an appendix. Copies of the ICT 1900 FORTRAN programme are available on request (to a limited extent). DESCRIPTIOIS
OF THE
PROGRAMME
Basic parts (1) Data input. This part has been subdivided according.to (see appendix).
the, different types of data
Calculation of equilibrium concentrations
1011
(2) Single calculation of equilibrium concentrations and print-out of results. (3) Calculation of equilibrium concentrations with medium values of formation constants and print-out of results, followed by calculation of concentration differences produced by variation of formation constants. Let C, be a calculated free metal, free ligand or complex concentration (using the medium value for each formation constant). Let Ci be the concentration value whichis obtained if the ith formation constant is altered in a prescribed manner (i = 1, 2, . . .) KV, where KV = number of formation constants to be varied). Then the differences Ci - Co (i = I-KV) are calculated and printed for any of the above-mentioned concentrations. Additionally, from the relationship AZ = “c (Ci - C,)’ i=
1
A is determined, which gives an estimate of the expected overall uncertainty of C. (4) Instead of the differences Ci - C,, partial derivatives are determined for the calculated concentrations with respect to the logarithms of the formation constants, by solving the corresponding system of linear equations. The quantities D2 and D (see appendix, item 38) are calculated by using the linear terms of the Taylor series, giving the concentrations as functions of the formation constants. (5) Special values are assigned to the formation constants, which may be expected to produce large or small values of calculated concentrations. This part of the programme makes use of the results of parts (3) or (4). To produce a large value of a concentration C, the maximum value is given to the formation constant BK(k) if C is increased by increasing BK(L), and uice versa for producing small C. As long as partial derivatives X/aBK(k) do not change sign, this approach should yield the absolute maxima and minima of the concentrations within the BK-limits given. However, such changes of sign are possible and therefore the calculated concentrations are not necessarily identical with the absolute extrema. On the other hand, because this method requires only constancy of sign of the partial derivatives, it describes the variability of concentrations within a larger range of BK values than does estimation of variability with the aid of the partial derivatives (which requires constancy of the derivatives themselves). It may be used for computing minimum ranges of the concentrations. Transferring control to the dlrerent parts of the programme via punched card input
In each of the above-mentioned parts of the programme, the last statement causes a new punched card (later referred to as a steering-card) to be read in, which bears an integer number Il. By aid of 11 the computer learns which part of the programme it has to go to. It is clear that such a flexible programme organization permits the programme to be adapted to quite different sequences of readings and calculations by simply arranging the punched cards in the appropriate order. Repeating a calculation with only a few parameters altered is greatly facilitated in this way, because only those few parameters are read in by the computer, leaving all others unaltered. EXAMPLE
OF A BIOLOGICAL
APPLIdATION
In the human erythrocyte the most important constituents are ATP, 2,3-DPG, magnesium, potassium and sodium ions, and the 16 possible complex species formed from them. For biochemical purposes it is of interest to know not only their total concentrations, but also
1012
G. A. CUMME Table 1. Logarithms of cumulative formation constants, referring to concentrations. Complex species 1. 2. 3. 4. 5. 6. 7. 8. 9. 10. 11. 12. 13. 14. 15. 16.
HDFG H,DPG HsDPG KDPG HKDPG NaDPG HNaDPG MgDPG HMgDPG Mg,DPG MgATP HM@TP KATP NaATp HATP H,ATP
log WC) 8.21 14.84 18.36 1893 P-16 1.93 P-19 3.87 IO.92 6.12 4.96 9.80 l-15 1.18 6.95 10.88
the free concentrations, and the concentrations of the complex species. Values for the constants have been taken from O’Sullivan and Perrin” and from the references quoted by them, together with new results of Achilles et al.’ These constants are presented in Table 1. The results of the equilibrium calculations are shown in Table 2. Total concentrations of metals and ligands as found in human erythrocytes have been chosen (the values are given in column 1; they coincide exactly with those selected by Achilles et al.‘). P, = 7.2 From the equilibrium concentrations presented in column 2, the importance of the potassium complexes with 2,3-DPG and the magnesium complex with ATP becomes quite clear (cJ Achilles et aL3). For each single concentration value, the variation range due to the uncertainty of the association constants has been estimated according to basic part (5), allowing each association constant to be multiplied by 1.25 or O-8 (corresponding roughly to &O-l pK units). Whereas the importance of magnesium and potassium complex-formation with ATP and 2,3-DPG respectively is not removed by varying the constants, the single values are found to be altered remarkably (c$ columns 3a and 3b). Thus, in comparing biochemical model properties with observational results, the possibility that special features of the model are sensitive to changes within the uncertainty ranges of the complex formation constants should be taken into account. Column 4 of Table 2 gives some insight into the dependence of the calculated concentrations on formation constants. For each concentration C, the two complex formation constants causing (by their variation) the greatest alteration of C are given. The concentration of free 2,3-DPG depends most strongly on the formation constants of the complexes HKDPG and KDPG (numbers 5 and 4); these are the two DPG-complexes with the highest concentrations (cJ column 2). The concentration of free ATP is most strongly influenced by the formation constants of the complexes MgATP and MgDPG. Because ATP is bound nearly exclusively by formation of MgATP, the MgDPG complex influences the free ATP via Mg more strongly than the other ATP complexes do. For the free magnesium ion concentration, we find the formation constants of MgDPG and HMgDPG. Increasing the MgATP formation constant produces a smaller decrease in free magnesium concentration,
Calculation
of equilibrium
Table 2. Results of different equilibrium 1 Species (tot. concn.) mM_ DPG (7.2) ATP (20) Mg (35) K (130.0) Na (20.0) 1. HDPG 2. HzDPG 3. HjDPG 4. KDPG 5. HKDPG 6. NaDPG 7. HNaDPG 8. MgDPG 9. HMgDPG 10. MglDPG 11. MgATP 12. HMgATP 13. KATP 14. NaATP 15. HATP 16. HIATP
2 Concn. (med. formation constants) mM 0.154 0.0390 0.520 1265 19.4 1.57 0.424 8*86E-5 1.66 1.78 0.255 0.292 @593 0.420 0.0548 1.85 8*07E-3 0.0697 11.5E-3 2.19E-2 11.8E-6
calculations 4 Most important
3b
3a Concn. min. mM 0.123 0.0253 0.388 125.8 19.2 1.10 0.278 566E-5 1.17 1.26 0.166 0.191 0.423 0288 QO309 1.74 5*09E-3 0.0369 5*97E-3 0.0114 6.1lE-6
1013
concentrations
differences
max. mM
mM 0192 0.0612 0.674 127.2 19.6 2.16 0.638 13.8E-5 2.26 240 0.387 0442 @806 0.597 0.0887 1.91 12.7E-3 0131 22.3E- 3 0.0424 23*lE-6
5; 11; 8; 1; 7; ;: 3; 4; 5; 6; 7; 8; 9; 10; 11; 12; 13; 14; 15; 16;
-9*67E-3 -685E-3 -0.0353 0.195 -0.0642 0.280 0.0974 2.21E-5 0.289 0.301 0.0597 0.0679 00919 0.0761 0.0118 0.0265 2*01E-3 0.0159 2*83E-6 5*34E-3 2*95E-6
4; 8; 9; 5; 6; 5;
-9*07E-3 2.638 - 3 -0.0255 -0.192 - 0.0560 -0.0990
:; - -0.0266 -5*56E-6 :: * -0.107 5; 5; 9; 8; 8; 13; 11; 11; 11; 11; 11;
-0.0156 -0.0179 -0.0326 -0.0319 -7+iOE-3 -0.0146 -1.52E-3 -0.0122 -2*02E-6 -3.858-3 -207E-6
Note: uE - b means a. 10Tb Column 1: free and complex species considered. For metals and ligands, total concentrations (mM) are given in parentheses. Column 2: concentrations as calculated from total concentrations and cumulative complex formation constants (cjI Table 1). Columns 3a and 3b: each complex formation constant was allowed to be multiplied by either 1.25 or 0.8, the factors being chosen so as to minimize or maximize respectively the concentration of each species mentioned in column 1. Thus for each species another combination of values for the cumulative formation constants had to be used (cf. part (3) of basic programme parts). Column 4: in each line, “ i;c j;d” means : the concentration as given in column 2 is incremented by c or d, if the formation constant of the ith (or jth) complex species is multiplied by 1.25~ multiplying another formation constant by the same factor leads to a smaller increment than d.
because ATP is nearly saturated by magnesium and the.concentration of bound substrate (MgATP) depends only weakly on the formation constant within the saturation range. In the case of potassium, its large total concentration causes the free potassium concentration to be nearly independent of the choice of formation constants (c$ columns 3a and 3b). Increasing the formation constant of HDPG decreases the concentrations of KDPG and HKDPG and thus produces more free potassium. Increasing the constant for KDPG or HKDPG increases one complex at the expense of the other, but in total more potassium than before is bound and thus the free potassium concentration decreases. As to the complex species, each complex concentration depends most strongly upon the formation constant for that complex. Secondarily, many ATP-complexes depend on the formation constant of MgATP (which exerts the strongest influence on free ATP concentration) and many DPG-complex concentrations depend on the formation constant of HKDPG (which exerts the strongest influence on free DPG concentration). Acknowledgement-The
author is indebted to Dr. W. Achilles for his advice while preparing the text.
1014
1. 2. 3. 4. 5.
G. A.
CUMME
REFERENCES D. D. Perrin and I. G. Sayce. Talunru, 1967, 14, 833. W. I. O’Sullivan and D. D. P&n, Biochemistry, 1964, 3, 18. W. Achilles, G. A. Cumme, and H. Hoppe, Acm Biol. Med. Germ. (to be published). G. An&t-egg. Heltr. Chim. Actu, 1962, 45, 901. N. Ingri, W. Kakolowicz, L. G. Sill6n and B. Wamqvist, Talunru, 1967, 14, 1261. APPENDIX Block sckeme of the programme
1
Read in deck of punched cards, the last card of which bears the word ” ENDKARTE” in its columns nos. 1-8, and write a picture of each card on magnetic tape. Thereafter release-card-reader from the programme (this is for machines permitting multi-programming), and rewind magnetic tape. 2 Read in picture of a punched card from magnetic tape. If this picture by its structure turns out to be a picture of a steering-card*, go to 2b, if not, go to 2a. Brint appropriate message and finish actual job. 2 According to the integer number 11’ carried by the steering-card, go to the 11th item of the following list 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 20, 21, 22, 23, 24, 25, 26, 27, 39, 28, 29, 30, 31, 38, 18, 19 (e.g., a steering-card with 11 = 29 causes the computer to go to item 38 (see below), i.e. to calculate partial derivatives; whereas 11 = 24 would correspond to item 39, that is to finish the actual job. Thus, in each punched card deck the last steering-card has 11 = 24; it is followed immediately by a card with the word ‘ENDKARTE’). 3 Read from magnetic tape (hereafter abbreviated as m.t.) two card pictures containing an arbitrary text; write this text, time, and date; go to 2. 4 Read from m.t., and write EPS and ITLIM (in calculating free concentrations by iteration, the improvement of the free concentrations is ceased, if the calculated total concentrations coincide with their known values within relative errors of EPS, or if the number of improvements done exceeds ITLIM); go to 2. 5 Read from m.t. and write the integer numbers KEY 1, KEY 2, KEY 3 (the values of which are needed at items 31, 32, 35, 36 and 37); go to 2. 6 Read from m.t. and write the numbers of metals, ligands, and complex species (MM, ML, MK), the numbers MI and LI (for metals nos. l-MI, total concentrations are known, whereas for metals from MI + 1 to MM, the free concentrations are known, and vice uersa for ligands), and KV (for complex species nos. I-KV, the formation constants are varied automatically by the programme, see below); go to 2. 7 (This item is of minor importance and is therefore not explained here) 8-10 Read from m.t. and write the names of complex species, metals, and ligands respectively; go to 2. Read from m.t. and write NH(k), NM(m, k), NL(I, k), k = 1, . . . , MK; m = 1. . . . , MM; I = 1. . . . , 11 ML; go to 2. (This item is of minor importance and is therefore not explained here) :: Read from m.t. and write -log[H+] and -log([H+]. [OH-]); go to 2. 14-17 Read from m.t. and write given values for total metal, free metal, total ligand, and free ligand concentrations respectively; go to 2. 18 and (These items are of minor importance and are therefore not explained here) 19 Read from m.t. and write decadic logarithms of the cumulative formation constants; go to 2. 20 Calculate cumulative formation constants BK(k), k = 1, . . . , MK, according to the NH sign conven21 tion (c/. introductory section of oaoer) from the data read in at items 11. 13, and 20: go to 2. Read &om m.t. VA; the amouni of relative variation of cumulative formation constants; go to 2. 22 Write VA: calculate maximum values BK(k)*(l + VA) and minimum values BKQ)/(l 23 ., + VA1 for cumulative formation constants nos. l-KV; b’to 2. 24 Raad from m.t. and write individual maximum and minimum logarithms of cumulative formation constants nos. l-KV; go to 2. Calculate individual maximum and minimum values of cumulative formation constants iios. l-KV 25 according to the NH sign convention from the data read in at items 11, 13, and 24; go to 2. 26 and (These items are of minor importance and are not explained here) 27 * CJ section on transfer of control via punched card input.
Calculation of equilibrium concentrations
lOI
Look for so-called sub-complexes. A combination of metabolites A,. ..Z,is said to be a sub-complex, if these metabolites do appear in an arbitrary complex species exclusively as multiples of this combination (e.g., if the complex species present are A 2, ABC,, AD, DB&, DE, then BC2 is a subcomplex). If a sub-complex A,. ..Z, is present, and the equilibrium free concentrations of A-Z are much smaller than the total ones, then the calculated equilibrium concentrations are to be looked at with caution and may need improvement.* This case, however, occurred very seldom in our studies. Thus no improvement has been provided as yet. After printing an appropriate message; go to 2. 29 and (These items are’of minor importance and are not explained here) 30 If KEY 1 = 1, assign starting values for subsequent iterative improvement to free metal and ligand 31 concentrations; go to 32. If KEY I = 2, go immediately to 32. If KEY 3 = 1, 2, 3, 4 or 5, go to 33, 34, 35, 36 or 37 respectively. 32 Calculate free and complex equilibrium concentrations and print results; go to 2. Replace each cumulative formation constant by its medium value. Calculate equilibrium free and :: complex concentrations and prim results. Go to 2. (This branch of the programme is used, if after variation of formation constants a further calculation with medium values of formation constants but, e.g.: altered total concentrations is wanted.) Perform computations according to basic part (3) in or&r to calculate di&rences of free and 35 complex concentrations caused by variation of each formation constant the no. of which lies within l... KV. The alteration of each formation constant is done by assigning to it its maximum value (which has been calculated at item 23 or 25). After each equilibrium calculation, print out calculated concentrations, if KEY 2 = 1, and omit printing, if KEY 2 = 2. After varying the last formation constant, print list of concentration differences AC as caused by variation of formation constants BK(k) (one difference for each calculated concentration, and each varied formation constant). Calculate and print out A [see basic part (3)] for each calculated free metal, free ligand, and complex concentration. 36 and Perform computations according to basic part (5) in order to estimate ranges of calculated free 37 concentrations (36) or complex concentrations (37) due to variation of cumulative formation constants nos. l... KV. (These calculations require the calculations at item 35 or item 38 to have been done.) Print out calculated concentrations for each combination of formation constants if KEY 2 = 1, and omit printing if KEY 2 = 2. At any case print list of ranges for all calculated free or complex conce-ntrations. Then go to 2. 38 Compute partial derivatives X/a log BK(k) of the calculated free and complex concentrations C with respect to the logarithms of the formation constants log BK(&) by solving the corresponding system of linear equations (this calculation requires the calculations at item 33 or 34 to have been done). Compute for each calculated free and complex concentration C 28
Dz = k$, [aC/aBK(k) . ABK(/#
39
with ABK = maximum minus medium value of BK (cf. items 21-25). Calculate D as well, D giving an estimation for the expected overall uncertainty of C due to the uncertainty of the formation constants. Then go to 2. Write “Job done, ” time and date. Finish actual job. Zusammenfassung-Metallkomplexe tfeten oft als Substrate in enzymatischen Reaktionen auf. Da ihre Konzentrationen nicht direkt gemessen werden k&men, miiasen sie aus den Bildungskonstanten und den Gesamtkonzentrationen der verschiedenen Metall- und Ligandenspezies berechnet werden. Ein von Sayce erstellter Iterations-Algorithmus wurde dazu verwendet, ein Fortran-(ICT 1900) Programm zu formulieren, das folgende Anforderungen erfilllt : 1) Mit dem Programm kann man die Gesamtkonzentrationen der Ausgangsspezies berechnen, wenn die Konzentrationen allerfreien Metall- und Ligandenspezies bekannt sind, oder die Konzcntrationen aller freien Spezies zu berechnen, von denen die Gesamtkonzentration bekannt ist. Diese Eigenschaft des Programms ist dann von Nutzen, wenn man die Gesamtmengen berechnen will, die man in eine Lasung einsetzen mu& urn vorgeschriebene Konzentrationen an freien Spezies zu erhalten. 2) Urn die Unsicherheit da Bildungskonstanten der Komplexe zu bertlcksichtigen, erlaubt das Programm, die eingespeisten Werte der Bildungskonstanten automat&h zu variieren und die sich ergebende Xnderung der berechneten Konzentrationen zu ermitteln.
l This situation may also a& if complex species containing metabolites A-Z in a proportion different from that of the sub-complex A,. ..Z, are present, but with negligible concentrations (c$ G. Anderegg).
1016
G. A.
CUMME
R&m&-Les complexes mttalliques a&sent souvent comme substrats dans des riactions enzymatiques. Parce que leurs concentrations ne peuvent Btre mesurkes directement, elles doivent 6tre calcukes h partir des constantes de formation et des soncentrations totales de diffhentes es&es de m&w et coordinats. Un algorithme itkratif babli par Sayce a ktt utilisk pour formuler u11programme FORTRAN (In BOO),qui r&pond aux exigence suivantes. (1) Le programme peut &re fait pour calculer les concentrations totales pour les espkces apparentks si la concentration librede n’importe quelle espkce de m&al ou de coordinat est connue ou pour calculer les concentrations libres de toutes les espkes pour lesquelles les concentrations totales sont connues. Cette caractkistique du programme est utile pour le calcul des quantitis totales qui doivent 2tre introduites dans une solution afin d’obtenir les concentrations libres prescrites. (2) Pour tenir compte de l’incertitude dans les constantes de formation de l’espkce complexe le programme permet de faire varier automatiquement les valeurs d’entrtk des constantes de formation et d’bvaluer la variation rksultante des concentrations calcuks.
Pergamon
Press, 1913. Printed in Great Britain.
CALCULATION CONCENTRATIONS
OF CHEMICAL EQUILIBRIUM OF COMPLEXING LIGANDS AND METALS
A FLEXIBLE COMPUTER PROGRAMME TAKING UNCERTAINTY IN FORMATION CONSTANT
ACCOUNT VALUES*
OF
G. A. CUMME Institut fiir PhysiologischeChemie, Friedrich-Schiller-Universitat,Jena, G.D.R. (Received 6 September 1912. Accepted 14 February 1973)
Summary-Metal complexesoften act as substratesin enzymicreactions. Becausetheir concentrations cannot be measured directly, they must be calculated from the formation constants and the total concentrations of the different metal and ligand species. An iteration algorithm devised by Sayce has been used to formulate a FORTRAN (ICT 1900) programme, which meets the following requirements. (I) The programme may be made to calculate the total concentrations for the parent species if the free concentration of any metal or ligand species is known, or to calculate the free concentrations of all species for which the total concentrations are known. This feature of the programme is useful for calculation of the total amounts which must be put into a solution in order to achieve prescribed free concentrations. (2) To allow for the uncertainty in the formation constants of the complex species the programme permits the input values of the formation constants to be varied automatically and the resulting variation of the calculated concentrations to be evaluated. Consider a solution which contains given total molar concentrations of metals 1, 2, . . . M, and ligands 1, 2, . . . L, and assume the pH to have been measured. The metals, ligands, Hf and OH- may form complex species, each of which may be described uniquely by the number of metal, ligand, H+ and OH- ions contained in it. At chemical equilibrium and given ionic strength, CK(k) = BK(k) . CHNHck)*4
CM(m)NM@“3 k, *,B CL(I)NL(f*k,
(1)
where CK(k) BK(k) CH NH(k)
= = = =
CM(m) = NM(m, k) = CL(Z) = NL(1, k) =
concentration of the kth complex species cumulative formation constant concentration of H+ number of H+ ions contained in the kth complex species [NH(k) < 0 means that the kth species contains -NH(k) OH- ions] free concentration of mth metal ion species number of ions of the mth metal contained in the kth complex species free concentration of the Ith ligand species number of particles of the Ith ligand contained in the kth complex species
* Part of this paper (without/the Information nationales 1970.
Berliner Symposium
iiber Struktur
on programme structure) was reported at the “VI. Interund Funktion der Erythroryten,” Berlin, G.D.R., August
1010
G. A.
tiMME
Here, the term “complex species ” is used in a generalized manner, so as to include protonated ligands and metal hydroxo-complexes. The concentrations of the complex species, CK, can be calculated easily according to equation (1) from BK, CH, CM and CL. By summing the products NM(m, k) * CK(k) for all k values, and adding CM(m), the total concentration CTM(m) of the mth metal ion is obtained. The total ligand concentration CTL(1) may be calculated similarly. If the free concentrations are to be calculated from given total concentrations, as is usually the case, a system of (M + L) non-linear equations giving CTM(l), . . . , CTM(M), CTL(l), . . . , CTL(L) as functions of the free concentrations has to be solved for the (M + L) unknown free concentrations. Computer programmes capable of solving this problem have been published by various authors, e.g., Perrin and Sayce,’ and Sill& et al.’ The algorithm devised by Sayce, has proved satisfactory in our programme for the largest system so far treated, with 3 metals, 22 ligands, and 93 complex species. Sillen’s more complicated programme’ permits the treatment of solid phases and the use of the mass-balance equation for the proton (Sayce’s very simple algorithm needs the pH to be given, and his programme does not consider solid phases), but these features are not necessary for the biochemical application in question. The algorithm of Sayce has been incorporated into a main programme which had to meet the following requirements. (1) In theoretical studies the question often arises as to how the concentration of one free or complex species depends on the free rather than on the total concentration of another given species. Thus, the programme should permit the free concentrations to be given (and the total concentrations to be calculated) for any subgroup of metals and ligands. (2) Often it is desired to repeat the calculation of free and complex species with only a few parameters altered (e.g., a single total concentration increased). In such cases, it should not be necessary to read in all input parameters anew. (3) The cumulative formation constants available are usually subject to varying degrees of uncertainty. Therefore, the programme should be able to vary the constants within their given limits and to calculate the deviations of free and complex concentrations produced thereby. Because for given total concentrations the free and complex concentrations are non-linear functions of the complex formation constants, their partial derivatives may change their sign. If they do not change sign, the largest and smallest possible free and complex concentrations occur at the BK limits. If, however, the partial derivatives do change in sign, maxima and minima of the free and complex concentrations may occur between the BK limits as well as at them, and it is more complicated to find them. The task of deciding of whether or not this is the case and of determining the absolute extrema has not been solved by our programme. Instead, different estimates of the variation of free and complex concentrations are made, as described later. In the second part of the paper, a description of the programme is given. An application is presented in the third section. The block scheme of the programme, giving detailed information about programme structure, is presented in an appendix. Copies of the ICT 1900 FORTRAN programme are available on request (to a limited extent). DESCRIPTIOIS
OF THE
PROGRAMME
Basic parts (1) Data input. This part has been subdivided according.to (see appendix).
the, different types of data
Calculation of equilibrium concentrations
1011
(2) Single calculation of equilibrium concentrations and print-out of results. (3) Calculation of equilibrium concentrations with medium values of formation constants and print-out of results, followed by calculation of concentration differences produced by variation of formation constants. Let C, be a calculated free metal, free ligand or complex concentration (using the medium value for each formation constant). Let Ci be the concentration value whichis obtained if the ith formation constant is altered in a prescribed manner (i = 1, 2, . . .) KV, where KV = number of formation constants to be varied). Then the differences Ci - Co (i = I-KV) are calculated and printed for any of the above-mentioned concentrations. Additionally, from the relationship AZ = “c (Ci - C,)’ i=
1
A is determined, which gives an estimate of the expected overall uncertainty of C. (4) Instead of the differences Ci - C,, partial derivatives are determined for the calculated concentrations with respect to the logarithms of the formation constants, by solving the corresponding system of linear equations. The quantities D2 and D (see appendix, item 38) are calculated by using the linear terms of the Taylor series, giving the concentrations as functions of the formation constants. (5) Special values are assigned to the formation constants, which may be expected to produce large or small values of calculated concentrations. This part of the programme makes use of the results of parts (3) or (4). To produce a large value of a concentration C, the maximum value is given to the formation constant BK(k) if C is increased by increasing BK(L), and uice versa for producing small C. As long as partial derivatives X/aBK(k) do not change sign, this approach should yield the absolute maxima and minima of the concentrations within the BK-limits given. However, such changes of sign are possible and therefore the calculated concentrations are not necessarily identical with the absolute extrema. On the other hand, because this method requires only constancy of sign of the partial derivatives, it describes the variability of concentrations within a larger range of BK values than does estimation of variability with the aid of the partial derivatives (which requires constancy of the derivatives themselves). It may be used for computing minimum ranges of the concentrations. Transferring control to the dlrerent parts of the programme via punched card input
In each of the above-mentioned parts of the programme, the last statement causes a new punched card (later referred to as a steering-card) to be read in, which bears an integer number Il. By aid of 11 the computer learns which part of the programme it has to go to. It is clear that such a flexible programme organization permits the programme to be adapted to quite different sequences of readings and calculations by simply arranging the punched cards in the appropriate order. Repeating a calculation with only a few parameters altered is greatly facilitated in this way, because only those few parameters are read in by the computer, leaving all others unaltered. EXAMPLE
OF A BIOLOGICAL
APPLIdATION
In the human erythrocyte the most important constituents are ATP, 2,3-DPG, magnesium, potassium and sodium ions, and the 16 possible complex species formed from them. For biochemical purposes it is of interest to know not only their total concentrations, but also
1012
G. A. CUMME Table 1. Logarithms of cumulative formation constants, referring to concentrations. Complex species 1. 2. 3. 4. 5. 6. 7. 8. 9. 10. 11. 12. 13. 14. 15. 16.
HDFG H,DPG HsDPG KDPG HKDPG NaDPG HNaDPG MgDPG HMgDPG Mg,DPG MgATP HM@TP KATP NaATp HATP H,ATP
log WC) 8.21 14.84 18.36 1893 P-16 1.93 P-19 3.87 IO.92 6.12 4.96 9.80 l-15 1.18 6.95 10.88
the free concentrations, and the concentrations of the complex species. Values for the constants have been taken from O’Sullivan and Perrin” and from the references quoted by them, together with new results of Achilles et al.’ These constants are presented in Table 1. The results of the equilibrium calculations are shown in Table 2. Total concentrations of metals and ligands as found in human erythrocytes have been chosen (the values are given in column 1; they coincide exactly with those selected by Achilles et al.‘). P, = 7.2 From the equilibrium concentrations presented in column 2, the importance of the potassium complexes with 2,3-DPG and the magnesium complex with ATP becomes quite clear (cJ Achilles et aL3). For each single concentration value, the variation range due to the uncertainty of the association constants has been estimated according to basic part (5), allowing each association constant to be multiplied by 1.25 or O-8 (corresponding roughly to &O-l pK units). Whereas the importance of magnesium and potassium complex-formation with ATP and 2,3-DPG respectively is not removed by varying the constants, the single values are found to be altered remarkably (c$ columns 3a and 3b). Thus, in comparing biochemical model properties with observational results, the possibility that special features of the model are sensitive to changes within the uncertainty ranges of the complex formation constants should be taken into account. Column 4 of Table 2 gives some insight into the dependence of the calculated concentrations on formation constants. For each concentration C, the two complex formation constants causing (by their variation) the greatest alteration of C are given. The concentration of free 2,3-DPG depends most strongly on the formation constants of the complexes HKDPG and KDPG (numbers 5 and 4); these are the two DPG-complexes with the highest concentrations (cJ column 2). The concentration of free ATP is most strongly influenced by the formation constants of the complexes MgATP and MgDPG. Because ATP is bound nearly exclusively by formation of MgATP, the MgDPG complex influences the free ATP via Mg more strongly than the other ATP complexes do. For the free magnesium ion concentration, we find the formation constants of MgDPG and HMgDPG. Increasing the MgATP formation constant produces a smaller decrease in free magnesium concentration,
Calculation
of equilibrium
Table 2. Results of different equilibrium 1 Species (tot. concn.) mM_ DPG (7.2) ATP (20) Mg (35) K (130.0) Na (20.0) 1. HDPG 2. HzDPG 3. HjDPG 4. KDPG 5. HKDPG 6. NaDPG 7. HNaDPG 8. MgDPG 9. HMgDPG 10. MglDPG 11. MgATP 12. HMgATP 13. KATP 14. NaATP 15. HATP 16. HIATP
2 Concn. (med. formation constants) mM 0.154 0.0390 0.520 1265 19.4 1.57 0.424 8*86E-5 1.66 1.78 0.255 0.292 @593 0.420 0.0548 1.85 8*07E-3 0.0697 11.5E-3 2.19E-2 11.8E-6
calculations 4 Most important
3b
3a Concn. min. mM 0.123 0.0253 0.388 125.8 19.2 1.10 0.278 566E-5 1.17 1.26 0.166 0.191 0.423 0288 QO309 1.74 5*09E-3 0.0369 5*97E-3 0.0114 6.1lE-6
1013
concentrations
differences
max. mM
mM 0192 0.0612 0.674 127.2 19.6 2.16 0.638 13.8E-5 2.26 240 0.387 0442 @806 0.597 0.0887 1.91 12.7E-3 0131 22.3E- 3 0.0424 23*lE-6
5; 11; 8; 1; 7; ;: 3; 4; 5; 6; 7; 8; 9; 10; 11; 12; 13; 14; 15; 16;
-9*67E-3 -685E-3 -0.0353 0.195 -0.0642 0.280 0.0974 2.21E-5 0.289 0.301 0.0597 0.0679 00919 0.0761 0.0118 0.0265 2*01E-3 0.0159 2*83E-6 5*34E-3 2*95E-6
4; 8; 9; 5; 6; 5;
-9*07E-3 2.638 - 3 -0.0255 -0.192 - 0.0560 -0.0990
:; - -0.0266 -5*56E-6 :: * -0.107 5; 5; 9; 8; 8; 13; 11; 11; 11; 11; 11;
-0.0156 -0.0179 -0.0326 -0.0319 -7+iOE-3 -0.0146 -1.52E-3 -0.0122 -2*02E-6 -3.858-3 -207E-6
Note: uE - b means a. 10Tb Column 1: free and complex species considered. For metals and ligands, total concentrations (mM) are given in parentheses. Column 2: concentrations as calculated from total concentrations and cumulative complex formation constants (cjI Table 1). Columns 3a and 3b: each complex formation constant was allowed to be multiplied by either 1.25 or 0.8, the factors being chosen so as to minimize or maximize respectively the concentration of each species mentioned in column 1. Thus for each species another combination of values for the cumulative formation constants had to be used (cf. part (3) of basic programme parts). Column 4: in each line, “ i;c j;d” means : the concentration as given in column 2 is incremented by c or d, if the formation constant of the ith (or jth) complex species is multiplied by 1.25~ multiplying another formation constant by the same factor leads to a smaller increment than d.
because ATP is nearly saturated by magnesium and the.concentration of bound substrate (MgATP) depends only weakly on the formation constant within the saturation range. In the case of potassium, its large total concentration causes the free potassium concentration to be nearly independent of the choice of formation constants (c$ columns 3a and 3b). Increasing the formation constant of HDPG decreases the concentrations of KDPG and HKDPG and thus produces more free potassium. Increasing the constant for KDPG or HKDPG increases one complex at the expense of the other, but in total more potassium than before is bound and thus the free potassium concentration decreases. As to the complex species, each complex concentration depends most strongly upon the formation constant for that complex. Secondarily, many ATP-complexes depend on the formation constant of MgATP (which exerts the strongest influence on free ATP concentration) and many DPG-complex concentrations depend on the formation constant of HKDPG (which exerts the strongest influence on free DPG concentration). Acknowledgement-The
author is indebted to Dr. W. Achilles for his advice while preparing the text.
1014
1. 2. 3. 4. 5.
G. A.
CUMME
REFERENCES D. D. Perrin and I. G. Sayce. Talunru, 1967, 14, 833. W. I. O’Sullivan and D. D. P&n, Biochemistry, 1964, 3, 18. W. Achilles, G. A. Cumme, and H. Hoppe, Acm Biol. Med. Germ. (to be published). G. An&t-egg. Heltr. Chim. Actu, 1962, 45, 901. N. Ingri, W. Kakolowicz, L. G. Sill6n and B. Wamqvist, Talunru, 1967, 14, 1261. APPENDIX Block sckeme of the programme
1
Read in deck of punched cards, the last card of which bears the word ” ENDKARTE” in its columns nos. 1-8, and write a picture of each card on magnetic tape. Thereafter release-card-reader from the programme (this is for machines permitting multi-programming), and rewind magnetic tape. 2 Read in picture of a punched card from magnetic tape. If this picture by its structure turns out to be a picture of a steering-card*, go to 2b, if not, go to 2a. Brint appropriate message and finish actual job. 2 According to the integer number 11’ carried by the steering-card, go to the 11th item of the following list 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 20, 21, 22, 23, 24, 25, 26, 27, 39, 28, 29, 30, 31, 38, 18, 19 (e.g., a steering-card with 11 = 29 causes the computer to go to item 38 (see below), i.e. to calculate partial derivatives; whereas 11 = 24 would correspond to item 39, that is to finish the actual job. Thus, in each punched card deck the last steering-card has 11 = 24; it is followed immediately by a card with the word ‘ENDKARTE’). 3 Read from magnetic tape (hereafter abbreviated as m.t.) two card pictures containing an arbitrary text; write this text, time, and date; go to 2. 4 Read from m.t., and write EPS and ITLIM (in calculating free concentrations by iteration, the improvement of the free concentrations is ceased, if the calculated total concentrations coincide with their known values within relative errors of EPS, or if the number of improvements done exceeds ITLIM); go to 2. 5 Read from m.t. and write the integer numbers KEY 1, KEY 2, KEY 3 (the values of which are needed at items 31, 32, 35, 36 and 37); go to 2. 6 Read from m.t. and write the numbers of metals, ligands, and complex species (MM, ML, MK), the numbers MI and LI (for metals nos. l-MI, total concentrations are known, whereas for metals from MI + 1 to MM, the free concentrations are known, and vice uersa for ligands), and KV (for complex species nos. I-KV, the formation constants are varied automatically by the programme, see below); go to 2. 7 (This item is of minor importance and is therefore not explained here) 8-10 Read from m.t. and write the names of complex species, metals, and ligands respectively; go to 2. Read from m.t. and write NH(k), NM(m, k), NL(I, k), k = 1, . . . , MK; m = 1. . . . , MM; I = 1. . . . , 11 ML; go to 2. (This item is of minor importance and is therefore not explained here) :: Read from m.t. and write -log[H+] and -log([H+]. [OH-]); go to 2. 14-17 Read from m.t. and write given values for total metal, free metal, total ligand, and free ligand concentrations respectively; go to 2. 18 and (These items are of minor importance and are therefore not explained here) 19 Read from m.t. and write decadic logarithms of the cumulative formation constants; go to 2. 20 Calculate cumulative formation constants BK(k), k = 1, . . . , MK, according to the NH sign conven21 tion (c/. introductory section of oaoer) from the data read in at items 11. 13, and 20: go to 2. Read &om m.t. VA; the amouni of relative variation of cumulative formation constants; go to 2. 22 Write VA: calculate maximum values BK(k)*(l + VA) and minimum values BKQ)/(l 23 ., + VA1 for cumulative formation constants nos. l-KV; b’to 2. 24 Raad from m.t. and write individual maximum and minimum logarithms of cumulative formation constants nos. l-KV; go to 2. Calculate individual maximum and minimum values of cumulative formation constants iios. l-KV 25 according to the NH sign convention from the data read in at items 11, 13, and 24; go to 2. 26 and (These items are of minor importance and are not explained here) 27 * CJ section on transfer of control via punched card input.
Calculation of equilibrium concentrations
lOI
Look for so-called sub-complexes. A combination of metabolites A,. ..Z,is said to be a sub-complex, if these metabolites do appear in an arbitrary complex species exclusively as multiples of this combination (e.g., if the complex species present are A 2, ABC,, AD, DB&, DE, then BC2 is a subcomplex). If a sub-complex A,. ..Z, is present, and the equilibrium free concentrations of A-Z are much smaller than the total ones, then the calculated equilibrium concentrations are to be looked at with caution and may need improvement.* This case, however, occurred very seldom in our studies. Thus no improvement has been provided as yet. After printing an appropriate message; go to 2. 29 and (These items are’of minor importance and are not explained here) 30 If KEY 1 = 1, assign starting values for subsequent iterative improvement to free metal and ligand 31 concentrations; go to 32. If KEY I = 2, go immediately to 32. If KEY 3 = 1, 2, 3, 4 or 5, go to 33, 34, 35, 36 or 37 respectively. 32 Calculate free and complex equilibrium concentrations and print results; go to 2. Replace each cumulative formation constant by its medium value. Calculate equilibrium free and :: complex concentrations and prim results. Go to 2. (This branch of the programme is used, if after variation of formation constants a further calculation with medium values of formation constants but, e.g.: altered total concentrations is wanted.) Perform computations according to basic part (3) in or&r to calculate di&rences of free and 35 complex concentrations caused by variation of each formation constant the no. of which lies within l... KV. The alteration of each formation constant is done by assigning to it its maximum value (which has been calculated at item 23 or 25). After each equilibrium calculation, print out calculated concentrations, if KEY 2 = 1, and omit printing, if KEY 2 = 2. After varying the last formation constant, print list of concentration differences AC as caused by variation of formation constants BK(k) (one difference for each calculated concentration, and each varied formation constant). Calculate and print out A [see basic part (3)] for each calculated free metal, free ligand, and complex concentration. 36 and Perform computations according to basic part (5) in order to estimate ranges of calculated free 37 concentrations (36) or complex concentrations (37) due to variation of cumulative formation constants nos. l... KV. (These calculations require the calculations at item 35 or item 38 to have been done.) Print out calculated concentrations for each combination of formation constants if KEY 2 = 1, and omit printing if KEY 2 = 2. At any case print list of ranges for all calculated free or complex conce-ntrations. Then go to 2. 38 Compute partial derivatives X/a log BK(k) of the calculated free and complex concentrations C with respect to the logarithms of the formation constants log BK(&) by solving the corresponding system of linear equations (this calculation requires the calculations at item 33 or 34 to have been done). Compute for each calculated free and complex concentration C 28
Dz = k$, [aC/aBK(k) . ABK(/#
39
with ABK = maximum minus medium value of BK (cf. items 21-25). Calculate D as well, D giving an estimation for the expected overall uncertainty of C due to the uncertainty of the formation constants. Then go to 2. Write “Job done, ” time and date. Finish actual job. Zusammenfassung-Metallkomplexe tfeten oft als Substrate in enzymatischen Reaktionen auf. Da ihre Konzentrationen nicht direkt gemessen werden k&men, miiasen sie aus den Bildungskonstanten und den Gesamtkonzentrationen der verschiedenen Metall- und Ligandenspezies berechnet werden. Ein von Sayce erstellter Iterations-Algorithmus wurde dazu verwendet, ein Fortran-(ICT 1900) Programm zu formulieren, das folgende Anforderungen erfilllt : 1) Mit dem Programm kann man die Gesamtkonzentrationen der Ausgangsspezies berechnen, wenn die Konzentrationen allerfreien Metall- und Ligandenspezies bekannt sind, oder die Konzcntrationen aller freien Spezies zu berechnen, von denen die Gesamtkonzentration bekannt ist. Diese Eigenschaft des Programms ist dann von Nutzen, wenn man die Gesamtmengen berechnen will, die man in eine Lasung einsetzen mu& urn vorgeschriebene Konzentrationen an freien Spezies zu erhalten. 2) Urn die Unsicherheit da Bildungskonstanten der Komplexe zu bertlcksichtigen, erlaubt das Programm, die eingespeisten Werte der Bildungskonstanten automat&h zu variieren und die sich ergebende Xnderung der berechneten Konzentrationen zu ermitteln.
l This situation may also a& if complex species containing metabolites A-Z in a proportion different from that of the sub-complex A,. ..Z, are present, but with negligible concentrations (c$ G. Anderegg).
1016
G. A.
CUMME
R&m&-Les complexes mttalliques a&sent souvent comme substrats dans des riactions enzymatiques. Parce que leurs concentrations ne peuvent Btre mesurkes directement, elles doivent 6tre calcukes h partir des constantes de formation et des soncentrations totales de diffhentes es&es de m&w et coordinats. Un algorithme itkratif babli par Sayce a ktt utilisk pour formuler u11programme FORTRAN (In BOO),qui r&pond aux exigence suivantes. (1) Le programme peut &re fait pour calculer les concentrations totales pour les espkces apparentks si la concentration librede n’importe quelle espkce de m&al ou de coordinat est connue ou pour calculer les concentrations libres de toutes les espkes pour lesquelles les concentrations totales sont connues. Cette caractkistique du programme est utile pour le calcul des quantitis totales qui doivent 2tre introduites dans une solution afin d’obtenir les concentrations libres prescrites. (2) Pour tenir compte de l’incertitude dans les constantes de formation de l’espkce complexe le programme permet de faire varier automatiquement les valeurs d’entrtk des constantes de formation et d’bvaluer la variation rksultante des concentrations calcuks.