Mass Action Law Conjugate Representation for General Chemical Mechanisms
aa r X i v : . [ phy s i c s . c h e m - ph ] O c t Mass action law conjugate representation forgeneral chemical mechanisms
V´ıctor Fair´en ∗ Benito Hern´andez–Bermejo
Departamento de F´ısica Fundamental, Universidad Nacional de Educaci´ona Distancia. Apartado 60.141, 28080 Madrid (Spain). Email: [email protected]
Abstract
Power-law rates constitute a common approximation to the generalanalysis of the stability properties of complex reaction networks. Wepoint out in this paper that this form for the rates does not need tobe assumed as an approximation for general rate-laws. On the contrary,any functional form for a rate law can be represented exactly in termsof power-laws. Moreover, we can uniquely associate to any set of kineticequations an equivalent ‘conjugate’ representation in terms of the well-known generalized Lotka-Volterra equations, standing for what we call per capita rates, which amounts to a great simplification in terms of thestructural form of the mathematical representation of a reaction network. ∗ Author tho whom all correspondence should be addressed.1 . Introduction
Mass action law chemical kinetics, and its corresponding mathematical mod-elling, has for long been considered as a prototype in nonlinear science. Wecan all recall how such archetypical schemes, as the Brusselator, Oregonator,Schl¨ogl model (to cite a few) which have constituted the vanguard in the pio-neering years. The simplicity of the stoichiometric rules and that of the algebraicstructure of the corresponding evolution equations has made chemical kineticsa traditional point of reference in modeling within such fields as populationbiology, quantitative sociology, prebiotic evolution and other biomathematicproblems, where a system is viewed as a collection of ‘species’ interacting asmolecules do. Moreover, as emphazised by ´Erdi and T´oth, even the algebraicstructure of the evolution equations from many other fields can be convertedinto ‘chemical language’, where a formal ‘analog’ in terms of a chemical reactionnetwork is defined.The interest of this common mathematical framework provided by chemicalkinetics is not only aesthetic. It has sparked the quest of theorems which connectthe structure of the chemical reaction network with the qualitative features ofthe solutions to the corresponding differential equations. In fact, we are talkingabout the search of theorems which would permit the knowledge of behaviorsopen to the system from an identification of certain patterns in the network,and that of the associated algebraic structure of the differential equations. Andconversely, mathematical propositions which should eventually point at whichof the properties of a chemical reaction network are to be selected for the ob-tainment of a given behavior. The ultimate goal being that of a classification ofnetworks, or at least of certain of their characteristics. The accomplishment ofthis purpose would certainly yield a tool of great practical importance in mod-eling. The zero deficiency theorem , and Vol’pert’s theorem are examples ofresults associating graph properties of the reaction network with the existenceof equilibrium points.In the context of classification of networks, a significant step forward hasbeen done by Clarke with the stoichiometric network analysis, with the help ofwhich he addresses the issue of connecting the topology of a given chemical net-work with: 1) The network stability problem (necessary and sufficient conditionsfor ensuring stability of steady states); 2)
The stability diagram problem (calcu-lation of the bifurcation set of an unstable network). Ross and collaborators i in reaction j , defined atsteady state X as: κ ij ( k ) = (cid:20) ∂ log v j ( X , k ) ∂ log X i (cid:21) X = X , (1)where v j ( X , k ) is the reaction rate, dependent upon a concentration vector X ,and a set of reaction constants, k . Then, the stability of the steady state, in anetwork involving n species and r reactions, is given by the solutions of:dd t δX i = n X m =1 r X j =1 ν ij κ mj (cid:20) v j X m (cid:21) X = X δX m (2)where ν ij are the elements of the net stoichiometric matrix.In principle, the rate laws v j ( X , k ) can have any form, but in (1) and (2)we effectively assume that in practice rate laws can be approximated locally,around the steady state, by the traditional power-law, v j ( X , k ) = k j n Y i =1 X κ ij i , (3)in which the reaction orders κ ij are not restricted to to the integer values givenby the law of mass action. These non integer values have already been found insituations where the mean-field approximation does not hold, as it happensin diffusion-limited kinetics, either in disordered media with fractal structures,or on regular lattices of dimensionality smaller than the critical value, d = 2. Here, the single elementary reaction2 A −→ productsdisplays a non-integer order of reaction, which may be even greater than 3 whenoccuring on fractal ‘dust’, with spectral dimension 0 ≤ d s ≤
1. Also, an ex-tensive use of the power-law approximation has been made by Savageau andcollaborators. Starting with the observation that enzyme-kinetic rates are wellrepresented by linear relations in logarithmic space, they have generalized thisstructural pattern to the analysis of many natural systems, encapsulating their3odeling in a systematic use of a version of the power-law formalism calledS-system approach. Finally, we can mention the different solutions suggestedto solve the ‘inversion problem’: the embedding of general differential equationsinto a unified formalism in terms of stoichiometric networks with power-lawkinetics; or, more specifically, mass-action kinetics. In this context, contribu-tions from Samardzija et al. , Poland and Kowalsky have tried differentroutes for producing stoichiometric network counterparts of well known proto-typical models, as Lorenz and R¨ossler systems, or the Van der Pol oscillator.The assumption of power-law rates is at the heart of most treatments try-ing to establish a unifying mathematical framework around the concept of astoichiometric network. This systematic approach opens, as Ross and collab-orators have shown, new horizons to the chemical dynamicists, inasmuch asthe structural analogies which might be discovered will be of help in configuringan association between the structure of a chemical mechanism and the expectedbehavior.We intend in the present article to stress how general functional forms forthe rate laws may be exactly encapsulated into a power-law formalism withoutresorting to a local approximation, as in (3). We will then rewrite the resultingkinetic equations into power-law rates equations, and that will permit to showthat the evolution equations for the per capita rates, defined as v j ( X , k ) X i , are always in the form of generalized Lotka-Volterra equations, no matter whatis the particular form of the original rate equations. An equivalent ‘conjugate’network may be associated to this generalized Lotka-Volterra representation,which involves unimolecular, bimolecular and pseudounimolecular steps. Weshall discuss the properties of this transformation and show that it leads to aunique generalized Lotka-Volterra representation for a given reaction network. II. Exact equivalence to power-law rates systems
Within the power-law formalism, the kinetic equations for a given speciesinvolved in a mechanism with r reactions are:˙ X i = r X j =1 k j ν ij n Y k =1 X κ kj k , (4)4o which we shall refer, from now on, as power-law rates systems.The question now is to demonstrate that general functional relations for therate laws are amenable to an equivalent power-law form without resorting toapproximations. The procedure to do so is well known , and can be illustratedwith a simple example. Assume the following mechanism of pseudoreactions A k −→ Y (5 .a ) B k −→ X (5 .b ) X k ′ −→ B (5 .c ) X + Y k −→ products (5 .d )which constitute an early model by Degn and Harrison , to account for theoscillations in the peroxidase-oxidase reaction:2 NADH + O + 2H + −→ + + 2H OReaction (5.d) is the peroxidase enzyme catalyzed oxidation of the NADH ( Y )by dilute oxygen ( X ), which was assumed in the model to be inhibited athigh concentrations of the latter. According to Degn and Harrison the cor-responding rate law was taken to follow a Michaelis-Menten form, suggestingfor model (5), when species A and B are held constant, the following equationsin dimensionless-form: ˙ X = B − X − XY qX (6 .a )˙ Y = A − XY qX (6 .b )We now introduce the auxiliary variable Z = (1 + qX ) − , which converts ther.h.s. of (6) into polynomial form, but which calls for a supplementary equationfor that same auxiliary variable. It will be˙ Z = ∂Z∂X ˙ X , which again is polynomial provided ∂Z/∂X is already in such form. This canbe proved to be the case for smooth functions, though we shall not discuss thedemonstration here (interested readers are referred to Kerner and Hern´andez–Bermejo and Fair´en ). After elementary algebra, we find for (6)˙ X = B − X − XY Z (7 .a )5 Y = A − XY Z (7 .b )˙ Z = − qBXZ + 2 qX Z + 2 qX Y Z , (7 .c )which is written in terms of power-law rates. In going from (5) to (7) thedimensionality of the kinetic equations has been increased. The equivalencebetween these two sets of equations will be ensured if the initial condition forvariable Z is taken to be Z (0) = (1 + qX (0) ) − (we again refer to Kerner for further details). According to (7) the step in (5.d) is substituted by the‘kinetically equivalent’ set of pseudoreactions X + Y + Z −→ Z + products (8 .a ) Q + B + X + 2 Z −→ X + products (8 .b ) Q + 2 X + 2 Z −→ X + 4 Z (8 .c ) Q + 2 X + Y + 3 Z −→ X + Y + 5 Z (8 .d )The high molecularity of the pseudoreactions (8) might eventually be reducedif additional ‘auxiliary variables’ were introduced in (7), i.e.: XZ = W ; with acorresponding increase in the number of items in the kinetically equivalent setof pseudoreactions.The previous procedure can be systematically carried out for any rate-lawrepresented by a smooth function. It amounts to a reduction in the degree ofnonlinearity to power-laws by the labelling under ‘auxiliary variables’ of func-tional expressions of the original independent variables. The initial conditionsfor these new variables are then automatically prescribed by the same functionalexpressions from which they are defined. Once the procedure is complete, thenew variables are understood to represent the concentrations of some ‘virtual’species, reacting in accordance to some appropriate mechanism, as in (8).Up to this point, and before proceeding any further, a comment is needed onthe procedure leading from (5) to (8). The previous method should be viewedas a ‘protocol’ for rewriting a system of ordinary nonlinear differential equationsinto a pattern formally identifiable as one describing the evolution of a collectionof ‘objects’ interacting according to the rules of the law of mass action. We thendo refer to ‘pseudospecies’, rather than to chemical species, because no actualchemical process has been found to obey such schemes (network). This is simplydue to the fact that the combinations and scenarios open to objects behaving6nder the simple rules of the law of mass action clearly outnumber the actu-ally known chemical processes, including those which are seriously consideredchemically plausible but have not been actually obseved. III. Conjugate representation in terms of generalized Lotka-Volterraequations
We can now return to (4) and write it in a slightly different way:˙ X i = X i r X j =1 k j ν ij n Y k =1 X κ kj − δ ik k = X i m X j =1 A ij n Y k =1 X B jk k (9)where i = 1 , . . . , n and δ ik is the Kronecker delta symbol.The kinetic equations are now written in terms of the per capita rates Y j = n Y k =1 X B jk k , j = 1 , . . . , m (10)In (9) we implicitly assume that there are actually m distingishable per capita rates in a network with r reactions: m is not necessarily equal to r , for differ-ent reactions might possess the same per capita rates, and conversely, a singlereaction will generate a specific per capita rate for each of the relevant speciesinvolved.We shall henceforth assign, in (9) and (10), the label j = 1 to the constant per capita rate (that with B k = 0 , k = 1 , . . . , n ). It will then be understoodthat Y , . . . , Y m will have at least one non-null B jk entry. Then, B will be an m × n matrix with a null first row, and A an n × m matrix with its first columnfilled in with the coefficients of Y (which will correspond to the linear rates inthe network).Let us now take time derivatives of both sides of (10). For Y we have˙ Y = 0 (11)to which we can assign, without loss of generality, the solution Y = 1 (12)On the other hand, for Y , . . . , Y m , we have˙ Y j = n X s =1 d Y j d X s ˙ X s = n X s =1 B js n Y k =1 X B jk − δ ks k =7 Y j m X q =1 n X s =1 B js A sq ! n Y p =1 X B qp p = Y j m X q =1 L jq Y q (13)Then, the resulting equations of motion for the per capita rates are in the formof generalized Lotka-Volterra equations. They might themselves, in turn, beassigned to the time evolution of a set of pseudospecies in a ‘chemical network’comprising, at most, bimolecular steps in the relevant pseudospecies. This resultis universal and applicable to any set of kinetic equations with power-law rates,and by extension, as demonstrated before, it is also valid for general functionalrelations for the rates.The generalized Lotka-Volterra equations (13) constitute a conjugate rep-resentation of the original kinetic equations (4) for a given chemical network.They provide a much stronger unifying structure than that associated to (4);a structure for which there exist several tools for studying the features of thesolutions, and which are straightforwardly translatable into a graph theorysetting. This is not the place to discourse upon the issues related to the gener-alized Lotka-Volterra equations, because most of the general results of interestto the chemical dynamicist (related to stability of steady states) are highlymathematical and fall outside the scope of this paper. We shall neverthelesstouch upon some aspects of interest later on, after giving an example of theprocedure of obtainment of (13). Example:
Let us take (7). According to the notation in (9) we can writefor matrices A and B : A = − B − A − − qB q q (14) B = − − , (15)from which we can easily calculate matrix L = B · A , in (13). B is simply theorder matrix for the per capita rates involved in (7), and can be systematically8ritten down once we display (7) in the following form:ddt ln( X ) = ˙ XX = − BX − − Y Z (16 .a )ddt ln( Y ) = ˙ YY = AY − − XZ (16 .b )ddt ln( Z ) = ˙ ZZ = − qBXZ + 2 qX Z + 2 qX Y Z , (16 .c )whilst matrix A is obtained from the corresponding coefficients. As indicatedbefore, the first row in (15) has zero entries. Correspondingly, the general-ized Lotka-Volterra matrix L = B · A will also have a zero first row, which isunderstandably assignable to (11).A mechanism (or network) may be associated to the generalized Lotka-Volterra equations in terms of unimolecular, pseudounimolecular and bimolec-ular steps which follow the law of mass-action. This mechanism is a sort of‘conjugate’ reactional scheme to the original network, and might be as wellbe used as a complementary (or alternative) representation in investigations di-rected torwards the classification of complex reaction networks. The structuralsimplicity of the generalized Lotka-Volterra equations (which are characterizedby a single algebraic object: matrix L , the properties of which are strongly re-lated to certain important features of the solutions), as well as their ubiquity inmany scientific disciplines, makes this approach particularly attractive.Upon examination of matrix L obtained from, say, (14) and (15), we inferthat its last three rows will be linear combinations of the three immediatelypreceding ones. By construction, this fact is generalizable to any matrix L ( m − n − n rows). This pattern determinesthe establishment of an associated modular mass-action law reaction network.A pseudoreaction template is ascribed to each entry of the independent rows inmatrix L . Those pseudoreactions templates constitute the modules, or buildingblocks, which the whole reaction network is made of.In connection to the example of the peroxidase-oxidase model (7), theseconstitutive units are, in view of (14) and (15), given by ± Y j ( − − BY + Y ) (17 .a ) ± Y j ( AY − (1 + 2 qB ) Y + 2 qY + 2 qY ) (17 .b )9 Y j ( − AY + Y ) (17 .c )for any j ≥ Y j −→ Y j + B + Y −→ B + Y Y j + Y −→ Y j + Y with (–) sign Y j −→ Y j Y j + B + Y −→ Y j + B + Y Y j + Y −→ Y IV. Properties of the transformation to the conjugate representation
In order to demonstrate some important results regarding the validity andscope of the previous manipulations, we will consider in this section the mostusual case in which m > n , that is, the number of per capita rates is greaterthan that of variables. For example, in the peroxidase equations (16) we have m = 7 and n = 3. We will also assume that the rank of matrix B is maximum:rank( B ) = n .A necessary condition for ensuring the equivalence between system (9) andthe generalized Lotka-Volterra equations (13) is that the transformation relat-ing them preserves the topological characteristics of the solutions. We shallprove that this is indeed the case here. A sufficient condition for demonstratingthis statement is the existence of a continuous, differentiable and invertibleapplication connecting the initial and final phase spaces. Since the dimensionof the generalized Lotka-Volterra system is greater than that of (9), such anapplication should connect the phase space of (9) and the n -dimensional subsetof R m into which it is mapped.We can write the transformation (10) relating the power-law rates systemvariables to the per capita rates as: Y j = m Y k =1 X ˜ B jk k , j = 1 , . . . , m , (18)10here X n +1 = . . . = X m = 1 and ˜ B is an m × m matrix, defined as:˜ B = B (cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12) B ′ ,n +1 . . . B ′ ,m ... ... B ′ m,n +1 . . . B ′ m,m ≡ ( B | B ′ ) (19)Here B is the already known m × n matrix of exponents of the per capita rates, and B ′ is a m × ( m − n ) matrix of arbitrary entries. These entries canalways be selected in such a way that ˜ B is invertible. Equation (18) is obviouslydifferentiable. Thus, we only need to prove that it is one to one and invertible.If we take logarithms in both sides of (18): ln( Y )...ln( Y m ) = ˜ B ln( X )...ln( X n )0...0 (20)Since rank( ˜ B )= m , for any two vectors X and X ′ we have ˜ B ln( X ) = ˜ B ln( X ′ ),unless X = X ′ . Thus the map (18) is one to one and invertible ( ˜ B invertible)and the topology is preserved by the transformation.The original variables of the power-law rates system (9) can be retrievedfrom those of the generalized Lotka-Volterra system by means of two differentprocedures. The first one is obtained by writing system (9) in the separableform: ˙ X i X i = m X j =1 A ij Y j ( t ) , i = 1 , . . . , n (21)Then the X i result from the formal integrations: X i ( t ) = X i (0) exp m X j =1 A ij Z t Y j ( t ′ )d t ′ (22)The second approach is purely algebraic and does not require any integration.Since B is an m × n matrix, with m > n , and rank( B ) is maximum, then thereexists an n × n invertible submatrix B n of B . Let B n = B i . . . B i n ... ... B i n . . . B i n n , { i , . . . , i n } ⊂ { , . . . , m } (23)11his implies that: ln( Y i )...ln( Y i n ) = B n ln( X )...ln( X n ) (24)Since B n is invertible, this finally leads to: X k ( t ) = n Y p =1 (cid:2) Y i p ( t ) (cid:3) ( B − n ) kp , k = 1 , . . . , n (25)The time evolution or stability properties of a given reaction network might beanalyzed in anyone of these two alternative descriptions ( X or Y ), for they arecompletely equivalent. However, as far as structural properties are concerned,the generalized Lotka-Volterra form (in terms of Y ) seems preferable for it ismathematically characterized by a single algebraic object: matrix L .As we have seen, to every power-law rates system (9) a single generalizedLotka-Volterra conjugate system can be associated. The question now is to whatextent this is also valid for general rate-laws:˙ X = F ( X ) (26)The way for finding the conjugate representation consists, as we saw in SectionsII and III in the peroxidase-oxidase example, in the introduction of auxiliaryvariables for functional rate-laws in the right hand side of (26) not comply-ing to the power-law rates system format. This always leads to a power-lawrates system from which the obtainment of the generalized Lotka-Volterra sys-tem is straightforward. We shall see that, to a great extent, the generalizedLotka-Volterra representative is unique for every system of the form (26), andis independent of the specific choice of auxiliary variables. Instead of a formalapproach, we shall consider in more detail the peroxidase example. However,the results that we shall display can be proved rigorously. , Let us generalize the procedure of Section II by introducing an auxiliaryvariable of the form: Z = X α Y β (1 + qX ) γ , (27)where α , β and γ are real parameters and γ = 0. After some algebra, theintroduction of this general variable leads to a family of ( α, β, γ )-dependent12ower-law rates systems with matrices: A ( α, β, γ ) = − B − A − − α Bα − α Aβ qBγ − β − qγ − qγ (28) B ( α, β, γ ) = − α/γ β/γ − /γ − α/γ β/γ − /γ α/γ β/γ − /γ α/γ β/γ − /γ (29)However, the product L = B ( α, β, γ ) · A ( α, β, γ ) is independent of α , β and γ . Since L is the matrix associated to the conjugate generalized Lotka-Volterrasystem, this means that such representation is unique, independently of thechoice of the auxiliary variables. Of course, this matrix L coincides with theone obtained from the product of (15) and (14), which are particular cases of(29) and (28), respectively, with α = 0 , β = 0 , γ = −
1. This is consistent withthe fact that the variables of the generalized Lotka-Volterra representative areindependent of ( α, β, γ ). From matrix (29), they are:1 ; X − ; X α/γ Y β/γ Z − /γ = Y qX ; Y − ; X α/γ Y β/γ Z − /γ = X qX ; X α/γ Y β/γ Z − /γ = X qX ; X α/γ Y β/γ Z − /γ = X Y (1 + qX ) (30)This implies that the initial conditions of the conjugate generalized Lotka-Volterra system will also be unique. In other words: To every general systemof the form (26) a single generalized Lotka-Volterra system can be associatedby means of this procedure.
Although the process leads to an infinite family ofintermediate power-law rates systems, all of them possess the same generalizedLotka-Volterra representative, irrespective of the parameters α , β and γ , andare thus all equivalent. This property supports our assertions in favor of thegeneralized Lotka-Volterra as a unifying format.13 . Conclusions We have stressed how the power-law formalism can be a referential formatfor general functional forms for chemical rate-laws. On encapsulating a chemicalmechanism (or network) under a power-law formalism, there is no need to resort,as we have shown, to any kind of local approximation in terms of that samepower-law formalism, even if it seems justified experimentally. Instead, simplemanipulations of elementary calculus convert non polynomial kinetic equationsinto power-law differential equations, completely equivalent to the original oneswhen appropriate initial conditions are assumed.For power-law rates an interesting universal relationship has been obtained.When these rate laws are considered as per capita rates (or, equivalently, interms of logarithmic derivatives) they obey a set of generalized Lotka-Volterraequations. The specific matrix characterizing this generalized Lotka-Volterrasystem is independent of the particular embedding procedure when transforminggeneral rate laws into a power-law formalism. Also, to each particular power-lawrates system corresponds a unique and mathematically equivalent generalizedLotka-Volterra system. The latter may then be considered a conjugate repre-sentation of any chemical network.Much attention has been devoted in the literature to the generalized Lotka-Volterra equations, a fact which is not independent of their structural simplicityand their ubiquity in many scientific disciplines, ranging from population biologyto laser physics. This is particularly attractive in the context of classificationof chemical networks.A conjugate chemical network , with at most bimolecular steps in the essentialspecies, may be associated to the generalized Lotka-Volterra equations. Thenetwork is purely conceptual and should not be thought of as the canonical reactional scheme of any chemical process. It should be regarded as an abstractequivalent representation of a model system in the familiar language of mass-action kinetics. Its inmediate interest in the modeling of actual chemical systemsmay be presently a subject of debate, for many critics argue that the fieldof chemical network dynamics has not yet produced any result of chemicalimportance. This point of view should be seriously reconsidered in the light ofrecent work by Ross and collaborators. cknowledgements: This work has been supported by the DGICYT(Spain), under grant PB94-0390. B. H. acknowledges a doctoral fellowship fromComunidad Aut´onoma de Madrid.
References and notes ´Erdi, P.; T´oth, J. Mathematical Models of Chemical Reactions;
ManchesterUniversity Press: Manchester, 1989; pp. 1-13. Pielou, E. C.
Mathematical Ecology;
John Wiley & Sons: New York, 1977. Weidlich, W.; Haag, G.
Concepts and Models of a Quantitative Sociology;
Springer-Verlag: Berlin, 1983. Hofbauer, J.; Sigmund, K.
The Theory of Evolution and Dynamical Systems;
Cambridge University Press: Cambridge, 1988. Murray, J.D.
Mathematical Biology, See the articles by Othmer, H.G. and Feinberg, M. in
Modelling of ChemicalReaction Systems;
Ebert, K.H., Deuflhard, P., J¨ager W., Eds.; Springer-Verlag:Berlin, 1981. Feinberg, M.
Arch. Ratl. Mech. Anal. , , 1. Horn, F.; Jackson, R.
Arch. Ratl. Mech. Anal. , , 81. See reference 1, pp. 45-48. Clarke, B. L.
Adv. Chem. Phys. , , 1. Eiswirth, M.; Freund, A.; Ross, J.
Adv. Chem. Phys. , , 127. Cheva-lier, T.; Schreiber, I.; Ross, J. J. Phys. Chem. , , 6776. Hung, Y. F.;Ross, J. J. Phys. Chem. , , 1974. Hung, Y. F.; Schreiber, I.; Ross, J. J.Phys. Chem. , , 80. Stemwedel, J. D.; Ross, J. J. Phys. Chem. , , 1988. Argyrakis, P. In
Fractals, Quasicrystals, Chaos, Knots and Algrebraic Quan-tum Mechanics;
Amann, A., Cederbaum, L., Gans, W., Eds.; Kluwer: NewYork, 1988; p. 53. Klymko, P. W.; Kopelman, R.
J. Phys. Chem. , , 3686. Klymko, P.W.; Kopelman, R. J. Phys. Chem. , , 4565. Anacker, L. W.; Kopelman,R. J. Chem. Phys. , , 6402. Kopelman, R. J. Stat. Phys. , ,185. Lin, A.; Kopelman, R.; Argyrakis, P. Phys. Rev. E , , 1502. See Anacker and Kopelman in reference 13.15 Savageau, M. A.
J. Theor. Biol. , , 365. Voit, E. O., Ed.
Canonical Nonlinear Modelling: S-system Approach to Un-derstanding Complexity;
Van Nostrand: New York, 1991. See reference 1, p. 64. Samardzija, N.; Greller, L. D.; Wasserman, E.
J. Chem. Phys. , ,2296. Poland, D.
Physica D , , 86. Kowalski, K.
Chem. Phys. Lett. , , 167. Kerner, E. H.
J. Math. Phys. , , 1366. Hern´andez–Bermejo, B.; Fair´en, V.
Phys. Lett. A , , 31. Degn, H.; Harrison, D. E. F.
J. Theoret. Biol. , , 238. Fair´en, V.; Velarde, M. G.
J. Math. Biol. , , 147. Takeuchi, Y.; Adachi, N.; Tokumaru, H.
Math. Biosci. , , 119. Peschel, M.; Mende, W.
The Predator–Prey Model. Do we live in a VolterraWorld?
Springer-Verlag: Wien–New York, 1986. Jackson, E. A.
Perspectives of Nonlinear Dynamics,
Vol. 1; 1st ed.; Cam-bridge University Press: Cambridge, 1994; pp. 21-23. Hern´andez–Bermejo, B.; Fair´en, V.