# Complex Reaction Kinetics in Chemistry: A Unified Picture Suggested by Mechanics in Physics

Authors: Elena Agliari; Adriano Barra; Giulio Landolfi; Sara Murciano; Sarah Perrone
Journal: Complexity
 (2018)

---

## Abstract

Complex biochemical pathways can be reduced to chains of elementary reactions, which can be described in terms of chemical kinetics. Among the elementary reactions so far extensively investigated, we recall the Michaelis-Menten and the Hill positive-cooperative kinetics, which apply to molecular binding and are characterized by the absence and the presence, respectively, of cooperative interactions between binding sites. However, there is evidence of reactions displaying a more complex pattern: these follow the positive-cooperative scenario at small substrate concentration, yet negative-cooperative effects emerge as the substrate concentration is increased. Here, we analyze the formal analogy between the mathematical backbone of (classical) reaction kinetics in Chemistry and that of (classical) mechanics in Physics. We first show that standard cooperative kinetics can be framed in terms of classical mechanics, where the emerging phenomenology can be obtained by applying the principle of least action of classical mechanics. Further, since the saturation function plays in Chemistry the same role played by velocity in Physics, we show that a relativistic scaffold naturally accounts for the kinetics of the above-mentioned complex reactions. The proposed formalism yields to a unique, consistent picture for cooperative-like reactions and to a stronger mathematical control.


---

## Body

## 1. Introduction

### 1.1. The Chemical Kinetics Background

The mathematical models that describe reaction kinetics provide chemists and chemical engineers with tools to better understand, depict, and possibly control a broad range of chemical processes (see, e.g., [1, 2]). These include applications to pharmacology, environmental pollution monitoring, and food industry. In particular, biological systems are often characterized by complex chemical pathways whose modeling is rather challenging and can not be recast in standard schemes [3–15] (see also [16–19] for a different perspective). In general, one tries to split such sophisticated systems into a set of elementary constituents, in mutual interaction, and for which a clear formalization is available [20–25].In this context, one of the best consolidated, elementary scheme is given by the Michaelis-Menten law. This was originally introduced by Leonor Michaelis and Maud Menten to describe enzyme kinetics and can be applied to systems made of two reactants, sayA (the binding molecule or, more generally, the binding sites of a molecule) and B (the free ligand, i.e., the substrate), which can bind (and unbind) to form the product AB. If we call S the concentration of free ligand, Y thesaturation function (orfractional occupancy), namely, the fraction of bound molecules (Y∈[0,1]), and, accordingly, 1-Y the fraction of the unbound molecules, under proper assumptions, one can write (1)S1-Y=kY,where k is the proportionality constant between response and occupancy (otherwise stated, it is the ratio between the dissociation and the association constants). In particular, as standard, it is assumed that (a) the reaction is in a steady state, with the product being formed and consumed at the same rate, (b) the free ligand concentration is in large excess over that of the binding molecules in such a way that it can be considered as constant along the reaction, and (c) all the binding molecules are equivalent and independent. Also, the derivation of the Michaelis-Menten law is based on the law of mass action.By reshuffling the previous equation we getY=S/(S+k) which allows stating that k is the concentration of free ligand at which 50% of the binding sites are occupied (i.e., when S=k, then Y=1/2). Thus, denoting with S0 the half-saturation ligand concentration, we get (2)Y=SS+S0.This equation represents a rectangular hyperbola with horizontal asymptote corresponding to full saturation; that is, Y=1; this is the typical outcome expected for systems where no interaction between binding sites is at work [26]. This model immediately settled down as the paradigm for Chemical Kinetics, somehow similarly to the perfect gas model (where atoms, or molecules, collisions apart, do not interact) of the Kinetic Theory in the early Statistical Physics [27]. Nevertheless, deviations from this behaviour were not late to arrive: the most common phenomenon was the occurrence of a positive cooperation among the binding sites of a multisite molecule. Actually, many polymers and proteins exhibit cooperativity, meaning that the ligand binds in a nonindependent way: if, upon a ligand binding, the probability of further binding (by other ligands) is enhanced, the system is said to showpositive cooperativity.To fix ideas, let us make a practical example and let us consider the case of a well-known protein, that is, the hemoglobin. This is responsible of oxygen transport throughout the body and it ultimately allows cellular respiration. Such features stem from hemoglobin’s ability to bind (and to dislodge as needed) up to four molecules of oxygen in a nonindependent way: if one of the four sites has captured an oxygen molecule, then the probability that the remaining three sites will capture further oxygen increases, and vice versa. As a result, if the protein is in an environment rich of oxygen (e.g., in the lungs), it readily binds up to four molecules of oxygen, and, as much readily, it gets rid of them when crossing an oxygen-deficient environment. To study quantitatively its behaviour one typically measures its characteristic input-output relation. This can be achieved by considering a set ofM elementary experiments where these proteins, in the same amount for each experiment, are prepared in a baker and allowed to bind oxygen, which is supplied at different concentrations Si for different experiments (e.g., S1<S2<⋯<SM). We can then construct a Cartesian plane, where on the abscissas we set the concentration of the ligand S (oxygen in this case, i.e., theinput) while on the y-axes we put the fraction of protein bound sites Y (the saturation function, i.e., theoutput). In this way, for each experiment, once reached the chemical equilibrium, we get a saturation level and we can draw a point in the considered Cartesian plane; interpolating between all the points a sigmoidal curve will emerge (see Figure 1). Archibald V. Hill formulated a description for the behavior of Y with respect to S: the so-called Hill equation empirically describes the fraction of molecules binding sites, occupied by the ligand, as a function of the ligand concentration [28–31]. This equation generalizes the Michaelis-Menten law (2) and reads as (3)Y=kSnHS0+SnH,where nH is referred to as Hill coefficient and can be interpreted as the effective number of binding sites that are interacting with each other. This number can be measured as the slope of the curve log[Y/(1-Y)] versus log(S), calculated at the half-saturation point. Of course, if nH=1 there is no cooperation at all and each binding site acts independently of the others (and, consistently, Michaelis-Menten kinetics is restored), and vice versa; if nH>1, the reaction is said to be cooperative (just like in hemoglobin), and if nH≫1 the cooperation among binding sites is so strong that the sigmoid becomes close to a step function and the kinetics is named ultrasensitive.Figure 1
These plots show comparison between data from experiments (symbols) and best-fits through (22) (lines). Data refer to noncooperative and positive-cooperative systems [33, 34] (a) and an ultrasensitive system [35] (b). For the latter we report two fits: dashed line is the result obtained by constraining the system to be cooperative but not ultrasensitive (i.e., J≤1), while solid line is the best fit (without constraints) which yields to J~1.1, hence a “first-order phase transition” in the language of statistical mechanics. The relative goodness of the fits is Rcoop2~0.85 and Rultra2~0.94, confirming an ultrasensitive behavior. The tables in the bottom present the value of J derived from the best fit and the resulting nH according to (24); the estimate of the Hill coefficient taken from the literature is also shown for comparison. This figure was presented in [28].
(a)
(b)The Michaelis-Menten law, together with the extension by Hill, provided a good description for a bulk of chemical reactions; however, things were not perfect yet. For instance, some yeast’s proteins (e.g., the glyceraldehyde 3-phosphate dehydrogenase [32]) produced novel (mild) deviations from the Hill curve: for these enzymes, the cooperativity of their binding sites decreases while increasing the ligand concentration. The following work by Daniel E. Koshland allowed understanding this kind of phenomenology by further enlarging the theoretical framework through the introduction of the concept of negative cooperativity. In fact, in the previous example, beyond the positive cooperation between the binding sites there are also negative-cooperative effects underlying. Their effective action is to diminish the overall binding capabilities of the enzyme and thus to reduce the magnitude of its Hill coefficient.


### 1.2. The Mechanics Background

The progressive enlargement of a theoretical scaffold to fit the always increasing amount of evidences is a common feature in the historical development of scientific disciplines [36, 37]. This is the case also for Mechanics and, as we will see, the analogy with Chemical Kinetics goes far beyond this feature.Beyond Kinematics, which describes the motion of systems without considering their mass or the forces that caused the motion, in the seventeenth century Newton gave a sharp description of Mechanics, in the form of laws describing how masses dynamically respond when stimulated by an external force (or moment). Here, the input is the force while the output is the motion of the body. The Newtonian dynamics has been ruling for centuries and, in fact, it was so well-consolidated that scientists, among which Giuseppe L. Lagrange, William R. Hamilton, and Carl G. J. Jacobi, later reformulated the entire theory in a powerful and elegant variational flavor. The theory was overall brilliant to explain the perceivable reality, but with exceptions emerging in the limit oftoo little ortoo fast.We will focus on the latter. In the Newtonian world, if an applied force is kept constant over a mass, this will constantly accelerate, eventually reaching diverging velocities. This was perfectly consistent with the general credo that the speed of light was infinite. However, this postulate broke down in 1887 when the famous experiment by Albert A. Michelson and Edward Morley proved that such a velocity is actually finite. The next years were dense of novel approaches and ideas by many scientists, as Hendrik Lorentz and Hermann Minkowski, and culminated with the special relativity by Albert Einstein in 1905. According to this theory, no mass can exist whose velocity may diverge, the limiting speed being the speed of light. The classical Hamilton-Jacobi equations and Galilean transformations left the place to the Klein-Gordon formulation and Lorentz covariances and contravariances (the natural metric being Minkowskian) [38]. Clearly, classical mechanics were still a good reference framework for the vast majority of the data collected (much like the positive cooperativity accounted for the bulk of the empirical data in the chemical counterpart); however, there were rare phenomena (e.g., a muon decay in atmosphere [39]) that required a broader scaffold which, in the opportune limits, could recover the classical one.Although this historical connection between Chemical Kinetics and Classical Mechanics may look weird at a first glance, as we will prove, there is a formal analogy between their mathematical representations. In the next section we will summarize the main results concerning the analogy at the classical level. More sharply, the saturation plot of classical (positive-cooperative) chemical kinetics (namely, the input-output relation between the saturation function and the concentration of the substrate) can be derived by a minimum action principle that is the same that holds in classical mechanics, when describing a mass motion in the Hamilton-Jacobi framework. In this parallelism, the saturation function in Chemistry plays as the velocity in Physics: thus, exactly as what happens in special relativity, the velocity of the mass is bounded (by definition, the saturation function can not exceed one). Indeed, we can follow this mathematical equivalence and verify that there is actually a natural broader formulation for chemical kinetics that is exactly through the Klein-Gordon setting (rather than its classical Hamilton-Jacobi counterpart) and the theory as a whole is Lorentz-invariant. Remarkably, when read with chemical glasses, this extended relativistic setting allows for the anticooperative corrections that Koshland revealed in the study of the yeast enzymes, resulting in a complex mixture of positive and negative cooperation among binding sites.


## 1.1. The Chemical Kinetics Background

The mathematical models that describe reaction kinetics provide chemists and chemical engineers with tools to better understand, depict, and possibly control a broad range of chemical processes (see, e.g., [1, 2]). These include applications to pharmacology, environmental pollution monitoring, and food industry. In particular, biological systems are often characterized by complex chemical pathways whose modeling is rather challenging and can not be recast in standard schemes [3–15] (see also [16–19] for a different perspective). In general, one tries to split such sophisticated systems into a set of elementary constituents, in mutual interaction, and for which a clear formalization is available [20–25].In this context, one of the best consolidated, elementary scheme is given by the Michaelis-Menten law. This was originally introduced by Leonor Michaelis and Maud Menten to describe enzyme kinetics and can be applied to systems made of two reactants, sayA (the binding molecule or, more generally, the binding sites of a molecule) and B (the free ligand, i.e., the substrate), which can bind (and unbind) to form the product AB. If we call S the concentration of free ligand, Y thesaturation function (orfractional occupancy), namely, the fraction of bound molecules (Y∈[0,1]), and, accordingly, 1-Y the fraction of the unbound molecules, under proper assumptions, one can write (1)S1-Y=kY,where k is the proportionality constant between response and occupancy (otherwise stated, it is the ratio between the dissociation and the association constants). In particular, as standard, it is assumed that (a) the reaction is in a steady state, with the product being formed and consumed at the same rate, (b) the free ligand concentration is in large excess over that of the binding molecules in such a way that it can be considered as constant along the reaction, and (c) all the binding molecules are equivalent and independent. Also, the derivation of the Michaelis-Menten law is based on the law of mass action.By reshuffling the previous equation we getY=S/(S+k) which allows stating that k is the concentration of free ligand at which 50% of the binding sites are occupied (i.e., when S=k, then Y=1/2). Thus, denoting with S0 the half-saturation ligand concentration, we get (2)Y=SS+S0.This equation represents a rectangular hyperbola with horizontal asymptote corresponding to full saturation; that is, Y=1; this is the typical outcome expected for systems where no interaction between binding sites is at work [26]. This model immediately settled down as the paradigm for Chemical Kinetics, somehow similarly to the perfect gas model (where atoms, or molecules, collisions apart, do not interact) of the Kinetic Theory in the early Statistical Physics [27]. Nevertheless, deviations from this behaviour were not late to arrive: the most common phenomenon was the occurrence of a positive cooperation among the binding sites of a multisite molecule. Actually, many polymers and proteins exhibit cooperativity, meaning that the ligand binds in a nonindependent way: if, upon a ligand binding, the probability of further binding (by other ligands) is enhanced, the system is said to showpositive cooperativity.To fix ideas, let us make a practical example and let us consider the case of a well-known protein, that is, the hemoglobin. This is responsible of oxygen transport throughout the body and it ultimately allows cellular respiration. Such features stem from hemoglobin’s ability to bind (and to dislodge as needed) up to four molecules of oxygen in a nonindependent way: if one of the four sites has captured an oxygen molecule, then the probability that the remaining three sites will capture further oxygen increases, and vice versa. As a result, if the protein is in an environment rich of oxygen (e.g., in the lungs), it readily binds up to four molecules of oxygen, and, as much readily, it gets rid of them when crossing an oxygen-deficient environment. To study quantitatively its behaviour one typically measures its characteristic input-output relation. This can be achieved by considering a set ofM elementary experiments where these proteins, in the same amount for each experiment, are prepared in a baker and allowed to bind oxygen, which is supplied at different concentrations Si for different experiments (e.g., S1<S2<⋯<SM). We can then construct a Cartesian plane, where on the abscissas we set the concentration of the ligand S (oxygen in this case, i.e., theinput) while on the y-axes we put the fraction of protein bound sites Y (the saturation function, i.e., theoutput). In this way, for each experiment, once reached the chemical equilibrium, we get a saturation level and we can draw a point in the considered Cartesian plane; interpolating between all the points a sigmoidal curve will emerge (see Figure 1). Archibald V. Hill formulated a description for the behavior of Y with respect to S: the so-called Hill equation empirically describes the fraction of molecules binding sites, occupied by the ligand, as a function of the ligand concentration [28–31]. This equation generalizes the Michaelis-Menten law (2) and reads as (3)Y=kSnHS0+SnH,where nH is referred to as Hill coefficient and can be interpreted as the effective number of binding sites that are interacting with each other. This number can be measured as the slope of the curve log[Y/(1-Y)] versus log(S), calculated at the half-saturation point. Of course, if nH=1 there is no cooperation at all and each binding site acts independently of the others (and, consistently, Michaelis-Menten kinetics is restored), and vice versa; if nH>1, the reaction is said to be cooperative (just like in hemoglobin), and if nH≫1 the cooperation among binding sites is so strong that the sigmoid becomes close to a step function and the kinetics is named ultrasensitive.Figure 1
These plots show comparison between data from experiments (symbols) and best-fits through (22) (lines). Data refer to noncooperative and positive-cooperative systems [33, 34] (a) and an ultrasensitive system [35] (b). For the latter we report two fits: dashed line is the result obtained by constraining the system to be cooperative but not ultrasensitive (i.e., J≤1), while solid line is the best fit (without constraints) which yields to J~1.1, hence a “first-order phase transition” in the language of statistical mechanics. The relative goodness of the fits is Rcoop2~0.85 and Rultra2~0.94, confirming an ultrasensitive behavior. The tables in the bottom present the value of J derived from the best fit and the resulting nH according to (24); the estimate of the Hill coefficient taken from the literature is also shown for comparison. This figure was presented in [28].
(a)
(b)The Michaelis-Menten law, together with the extension by Hill, provided a good description for a bulk of chemical reactions; however, things were not perfect yet. For instance, some yeast’s proteins (e.g., the glyceraldehyde 3-phosphate dehydrogenase [32]) produced novel (mild) deviations from the Hill curve: for these enzymes, the cooperativity of their binding sites decreases while increasing the ligand concentration. The following work by Daniel E. Koshland allowed understanding this kind of phenomenology by further enlarging the theoretical framework through the introduction of the concept of negative cooperativity. In fact, in the previous example, beyond the positive cooperation between the binding sites there are also negative-cooperative effects underlying. Their effective action is to diminish the overall binding capabilities of the enzyme and thus to reduce the magnitude of its Hill coefficient.


## 1.2. The Mechanics Background

The progressive enlargement of a theoretical scaffold to fit the always increasing amount of evidences is a common feature in the historical development of scientific disciplines [36, 37]. This is the case also for Mechanics and, as we will see, the analogy with Chemical Kinetics goes far beyond this feature.Beyond Kinematics, which describes the motion of systems without considering their mass or the forces that caused the motion, in the seventeenth century Newton gave a sharp description of Mechanics, in the form of laws describing how masses dynamically respond when stimulated by an external force (or moment). Here, the input is the force while the output is the motion of the body. The Newtonian dynamics has been ruling for centuries and, in fact, it was so well-consolidated that scientists, among which Giuseppe L. Lagrange, William R. Hamilton, and Carl G. J. Jacobi, later reformulated the entire theory in a powerful and elegant variational flavor. The theory was overall brilliant to explain the perceivable reality, but with exceptions emerging in the limit oftoo little ortoo fast.We will focus on the latter. In the Newtonian world, if an applied force is kept constant over a mass, this will constantly accelerate, eventually reaching diverging velocities. This was perfectly consistent with the general credo that the speed of light was infinite. However, this postulate broke down in 1887 when the famous experiment by Albert A. Michelson and Edward Morley proved that such a velocity is actually finite. The next years were dense of novel approaches and ideas by many scientists, as Hendrik Lorentz and Hermann Minkowski, and culminated with the special relativity by Albert Einstein in 1905. According to this theory, no mass can exist whose velocity may diverge, the limiting speed being the speed of light. The classical Hamilton-Jacobi equations and Galilean transformations left the place to the Klein-Gordon formulation and Lorentz covariances and contravariances (the natural metric being Minkowskian) [38]. Clearly, classical mechanics were still a good reference framework for the vast majority of the data collected (much like the positive cooperativity accounted for the bulk of the empirical data in the chemical counterpart); however, there were rare phenomena (e.g., a muon decay in atmosphere [39]) that required a broader scaffold which, in the opportune limits, could recover the classical one.Although this historical connection between Chemical Kinetics and Classical Mechanics may look weird at a first glance, as we will prove, there is a formal analogy between their mathematical representations. In the next section we will summarize the main results concerning the analogy at the classical level. More sharply, the saturation plot of classical (positive-cooperative) chemical kinetics (namely, the input-output relation between the saturation function and the concentration of the substrate) can be derived by a minimum action principle that is the same that holds in classical mechanics, when describing a mass motion in the Hamilton-Jacobi framework. In this parallelism, the saturation function in Chemistry plays as the velocity in Physics: thus, exactly as what happens in special relativity, the velocity of the mass is bounded (by definition, the saturation function can not exceed one). Indeed, we can follow this mathematical equivalence and verify that there is actually a natural broader formulation for chemical kinetics that is exactly through the Klein-Gordon setting (rather than its classical Hamilton-Jacobi counterpart) and the theory as a whole is Lorentz-invariant. Remarkably, when read with chemical glasses, this extended relativistic setting allows for the anticooperative corrections that Koshland revealed in the study of the yeast enzymes, resulting in a complex mixture of positive and negative cooperation among binding sites.


## 2. TheStandard Mathematical Scaffold for Classical Cooperativity

As anticipated in Section1.1, cooperativity is a widespread phenomenon in Chemistry and its underlying mechanisms can be multiple: for example, if the adjacent binding sites of a protein can accommodate charged ions, the attraction/repulsion between the ions themselves may result in a positive/negative kinetics; in most common cases, the bonds with the substrate modify the protein conformational structure, by influencing possible further links in an allosteric way [21, 40]. Whatever the origin, cooperativity in Chemistry is a typicalemergent property that directly relates the microscopic description of a system at the single binding-site level, with the macroscopic properties shown by its constituent molecules, cells, and organisms; thus the use of Statistical Physics for its investigation appears quite natural [26, 28]. Usually, in Statistical Physics one is provided with (inverse) temperature β and with Hamiltonian (i.e., a cost-function) H(σ,J,h) describing the model at the microscopic level, namely, in terms of elementary variables σi,σj, couplings among elementary variables Jij and external fields hi acting over these. The goal is to obtain thefree energyA(β,J,h) of the model, from which the average value of the macroscopic observables can be directly derived [26].


### 2.1. Formulation of the Problem: The Thermodynamical Free Energy

In the following we summarize the minimal assumptions needed when modelling chemical kinetics from the Statistical Physics perspective; for a more extensive treatment of this kind of modelling we refer to [21, 26, 28, 29, 41], while for a rigorous explanation of the underlying equivalence between Statistical Mechanics and Analytical Mechanics we refer to the seminal works by Guerra [42], dealing with the Sherrington-Kirkpatrick model (and then deepened in, e.g., [43–46]), and by Brankov and Zagrebnov in [47], dealing with the Husimi-Temperley model (and then deepened in, e.g., [48–51]).(i)
Each binding site may or may not be occupied by a ligand: this allows us to code its state (empty versus full) by a Boolean variable. For the genericith site, we will use an Ising spin σi=±1, where σi=-1 represents an empty ith site, and vice versa; σi=+1 means that the ith site is occupied. Clearly, if there are overall N binding sites, i∈(1,…,N).(ii)
It is rather inconvenient (and ultimately unnecessary) to deal with the whole setσi,…,σN if we are interested in the properties of large numbers of these variables (i.e., in the so-calledthermodynamic limit corresponding to N→∞). If we want to distinguish between a fully empty state σi=-1∀i∈(1,…,N) (ordered case), a fully occupied state σi=+1∀i∈(1,…,N) (ordered case), and a completely random case where σi=±1 with equal probability (disordered case), it is convenient to introduce theorder parameter for these variables as themagnetizationM (this term stems from the original application of the Statistical Mechanics model in the context of magnetism) that reads as the arithmetic average of the spin state, namely, (4)M=1N∑i=1Nσi∈-1,+1.There is a univocal relation between the magnetization in Physics and the saturation function in Chemistry, where, we recall, we denote withY∈[0,1] the fractional occupation of the binding sites. In fact, one has [28, 29] (5)Y=121+M.Equation (5) constitutes the first bridge between the Chemistry we aim to describe (via the saturation functionY) and the Physics that we want to use (via the magnetizationM).(iii)
All the binding sites interact with the ligand by the same strength. This is a standard assumption in Chemical Kinetics [29, 31, 52] and it means that the diffusion of the ligands is fast enough to ensure a homogeneous solution. The concentration of free ligands is mapped into a one-body contribution H1 in the cost-function. This term encodes for the action of anexternal magnetic field in such a way that, if the field acting on ith is positive, the spin will tend to align upwards (namely, this direction is energetically favored), and vice versa. This homogenous mixing assumption translates into a homogeneous external field h, and the related contribution reads as (6)H1σ,h=-h∑i=1Nσi.Notice thath plays as a chemical potential and, consistently, it can be related to the substrate concentration as (7)h=12log⁡SS0,S
0 being the value of the ligand concentration at half saturation.
Equation (7) constitutes the second bridge between the Chemistry we aim to describe (via the ligand concentration
S
) and the Physics that we want to use (via the magnetic field
h
).
(iv)
The binding sites can cooperate in a positive manner: this can be modelled by introducing a coupling between theσ variables. The simplest mathematical form is given by a two-body contribution H2 in the cost-function. This term encodes for the reciprocal interactions among binding sites and it reads as (8)H2σ,J=-JN∑i<jNσiσj,
whereJ≥0 is the interaction strength and the sum runs over all possible pairs; the normalization factor 1/N ensures the linear extensivity of the cost-function with respect to the system size. A positive value for J implies an imitative interaction among binding sites: configurations where spins tend to be aligned each others (namely, where sites tend to be either all occupied or all unoccupied) are energetically more favoured and will therefore be more likely.
(v)
Combining together the previous contributions we get the total Hamiltonian:(9)Hσ,J,h=H1σ,h+H2σ,J=-JN∑i<jNσiσj-h∑i=1Nσi=-NJ2NM2+hM.
It is possible to introduce the free energy associated with such a Hamiltonian as(10)Aβ,J,h=1Nlog∑σ1,…,σN2Nexp⁡-βHσ,J,h(11)=1Nlog∑σ1,…,σN2Nexp⁡βJ2N∑i,jNσiσj+βh∑iNσi,
whereβ is the inverse temperature in proper units and the sum runs over all possible configurations. The free energy is a key observable because it corresponds to the difference between the internal energy U and the entropy S (at given temperature), that is, A(β,J,h)=S(β,J,h)-βU(β,J,h). If we could obtain an explicit expression for A(β,J,h) in terms of the order parameter M, we could obtain an expression for the magnetization expected at equilibrium by imposing δMA(β,J,h)=0; in fact, this implies that we are simultaneously asking for the minimum energy and the maximum entropy.
Notice that, having stated the two bridges given by (5) and (7), other mappings between the two fields (e.g., the relation between the coupling strength J and the Hill coefficient nH; see (24) later on) emerge spontaneously as properties of the thermodynamic solutions of the problems.


### 2.2. Resolution of the Problem: The Mechanical Action

We want to find an explicit expression (in terms ofM) for the free energy defined in (10). To this task let us rename -βJ=t and βh=x and let us think of these fictitious variables as a time and a space, respectively. Thus, we can write the free energy as (12)At,x=1Nlog∑σ1,…,σN2Nexp⁡-t2N∑i,jNσiσj+x∑iNσi,where we also wrote ∑i<j‍σiσj~(1/2)∑i,j‍σiσj, which implies vanishing corrections in the thermodynamic limit. If we work out thespatial andtemporal derivatives of the free energy (12) we obtain(13)dAt,xdt=-12M2t,x,dAt,xdx=Mt,x,where the average 〈·〉t,x for a generic observable O depending on the spin configuration is defined as(14)Ot,x=∑σOexp⁡N-t·M2+x·M∑σexp⁡N-t·M2+x·M,and, posing t=-βJ and x=βh, the Boltzmann average for the original system (9) is recovered and this shall be simply denoted as 〈·〉If we now introduce a potentialV(t,x), defined as half the variance of the magnetization, that is, (15)Vt,x=12M2-M2,we see that, by construction, the free energy of this model obeys the following Hamilton-Jacobi equation: (16)dAt,xdt+12dAt,xdx2+Vt,x=0,and therefore A(t,x) is also anaction of Classical Mechanics. We can simplify the previous equation by noticing that, for large enough volumes, the magnetization is a self-averaging quantity [26, 43]; thus in the infinite volume limit the potential must vanish; that is, limN→∞V(t,x)=0. Here, we are restricting to large volumes and we are therefore left with a Hamilton-Jacobi equation describing a free propagation; since the potential is zero, the Lagrangian L coupled to the motion is just the kinetic term: (17)L=12M2,that is, the analogous of the classical formula L=1/2mv2, where the mass m is set unitary (i.e., m=1), and the role of the velocity v is played by the average magnetization M. Solving the Hamilton-Jacobi equation is then straightforward: the solution is formally written as (18)At,x=At=0,x=x0+∫0tLt′,xdt′.The evaluation of the Cauchy condition A(t=0,x=x0) is trivial because, at t=0, the coupling between variables disappears (see (10)), while the integral of the Lagrangian over time reduces to the Lagrangian times time (as the potential is zero). Pasting these two contributions together we obtain (19)At,x=ln⁡2+ln⁡cosh⁡x0+12M2·t.Finally, noticing that the equation of motion is a Galilean trajectory as x(t)=x0+〈M〉t (hence x0=x-〈M〉t) and recasting the solution back in the original variables, that is, t=-βJ and x=βh, we get the free energy associated with this general positive-cooperative reaction: (20)Aβ,J,h=ln⁡2+ln⁡cosh⁡βJM+βh-12βJM2.By extremizing A(β,J,h) with respect to M we get (21)dAβ,J,hdM=0⟹M=tanh⁡βJM+h.This result recovers the well-known self-consistency equation for the order parameter of the Curie-Weiss model in Statistical Mechanics [26, 43].


### 2.3. Chemical Properties of the Physical Solution

The self-consistent equation in (21) is an input-output relation for a general system of binary elements, possibly positively interacting, under the influence of an external field: the input in the system is the external field h and the output is the magnetization M. We can rewrite (21) in a chemical jargon by using the bridges coded in (5) and (7) and fixing, for the sake of simplicity, S0=1; that is, (22)YJ,S=121+tanh⁡J2Y-1+12ln⁡S=Se2J2Y-11+Se2J2Y-1.Before proceeding, we check that if cooperation disappears (i.e., binding sites are reciprocally independent), the Michaelis-Menten scenario is recovered. Posing J=0 in the equation above we get (23)YJ=0,S=S1+S,that is (apart a constant factor that can be reintroduced by taking S0=k, rather than 1), the Michaelis-Menten equation (see (2)).One step forward, we now take into account the couplingJ and relate it to the Hill coefficient nH. The latter is defined in Chemistry as the slope of Y(S) at half saturation (i.e., when Y=1/2), and we can obtain its expression following this prescription by using (22), namely, (24)nH=1Y1-YdYdSY=1/2=11-J.We note that as J→0 we get, as expected, nH→1: if there is no cooperation between binding sites, the Hill coefficient must be unitary; further, the stronger the coupling J, the (hyperbolically) larger the value of the Hill coefficient. In particular, for J→1 the kinetics get ultrasensitive and discontinuities emerge. We remark that, with simpler statistical mechanics model as linear chains of spins, phase transitions are not allowed; hence ultrasensitive behavior can not be taken into account: the present framework is the simplest nontrivial scheme where all these phenomena can be recovered at once (see Figure 1 and [28] for more details on ultrasensitive kinetics).Also, it is worth highlighting the full consistency between our treatment of ultrasensitive kinetics and more standard ones as, for instance, reported in [2] (see eq. 5.17 therein), where the expression for the Hill coefficient can be translated into our formulation as (25)nH=NM2-M21-M2.We see that for 〈M〉→±1 the Hill coefficient diverges, which is the signature of an ultrasensitive behavior: this is perfectly coherent with our approach where, in that limit, the input-output relation (see the hyperbolic tangent (22)) becomes a step function.However, as mentioned in the Introduction, this theory has its flaws, in Chemistry as well as in Mechanics. Regarding the former, the complex picture of yeast’s enzymes evidenced by Koshland [32, 53], where positive and negative cooperativity appear simultaneously (and with the anticooperativity effect getting more and more pronounced as the substrate concentration is raised), still escapes from this mathematical architecture. Further, from the mechanical point of view, two weird things happen: the velocity M is bounded by c=1, while in Classical Mechanics the velocity may diverge; further, if we look at the Boltzmann factor in the free energy (see (12)), this reads as expN-tM2/2+xM and, recalling that the kinetic energy in this mechanical analogy reads as M2/2 (the mass is unitary, thus velocity and momentum coincide), we are allowed to interpret A(β,J,h) as a real action. From this perspective, the exponent can be thought of as the coupling between the stress-energy tensor and the metric tensor: a glance at the form of the Boltzmann factor reveals that the natural underlying metric is (-1,+1) rather than (+1,+1) as in classical Euclidean frames, or in other words, it is of the Minkowskian type. All these details point toward the generalization of the equivalence including special relativity.Plan of the next section is to follow the mechanical path and extend the classical kinetic energy including relativistic corrections and then to investigate its implications. We will see that in the broader, relativistic framework for chemical kinetics the deviations that Koshland explained adding an anticooperative interactions, beyond the cooperative ones, at high ligand’s doses are the chemical analogies of the deviation from classical mechanics at high velocities observed in special relativity.


## 2.1. Formulation of the Problem: The Thermodynamical Free Energy

In the following we summarize the minimal assumptions needed when modelling chemical kinetics from the Statistical Physics perspective; for a more extensive treatment of this kind of modelling we refer to [21, 26, 28, 29, 41], while for a rigorous explanation of the underlying equivalence between Statistical Mechanics and Analytical Mechanics we refer to the seminal works by Guerra [42], dealing with the Sherrington-Kirkpatrick model (and then deepened in, e.g., [43–46]), and by Brankov and Zagrebnov in [47], dealing with the Husimi-Temperley model (and then deepened in, e.g., [48–51]).(i)
Each binding site may or may not be occupied by a ligand: this allows us to code its state (empty versus full) by a Boolean variable. For the genericith site, we will use an Ising spin σi=±1, where σi=-1 represents an empty ith site, and vice versa; σi=+1 means that the ith site is occupied. Clearly, if there are overall N binding sites, i∈(1,…,N).(ii)
It is rather inconvenient (and ultimately unnecessary) to deal with the whole setσi,…,σN if we are interested in the properties of large numbers of these variables (i.e., in the so-calledthermodynamic limit corresponding to N→∞). If we want to distinguish between a fully empty state σi=-1∀i∈(1,…,N) (ordered case), a fully occupied state σi=+1∀i∈(1,…,N) (ordered case), and a completely random case where σi=±1 with equal probability (disordered case), it is convenient to introduce theorder parameter for these variables as themagnetizationM (this term stems from the original application of the Statistical Mechanics model in the context of magnetism) that reads as the arithmetic average of the spin state, namely, (4)M=1N∑i=1Nσi∈-1,+1.There is a univocal relation between the magnetization in Physics and the saturation function in Chemistry, where, we recall, we denote withY∈[0,1] the fractional occupation of the binding sites. In fact, one has [28, 29] (5)Y=121+M.Equation (5) constitutes the first bridge between the Chemistry we aim to describe (via the saturation functionY) and the Physics that we want to use (via the magnetizationM).(iii)
All the binding sites interact with the ligand by the same strength. This is a standard assumption in Chemical Kinetics [29, 31, 52] and it means that the diffusion of the ligands is fast enough to ensure a homogeneous solution. The concentration of free ligands is mapped into a one-body contribution H1 in the cost-function. This term encodes for the action of anexternal magnetic field in such a way that, if the field acting on ith is positive, the spin will tend to align upwards (namely, this direction is energetically favored), and vice versa. This homogenous mixing assumption translates into a homogeneous external field h, and the related contribution reads as (6)H1σ,h=-h∑i=1Nσi.Notice thath plays as a chemical potential and, consistently, it can be related to the substrate concentration as (7)h=12log⁡SS0,S
0 being the value of the ligand concentration at half saturation.
Equation (7) constitutes the second bridge between the Chemistry we aim to describe (via the ligand concentration
S
) and the Physics that we want to use (via the magnetic field
h
).
(iv)
The binding sites can cooperate in a positive manner: this can be modelled by introducing a coupling between theσ variables. The simplest mathematical form is given by a two-body contribution H2 in the cost-function. This term encodes for the reciprocal interactions among binding sites and it reads as (8)H2σ,J=-JN∑i<jNσiσj,
whereJ≥0 is the interaction strength and the sum runs over all possible pairs; the normalization factor 1/N ensures the linear extensivity of the cost-function with respect to the system size. A positive value for J implies an imitative interaction among binding sites: configurations where spins tend to be aligned each others (namely, where sites tend to be either all occupied or all unoccupied) are energetically more favoured and will therefore be more likely.
(v)
Combining together the previous contributions we get the total Hamiltonian:(9)Hσ,J,h=H1σ,h+H2σ,J=-JN∑i<jNσiσj-h∑i=1Nσi=-NJ2NM2+hM.
It is possible to introduce the free energy associated with such a Hamiltonian as(10)Aβ,J,h=1Nlog∑σ1,…,σN2Nexp⁡-βHσ,J,h(11)=1Nlog∑σ1,…,σN2Nexp⁡βJ2N∑i,jNσiσj+βh∑iNσi,
whereβ is the inverse temperature in proper units and the sum runs over all possible configurations. The free energy is a key observable because it corresponds to the difference between the internal energy U and the entropy S (at given temperature), that is, A(β,J,h)=S(β,J,h)-βU(β,J,h). If we could obtain an explicit expression for A(β,J,h) in terms of the order parameter M, we could obtain an expression for the magnetization expected at equilibrium by imposing δMA(β,J,h)=0; in fact, this implies that we are simultaneously asking for the minimum energy and the maximum entropy.
Notice that, having stated the two bridges given by (5) and (7), other mappings between the two fields (e.g., the relation between the coupling strength J and the Hill coefficient nH; see (24) later on) emerge spontaneously as properties of the thermodynamic solutions of the problems.


## 2.2. Resolution of the Problem: The Mechanical Action

We want to find an explicit expression (in terms ofM) for the free energy defined in (10). To this task let us rename -βJ=t and βh=x and let us think of these fictitious variables as a time and a space, respectively. Thus, we can write the free energy as (12)At,x=1Nlog∑σ1,…,σN2Nexp⁡-t2N∑i,jNσiσj+x∑iNσi,where we also wrote ∑i<j‍σiσj~(1/2)∑i,j‍σiσj, which implies vanishing corrections in the thermodynamic limit. If we work out thespatial andtemporal derivatives of the free energy (12) we obtain(13)dAt,xdt=-12M2t,x,dAt,xdx=Mt,x,where the average 〈·〉t,x for a generic observable O depending on the spin configuration is defined as(14)Ot,x=∑σOexp⁡N-t·M2+x·M∑σexp⁡N-t·M2+x·M,and, posing t=-βJ and x=βh, the Boltzmann average for the original system (9) is recovered and this shall be simply denoted as 〈·〉If we now introduce a potentialV(t,x), defined as half the variance of the magnetization, that is, (15)Vt,x=12M2-M2,we see that, by construction, the free energy of this model obeys the following Hamilton-Jacobi equation: (16)dAt,xdt+12dAt,xdx2+Vt,x=0,and therefore A(t,x) is also anaction of Classical Mechanics. We can simplify the previous equation by noticing that, for large enough volumes, the magnetization is a self-averaging quantity [26, 43]; thus in the infinite volume limit the potential must vanish; that is, limN→∞V(t,x)=0. Here, we are restricting to large volumes and we are therefore left with a Hamilton-Jacobi equation describing a free propagation; since the potential is zero, the Lagrangian L coupled to the motion is just the kinetic term: (17)L=12M2,that is, the analogous of the classical formula L=1/2mv2, where the mass m is set unitary (i.e., m=1), and the role of the velocity v is played by the average magnetization M. Solving the Hamilton-Jacobi equation is then straightforward: the solution is formally written as (18)At,x=At=0,x=x0+∫0tLt′,xdt′.The evaluation of the Cauchy condition A(t=0,x=x0) is trivial because, at t=0, the coupling between variables disappears (see (10)), while the integral of the Lagrangian over time reduces to the Lagrangian times time (as the potential is zero). Pasting these two contributions together we obtain (19)At,x=ln⁡2+ln⁡cosh⁡x0+12M2·t.Finally, noticing that the equation of motion is a Galilean trajectory as x(t)=x0+〈M〉t (hence x0=x-〈M〉t) and recasting the solution back in the original variables, that is, t=-βJ and x=βh, we get the free energy associated with this general positive-cooperative reaction: (20)Aβ,J,h=ln⁡2+ln⁡cosh⁡βJM+βh-12βJM2.By extremizing A(β,J,h) with respect to M we get (21)dAβ,J,hdM=0⟹M=tanh⁡βJM+h.This result recovers the well-known self-consistency equation for the order parameter of the Curie-Weiss model in Statistical Mechanics [26, 43].


## 2.3. Chemical Properties of the Physical Solution

The self-consistent equation in (21) is an input-output relation for a general system of binary elements, possibly positively interacting, under the influence of an external field: the input in the system is the external field h and the output is the magnetization M. We can rewrite (21) in a chemical jargon by using the bridges coded in (5) and (7) and fixing, for the sake of simplicity, S0=1; that is, (22)YJ,S=121+tanh⁡J2Y-1+12ln⁡S=Se2J2Y-11+Se2J2Y-1.Before proceeding, we check that if cooperation disappears (i.e., binding sites are reciprocally independent), the Michaelis-Menten scenario is recovered. Posing J=0 in the equation above we get (23)YJ=0,S=S1+S,that is (apart a constant factor that can be reintroduced by taking S0=k, rather than 1), the Michaelis-Menten equation (see (2)).One step forward, we now take into account the couplingJ and relate it to the Hill coefficient nH. The latter is defined in Chemistry as the slope of Y(S) at half saturation (i.e., when Y=1/2), and we can obtain its expression following this prescription by using (22), namely, (24)nH=1Y1-YdYdSY=1/2=11-J.We note that as J→0 we get, as expected, nH→1: if there is no cooperation between binding sites, the Hill coefficient must be unitary; further, the stronger the coupling J, the (hyperbolically) larger the value of the Hill coefficient. In particular, for J→1 the kinetics get ultrasensitive and discontinuities emerge. We remark that, with simpler statistical mechanics model as linear chains of spins, phase transitions are not allowed; hence ultrasensitive behavior can not be taken into account: the present framework is the simplest nontrivial scheme where all these phenomena can be recovered at once (see Figure 1 and [28] for more details on ultrasensitive kinetics).Also, it is worth highlighting the full consistency between our treatment of ultrasensitive kinetics and more standard ones as, for instance, reported in [2] (see eq. 5.17 therein), where the expression for the Hill coefficient can be translated into our formulation as (25)nH=NM2-M21-M2.We see that for 〈M〉→±1 the Hill coefficient diverges, which is the signature of an ultrasensitive behavior: this is perfectly coherent with our approach where, in that limit, the input-output relation (see the hyperbolic tangent (22)) becomes a step function.However, as mentioned in the Introduction, this theory has its flaws, in Chemistry as well as in Mechanics. Regarding the former, the complex picture of yeast’s enzymes evidenced by Koshland [32, 53], where positive and negative cooperativity appear simultaneously (and with the anticooperativity effect getting more and more pronounced as the substrate concentration is raised), still escapes from this mathematical architecture. Further, from the mechanical point of view, two weird things happen: the velocity M is bounded by c=1, while in Classical Mechanics the velocity may diverge; further, if we look at the Boltzmann factor in the free energy (see (12)), this reads as expN-tM2/2+xM and, recalling that the kinetic energy in this mechanical analogy reads as M2/2 (the mass is unitary, thus velocity and momentum coincide), we are allowed to interpret A(β,J,h) as a real action. From this perspective, the exponent can be thought of as the coupling between the stress-energy tensor and the metric tensor: a glance at the form of the Boltzmann factor reveals that the natural underlying metric is (-1,+1) rather than (+1,+1) as in classical Euclidean frames, or in other words, it is of the Minkowskian type. All these details point toward the generalization of the equivalence including special relativity.Plan of the next section is to follow the mechanical path and extend the classical kinetic energy including relativistic corrections and then to investigate its implications. We will see that in the broader, relativistic framework for chemical kinetics the deviations that Koshland explained adding an anticooperative interactions, beyond the cooperative ones, at high ligand’s doses are the chemical analogies of the deviation from classical mechanics at high velocities observed in special relativity.


## 3. TheGeneralized Mathematical Scaffold for Mixed Cooperativity

### 3.1. Relativistic Setting

The relativistic extension of the Hamiltonian (9) is defined by Hamiltonian of the form(26)Hσ,J,hN=-J1+M2-hM,where M=(1/N)∑iN‍σi as usual. Next, we have to insert (26) into the free energy (10):(27)At,x=1Nlog∑σ2Nexp⁡-t·N1+M2+x·M,where we already replaced t=-βJ and x=βh in order to work out their streaming that read as(28)∂At,x∂t=-1+M2t,x,∂At,x∂x=Mt,x,∂tt2At,x-∂xx2At,xN=1-1+M2t,x2+Mt,x2,where the Boltzmann averages ⋯t,x are defined as (using the magnetization as a trial function)(29)Mt,x=∑σMexp⁡N-t·1+M2+x·M∑σexp⁡N-t·1+M2+x·M.As before, the averages ⋯t,x will be denoted by ⋯ whenever evaluated in the sense of thermodynamics (i.e., for t=-βJ and x=βh). By a direct calculation, it is straightforward to see that expression (27) obeys the relativistic Hamilton-Jacobi equation: (30)∂tA2-∂xA2+VNt,x=1,Vt,x=1N□At,x,where the symbol □ represents the D’Alembert operator and V(t,x) is the potential whose expression is chosen in order to make the equation valid by construction and, this time, it is automatically Lorentz invariant. If the functional A(t,x) is sufficiently smooth (i.e., its derivatives are regular functions of t and x), in the thermodynamic limit, we have limN→∞V(t,x)=0; hence in this high-volume limit we are left with(31)∂μA2=1,which is the Klein-Gordon equation for a free relativistic particle with unitary mass in natural units (m0=1).In relativistic mechanics, the stress energy tensor of this particle is defined as(32)pμ=E,γv,where v is the classical velocity of the particle, γ=1/1-v2, and E=γm0=γ is the relativistic energy. In addition, the contravariant momentum is expressed through the action by the following equation:(33)pμ=-∂A∂xμ=1+M2t,x,Mt,x.Comparing (32) and (33), it is immediate to identify the magnetization as the relativistic dynamical variable:(34)Mt,x=v1-v2,while the Lorentz factor is(35)γ=1+Mt,x2.In the thermodynamic limit, the particle is free and its trajectories are the straight lines x=x0+vt. Since the relativistic Lagrangian L=-γ-1 is constant along these classical trajectories, the free energy can be computed as(36)At,x=A0,x0+∫0tdt′Lt′=A0,x0-tγ=log⁡2+log⁡cosh⁡x0-tγ=log⁡2+log⁡cosh⁡x-vt-t1+Mt,x2=log⁡2+log⁡cosh⁡x-Mt,xt1+Mt,x2-t1+Mt,x2.Setting t=-Jβ and x=βh, we finally get an explicit expression for the free energy:(37)Aβ=log⁡2+log⁡cosh⁡JβM1+M2+Jβ1+M2.Requiring that the free energy is extremal with respect to the magnetization (from a thermodynamical perspective this condition can be seen as the simultaneous effect of the minimum energy and the maximum entropy principles and from a mechanical perspective as the minimum action principle), the associated self-consistency equation becomes(38)M=tanh⁡JβM1+M2+βh.


### 3.2. The Classical Limit from a Chemical Perspective

Reading the self-consistency (38) in chemical terms, that is, using the bridges (5) and (7), we obtain(39)YS,J=121+tanh⁡Jβ2Y-11+2Y-12+β2log⁡SS0.We can now check whether, under suitable conditions, this broader theory recovers the classical limit. First, we notice that under the assumption of no interactions among binding sites (i.e., J=0) and replacing h=1/2log(S/S0), the Michaelis-Menten behaviour is recovered. This can be shown by rewriting (39) as(40)YS,J=Se2Jβ2Y-1/1+2Y-121+Se2Jβ2Y-1/1+2Y-12,where we also shifted S/S0→S for simplicity. For J=0 the previous equation reduces to Y(S)=S/(1+S). Further, taking the classical limit, at the lowest order, we have the following expansions:(41)11+M2=1-M22+OM3,M1+M2=M+OM3,such that (38) reduces to (21), in the physical context, and to (22), in the chemical context. Clearly, also the slope at Y=1/2 is preserved; hence, in the classical limit, we recover the expected expression for Hill coefficient (see (24)), namely,(42)nH=1Y1-Y∂Y∂SY=1/2=11-J.


### 3.3. Beyond the Classical Limit

To understand why we expect variations with respect to the Hill paradigm at relatively large values of the substrate concentration, we must check carefully the relativistic self-consistency (38). Let us assume we are working at not too high velocities (i.e., 〈M〉<1) and we can expand the argument inside the hyperbolic tangent; in particular, approximating 1/(1+x2)~1-x2/2, we get (43)M=tanh⁡βJM-βJ2M3+βh.The relativistic effects in chemical kinetics become transparent in this way: if we look at the field felt by the binding sites (i.e., the argument inside the hyperbolic tangent), we see that, beyond the standard Hill term βJ〈M〉 (that positively pairs binding sites together), another term appears that, this time, negatively pairs binding sites together. Retaining this level of approximation, we could write an effective Hamiltonian to generate (43) that reads as (44)Hσ,J,h=-JN∑i<jNσiσj+J2N3∑i<j<k<lσiσjσkσk;hence, beyond the two-body positive coupling coded by the first term, another four-body negative coupling appears. The latter is responsible for the deviation from the classical paradigm and these deviations are in full agreement with the Koshland generalization toward the concept of mixed positive and negative cooperativity [32].In particular, we can see at work the entire reasoning of Koshland who pointed out how, at large enough substrate concentration, the positivity of the reaction diminishes. In fact, for〈M〉~0 no relativistic effect can be noted. By increasing S (the input in the system), we get a growth in 〈M〉 (the output in the system): the latter raises in response of S and it is enhanced because of the two-body term in the effective Hamiltonian (44), the four-body term still being negligible. As S keeps on growing, 〈M〉 increases as well, up to a point where it reaches high enough values such that, from now on, also the four-body term inside the effective Hamiltonian (44) becomes relevant. At this point, a novel, anticooperative effect is naturally induced in the reaction and it yields to a reduction of the Hill coefficient. In the next analysis these qualitative remarks shall be addressed in more details.We focus on the definition of the Hill coefficient based on the Hill equation:(45)Y=SnHk+SnH.This equation accounts for the possibility that not all receptor sites are independent: here nH is the average number of interacting sites and the slope of the Hill plot. The latter is based on a linear transformation made by rearranging (45) and taking the logarithm: (46)logY1-Y=nlog⁡S-log⁡k.Thus, one plots log⁡Y/(1-Y) versus log⁡S, fits with a linear function and the resulting slope, calculated at the half-saturation point, and provides the Hill coefficient. As already underlined, the Michaelis-Menten theory corresponds to nH=1 and any deviations from a slope of 1 tell us about deviation from that model.For the (classical and relativistic) models analyzed here (coded in the Hamiltonians (9) and (26)) we can estimate the slope nH directly from the self-consistency equations (22) and (39). Let us start with the classical model. We preliminary notice that (47)ddlog⁡Slog⁡Y1-YY=1/2=1Y1-YdYdlog⁡SY=1/2=4dYdlog⁡SY=1/2.Therefore, we just need to evaluate dY/dlog(S) in Y=1/2, which reads as (48)dYdlog⁡SY=1/2=12sech2J2Y-1+12log⁡S2JdYdlog⁡S+12Y=1/2.Posing x=dY/dlog⁡SY=1/2 and noticing that S=1 when Y=1/2, we have (49)x=122Jx+12⟹x=1411-J.By plugging this result in (47), we finally have (50)nHclass=11-J.One can see that when J=0 the Hill coefficient is unitary as expected for noncooperative systems, when J>0 the coefficient is larger than 1, indicating that receptors are interacting, and when J<0 the coefficient is smaller than 1, as expected for negative cooperativity.Let us now move to the relativistic model. Again, we just need to evaluatedY/dlog(S) in Y=1/2, which, recalling (39), reads as(51)dYdlog⁡SY=1/2=12sech2J2Y-11+2Y-12+12log⁡S(52)×2J2+4YY-13/2dYdlog⁡S+12Y=1/2.Exploiting the fact that S=1 when Y=1/2, the previous expression simplifies as (53)dYdlog⁡SY=1/21411-J/27.Thus, we can write (54)nHrel=ddlog⁡Slog⁡Y1-YY=1/2=4dYdlog⁡SY=1/2=11-J/27.Note that nHclass/nHrel<1, confirming that the relativistic correction weakens the emerging cooperativity.


### 3.4. Further Robustness Checks

As stressed above, for a fixed interaction couplingJ, the relativistic model is expected to exhibit a lower cooperativity with respect to the classical model. In order to quantify this point we considered different quantifiers for cooperativity and we compared the outcomes for the relativistic and the classical models set at the same value of J. Let us start with the Koshland measure of cooperativity which is defined as the ratio (notice that the Koshland index κ is actually strongly related to the Hill coefficient (see, e.g., [2])) (55)κ=R0.9R0.1,where R0.9 denotes the substrate concentration corresponding to a 90% saturation, while R0.1 denotes the substrate concentration corresponding to a 10% saturation; that is, Y(R0.9)=0.9 and Y(R0.1)=0.1. In the noncooperative case one has R0.9/R0.1=81 and, accordingly, if the ratio is smaller than 81 (meaning that the saturation curve is relatively steep) one has positive cooperativity, while if the ratio is larger than 81 one has negative cooperativity. The advantage in using the index κ is that it can be easily measured, yet it ignores all information that can be derived from the shape of Y(S). In particular, this quantifier can be estimated starting from a Klotz plot (see, e.g., Figure 2(a)) where the saturation function is shown versus the logarithm of the (free) ligand concentration; in the presence of positive cooperativity this plot yields to a characteristic sigmoidal curve. For the models analyzed here we can estimate R0.9/R0.1 directly from the self-consistency equations (22)–(24), (26)–(39). Starting from the classical model and posing Y=0.9 and Y=0.1 we get, respectively,(56)910=121+tanh⁡J2×910-1+12log⁡S0.9,110=121+tanh⁡J2×110-1+12log⁡S0.1,and, with some algebra (recalling 2atanh⁡x=log⁡[(1+x)/(1-x)]),(57)log⁡S0.9=2atanh⁡45-85J=log⁡9-85J⟹S0.9=9e-8J/5,log⁡S0.1=2atanh⁡-45+85J=-log⁡9+85J⟹S0.1=19e8J/5;that is, (58)κclass=S0.9S0.1=81e-16J/5.Of course, when J=0 we recover the value 81, when J>0 we get Rclass<81, and when J<0 we get Rclass>81.Figure 2
Summary of analysis on cooperativity. (a) Klotz plot showing the saturation function Y versus the substrate concentration S (notice the logarithmic scale on the x-axis) for the relativistic (dashed line) and the classical (solid line) models. Data for Y are collected by solving numerically the self-consistent equations ((39) and (22), resp.) for J=0.9 and different values of S. Both models exhibit the sigmoidal shape typical of cooperative systems; however, the former displays a slower saturation. Analogous results are obtained for different values of J>0. (b) A Scatchard plot is built with the same collection of data by showing the ratio Y/S versus Y. Both models exhibit the concave-down shape typical of cooperative systems. However, for relatively small values of J the plot for the relativistic model is monotonically decreasing (see also Figure 3). (c) A Hill plot is built with the same collection of data by showing θ=Y/(1-Y) versus S; both observables are taken under the logarithm. When S is close to one (here S0=1) the relativistic and the classical model give overlapped curves, while when S is either very large or very small the two curves are shifted. (d) By further analyzing the plots in the previous panels we can derive estimates for the extent of cooperativity characterizing the systems. As explained in the main text, starting from data in (a) we measured the Kloshand quantifier κ=S0.9/S0.10 (⋄), by extrapolating the maximum value for data in (b) we get σ (∇), and by fitting the data in (c) at the half-saturation point we get nH (□). These estimates are obtained for both the relativistic (white symbol) and the classic (black symbols) models.
(a)
(b)
(c)
(d)Figure 3
The role of the interaction coupling
J. We resume the plots presented in (a)–(c) of Figure 2 and we show how they are affected by the interaction coupling J. In particular, we compare the outputs for J=0.1 (black), J=0.5 (blue), and J=0.9 (bright blue). Again, the relativistic model (dashed line) and the classical model (solid line) are compared. Notice that the gap between relativistic and classical model is larger when J is relatively large.
(a)
(b)
(c)Repeating analogous calculations for the relativistic model we get(59)910=121+tanhJ2×9/10-11+2×9/10-12+12log⁡S0.9,110=121+tanhJ2×1/10-11+2×1/10-12+12log⁡S0.1,and, with some algebra,(60)log⁡S0.9=2atanh⁡45-841J=log⁡9-841J⟹S0.9=9e-8J/41,log⁡S0.1=2atanh⁡-45+841J=-log⁡9+841J⟹S0.1=19e8J/41;that is, (61)κrel=S0.9S0.1=81e-16J/41.Again, one can check that when J=0 we recover the value 81, when J>0 we get κrel<81, and when J<0 we get κrel>81. Also, κrel/κclass=e-16J/41+16J/5>1. This means that, even with this quantifier, when fixing the same coupling constant J, the emerging cooperativity is weaker for the relativistic model, as expected.Next, let us consider the cooperativity quantifier derived from the Scatchard plot. We recall that this plot is built by showing the behavior ofY/S with respect to Y. In fact, according to the simplest scenario (this corresponds to the Michaelis-Menten theory and to Clark’s theory and it requires a set of simplifying assumptions, among which the interaction is reversible; all the binding molecules are equivalent and independent; the biological response is proportional to the number of occupied binding sites; the substrate only exists in either a free (i.e., unbound) form or bound to the receptor), at equilibrium, one can write (62)S1-YY=k,where k is the proportionality constant between response and occupancy (i.e., it is the ratio between the dissociation and the association constants), and rearranging (62) we have (63)YS=-Yk+1k.The previous expression fits the equation of a line for Y/S versus Y, whose slope is -1/k. The advantages in using the Scatchard plot is that it is a very powerful tool for identifying deviations from the simple model, which, without deviations, is represented by a straight line. In particular, a concave-up curve may indicate the presence of negative cooperativity between binding sites, while a concave-down curve is indicative of positive cooperativity. Also, in the latter case, the maxima occurs at the fractional occupancy Y∗ which fulfills (64)Y∗=σ-1σ,where σ provides another quantifier for cooperativity.Starting from the classical model, we can build the functionY/S, by first getting S as a function of Y, and can be obtained by inverting formula (22); namely, (65)SY=exp⁡2atanh⁡2Y-1-2J2Y-1.By deriving Y/S with respect to Y we get (66)ddYYSY=-e2J2Y-11-4J1-Y,which is identically equal to -1 when J=0, monotonically decreasing with Y when J>0 and monotonically increasing with Y when J<0. The (possible) root therefore provides the extremal point; that is, (67)Y∗=4J-14J,and, comparing with (64), we get (68)σclass=4J.We now repeat analogous calculations for the relativistic model. First, we getS as a function of Y, by inverting formula (39), namely, (69)S=exp2atanh⁡2Y-1-2J2Y-11+2Y-12.By deriving Y/S with respect to Y we get (70)ddYYS=-eJ2Y-1/YY-1+1/21-4J1-Y2-4Y1-Y3/2,which is again identically equal to -1 when J=0, but it is no longer monotonic when J≠0. More precisely, by studying d(Y/S)/dY we can derive that when J is relatively small, Y/S does not exhibit any extremal points, but there is a flex at intermediate values of Y; for intermediate values of J there is a minimum at small values of Y and a maximum at larger values of Y; for large values of J there is a maximum. The extremal points can be found as roots of a 6th degree function of Y. We can obtain an estimate of the value Y∗ corresponding to the maximum by recalling Y≤1 and neglecting high-order terms. In this way we get (71)Y∗≈-3+2J2--9+26J22-9+J2,and, comparing with (64), we get (72)σrel=2-9+J2-15+-9+26J2.The three plots considered here (i.e., Klotz, Scatchard, and Hill) and the related estimates for the extent of cooperativity are presented in Figure2. In particular, in (d) we compare the cooperativity quantifiers for several values of J: as anticipated, in general, for a given value of J, the relativistic model gives rise to a weaker cooperativity.We proceed our analysis by deepening the role of the coupling constantJ in the binding curves related to the two models. In Figure 3 we present Klotz’s plot (a), the Scatchard plot (b), and the Hill plot (c) for the relativistic and the classic models, comparing the outcomes for different values of J. As expected, the point corresponding to S=1 and Y=1/2 is a fixed point in each plot and, in general, the gap between the two models is enhanced when J is larger (i.e., when J is closer to 1). Also, when J is not too small, the Scatchard plot for the relativistic model displays a flex at small values of Y suggesting that when the saturation is small, the cooperativity is not truly positive.In the final part of this section we want to get deeper in the comparison between the classical and the relativistic models. To this aim, we solved numerically (39), for different values of S and of J, getting a set of data Y(S,J). We can think of this set of data as the result of a set of measurements where we collect the saturation value at a given substrate concentration. Now, assuming that in this experiment we have no hints about the underlying cooperative mechanisms, we may apply the formulas for the plain positive cooperativity and infer the value of J. More practically, we calculate numerically Y from the relativistic model for different values of S and of the coupling strength, referred to as Jrel for clarity. Next, we manipulate the set of data Y(S,Jrel) by inverting the formula in (22): as the value of S is assumed to be known, we can derive the coupling strength, referred to as Jclass, expected within a classical framework. In this way, we can compare the original coupling constant Jrel with the inferred one Jclass. We can translate these procedures in formulas as follows:(73)Jclass=atanh⁡2Y-1-1/2log⁡S2Y-12Y-1=tanh⁡Jrel2Y-11+2Y-12+12log⁡S⇓Jclass=Jrel1+2Y-12≤Jrel,with equality holding only when Y=1/2.In Figure4(a) we plot Jclass versus Jrel, for different values of S. Notice that the two parameters are related by a linear law, whose slope is smaller than 1 and decreases with S. This confirms that the relativistic model yields to a weak cooperativity. The negative contributions in the relativistic model get more effective when Jrel and S are large, as further highlighted in Figure 4(b).Figure 4
Comparison between relativistic and classical model. We performed numerical experiments where we obtained Y(S,Jrel) for the relativistic model by solving the self-consistent equation (39). From this data we inferred the expected classical coupling Jclass by inverting the self-consistent equation (22). We repeated the same operations for several values of S and Jrel. In (a) we show the inferred Jclass versus the fixed Jrel: different colors represent different values of S and the identity function is also shown for reference (dashed, black curve). Notice that, in general Jclass<Jrel and the inequality is enhanced as S grows. In (b) we show a contour plot for the ratio Jclass/Jrel versus h=log(S)/2 and Jrel. Again, one can notice that, in general, Jclass/Jrel<1 and this inequality is enhanced for relatively large values of S.
(a)
(b)


## 3.1. Relativistic Setting

The relativistic extension of the Hamiltonian (9) is defined by Hamiltonian of the form(26)Hσ,J,hN=-J1+M2-hM,where M=(1/N)∑iN‍σi as usual. Next, we have to insert (26) into the free energy (10):(27)At,x=1Nlog∑σ2Nexp⁡-t·N1+M2+x·M,where we already replaced t=-βJ and x=βh in order to work out their streaming that read as(28)∂At,x∂t=-1+M2t,x,∂At,x∂x=Mt,x,∂tt2At,x-∂xx2At,xN=1-1+M2t,x2+Mt,x2,where the Boltzmann averages ⋯t,x are defined as (using the magnetization as a trial function)(29)Mt,x=∑σMexp⁡N-t·1+M2+x·M∑σexp⁡N-t·1+M2+x·M.As before, the averages ⋯t,x will be denoted by ⋯ whenever evaluated in the sense of thermodynamics (i.e., for t=-βJ and x=βh). By a direct calculation, it is straightforward to see that expression (27) obeys the relativistic Hamilton-Jacobi equation: (30)∂tA2-∂xA2+VNt,x=1,Vt,x=1N□At,x,where the symbol □ represents the D’Alembert operator and V(t,x) is the potential whose expression is chosen in order to make the equation valid by construction and, this time, it is automatically Lorentz invariant. If the functional A(t,x) is sufficiently smooth (i.e., its derivatives are regular functions of t and x), in the thermodynamic limit, we have limN→∞V(t,x)=0; hence in this high-volume limit we are left with(31)∂μA2=1,which is the Klein-Gordon equation for a free relativistic particle with unitary mass in natural units (m0=1).In relativistic mechanics, the stress energy tensor of this particle is defined as(32)pμ=E,γv,where v is the classical velocity of the particle, γ=1/1-v2, and E=γm0=γ is the relativistic energy. In addition, the contravariant momentum is expressed through the action by the following equation:(33)pμ=-∂A∂xμ=1+M2t,x,Mt,x.Comparing (32) and (33), it is immediate to identify the magnetization as the relativistic dynamical variable:(34)Mt,x=v1-v2,while the Lorentz factor is(35)γ=1+Mt,x2.In the thermodynamic limit, the particle is free and its trajectories are the straight lines x=x0+vt. Since the relativistic Lagrangian L=-γ-1 is constant along these classical trajectories, the free energy can be computed as(36)At,x=A0,x0+∫0tdt′Lt′=A0,x0-tγ=log⁡2+log⁡cosh⁡x0-tγ=log⁡2+log⁡cosh⁡x-vt-t1+Mt,x2=log⁡2+log⁡cosh⁡x-Mt,xt1+Mt,x2-t1+Mt,x2.Setting t=-Jβ and x=βh, we finally get an explicit expression for the free energy:(37)Aβ=log⁡2+log⁡cosh⁡JβM1+M2+Jβ1+M2.Requiring that the free energy is extremal with respect to the magnetization (from a thermodynamical perspective this condition can be seen as the simultaneous effect of the minimum energy and the maximum entropy principles and from a mechanical perspective as the minimum action principle), the associated self-consistency equation becomes(38)M=tanh⁡JβM1+M2+βh.


## 3.2. The Classical Limit from a Chemical Perspective

Reading the self-consistency (38) in chemical terms, that is, using the bridges (5) and (7), we obtain(39)YS,J=121+tanh⁡Jβ2Y-11+2Y-12+β2log⁡SS0.We can now check whether, under suitable conditions, this broader theory recovers the classical limit. First, we notice that under the assumption of no interactions among binding sites (i.e., J=0) and replacing h=1/2log(S/S0), the Michaelis-Menten behaviour is recovered. This can be shown by rewriting (39) as(40)YS,J=Se2Jβ2Y-1/1+2Y-121+Se2Jβ2Y-1/1+2Y-12,where we also shifted S/S0→S for simplicity. For J=0 the previous equation reduces to Y(S)=S/(1+S). Further, taking the classical limit, at the lowest order, we have the following expansions:(41)11+M2=1-M22+OM3,M1+M2=M+OM3,such that (38) reduces to (21), in the physical context, and to (22), in the chemical context. Clearly, also the slope at Y=1/2 is preserved; hence, in the classical limit, we recover the expected expression for Hill coefficient (see (24)), namely,(42)nH=1Y1-Y∂Y∂SY=1/2=11-J.


## 3.3. Beyond the Classical Limit

To understand why we expect variations with respect to the Hill paradigm at relatively large values of the substrate concentration, we must check carefully the relativistic self-consistency (38). Let us assume we are working at not too high velocities (i.e., 〈M〉<1) and we can expand the argument inside the hyperbolic tangent; in particular, approximating 1/(1+x2)~1-x2/2, we get (43)M=tanh⁡βJM-βJ2M3+βh.The relativistic effects in chemical kinetics become transparent in this way: if we look at the field felt by the binding sites (i.e., the argument inside the hyperbolic tangent), we see that, beyond the standard Hill term βJ〈M〉 (that positively pairs binding sites together), another term appears that, this time, negatively pairs binding sites together. Retaining this level of approximation, we could write an effective Hamiltonian to generate (43) that reads as (44)Hσ,J,h=-JN∑i<jNσiσj+J2N3∑i<j<k<lσiσjσkσk;hence, beyond the two-body positive coupling coded by the first term, another four-body negative coupling appears. The latter is responsible for the deviation from the classical paradigm and these deviations are in full agreement with the Koshland generalization toward the concept of mixed positive and negative cooperativity [32].In particular, we can see at work the entire reasoning of Koshland who pointed out how, at large enough substrate concentration, the positivity of the reaction diminishes. In fact, for〈M〉~0 no relativistic effect can be noted. By increasing S (the input in the system), we get a growth in 〈M〉 (the output in the system): the latter raises in response of S and it is enhanced because of the two-body term in the effective Hamiltonian (44), the four-body term still being negligible. As S keeps on growing, 〈M〉 increases as well, up to a point where it reaches high enough values such that, from now on, also the four-body term inside the effective Hamiltonian (44) becomes relevant. At this point, a novel, anticooperative effect is naturally induced in the reaction and it yields to a reduction of the Hill coefficient. In the next analysis these qualitative remarks shall be addressed in more details.We focus on the definition of the Hill coefficient based on the Hill equation:(45)Y=SnHk+SnH.This equation accounts for the possibility that not all receptor sites are independent: here nH is the average number of interacting sites and the slope of the Hill plot. The latter is based on a linear transformation made by rearranging (45) and taking the logarithm: (46)logY1-Y=nlog⁡S-log⁡k.Thus, one plots log⁡Y/(1-Y) versus log⁡S, fits with a linear function and the resulting slope, calculated at the half-saturation point, and provides the Hill coefficient. As already underlined, the Michaelis-Menten theory corresponds to nH=1 and any deviations from a slope of 1 tell us about deviation from that model.For the (classical and relativistic) models analyzed here (coded in the Hamiltonians (9) and (26)) we can estimate the slope nH directly from the self-consistency equations (22) and (39). Let us start with the classical model. We preliminary notice that (47)ddlog⁡Slog⁡Y1-YY=1/2=1Y1-YdYdlog⁡SY=1/2=4dYdlog⁡SY=1/2.Therefore, we just need to evaluate dY/dlog(S) in Y=1/2, which reads as (48)dYdlog⁡SY=1/2=12sech2J2Y-1+12log⁡S2JdYdlog⁡S+12Y=1/2.Posing x=dY/dlog⁡SY=1/2 and noticing that S=1 when Y=1/2, we have (49)x=122Jx+12⟹x=1411-J.By plugging this result in (47), we finally have (50)nHclass=11-J.One can see that when J=0 the Hill coefficient is unitary as expected for noncooperative systems, when J>0 the coefficient is larger than 1, indicating that receptors are interacting, and when J<0 the coefficient is smaller than 1, as expected for negative cooperativity.Let us now move to the relativistic model. Again, we just need to evaluatedY/dlog(S) in Y=1/2, which, recalling (39), reads as(51)dYdlog⁡SY=1/2=12sech2J2Y-11+2Y-12+12log⁡S(52)×2J2+4YY-13/2dYdlog⁡S+12Y=1/2.Exploiting the fact that S=1 when Y=1/2, the previous expression simplifies as (53)dYdlog⁡SY=1/21411-J/27.Thus, we can write (54)nHrel=ddlog⁡Slog⁡Y1-YY=1/2=4dYdlog⁡SY=1/2=11-J/27.Note that nHclass/nHrel<1, confirming that the relativistic correction weakens the emerging cooperativity.


## 3.4. Further Robustness Checks

As stressed above, for a fixed interaction couplingJ, the relativistic model is expected to exhibit a lower cooperativity with respect to the classical model. In order to quantify this point we considered different quantifiers for cooperativity and we compared the outcomes for the relativistic and the classical models set at the same value of J. Let us start with the Koshland measure of cooperativity which is defined as the ratio (notice that the Koshland index κ is actually strongly related to the Hill coefficient (see, e.g., [2])) (55)κ=R0.9R0.1,where R0.9 denotes the substrate concentration corresponding to a 90% saturation, while R0.1 denotes the substrate concentration corresponding to a 10% saturation; that is, Y(R0.9)=0.9 and Y(R0.1)=0.1. In the noncooperative case one has R0.9/R0.1=81 and, accordingly, if the ratio is smaller than 81 (meaning that the saturation curve is relatively steep) one has positive cooperativity, while if the ratio is larger than 81 one has negative cooperativity. The advantage in using the index κ is that it can be easily measured, yet it ignores all information that can be derived from the shape of Y(S). In particular, this quantifier can be estimated starting from a Klotz plot (see, e.g., Figure 2(a)) where the saturation function is shown versus the logarithm of the (free) ligand concentration; in the presence of positive cooperativity this plot yields to a characteristic sigmoidal curve. For the models analyzed here we can estimate R0.9/R0.1 directly from the self-consistency equations (22)–(24), (26)–(39). Starting from the classical model and posing Y=0.9 and Y=0.1 we get, respectively,(56)910=121+tanh⁡J2×910-1+12log⁡S0.9,110=121+tanh⁡J2×110-1+12log⁡S0.1,and, with some algebra (recalling 2atanh⁡x=log⁡[(1+x)/(1-x)]),(57)log⁡S0.9=2atanh⁡45-85J=log⁡9-85J⟹S0.9=9e-8J/5,log⁡S0.1=2atanh⁡-45+85J=-log⁡9+85J⟹S0.1=19e8J/5;that is, (58)κclass=S0.9S0.1=81e-16J/5.Of course, when J=0 we recover the value 81, when J>0 we get Rclass<81, and when J<0 we get Rclass>81.Figure 2
Summary of analysis on cooperativity. (a) Klotz plot showing the saturation function Y versus the substrate concentration S (notice the logarithmic scale on the x-axis) for the relativistic (dashed line) and the classical (solid line) models. Data for Y are collected by solving numerically the self-consistent equations ((39) and (22), resp.) for J=0.9 and different values of S. Both models exhibit the sigmoidal shape typical of cooperative systems; however, the former displays a slower saturation. Analogous results are obtained for different values of J>0. (b) A Scatchard plot is built with the same collection of data by showing the ratio Y/S versus Y. Both models exhibit the concave-down shape typical of cooperative systems. However, for relatively small values of J the plot for the relativistic model is monotonically decreasing (see also Figure 3). (c) A Hill plot is built with the same collection of data by showing θ=Y/(1-Y) versus S; both observables are taken under the logarithm. When S is close to one (here S0=1) the relativistic and the classical model give overlapped curves, while when S is either very large or very small the two curves are shifted. (d) By further analyzing the plots in the previous panels we can derive estimates for the extent of cooperativity characterizing the systems. As explained in the main text, starting from data in (a) we measured the Kloshand quantifier κ=S0.9/S0.10 (⋄), by extrapolating the maximum value for data in (b) we get σ (∇), and by fitting the data in (c) at the half-saturation point we get nH (□). These estimates are obtained for both the relativistic (white symbol) and the classic (black symbols) models.
(a)
(b)
(c)
(d)Figure 3
The role of the interaction coupling
J. We resume the plots presented in (a)–(c) of Figure 2 and we show how they are affected by the interaction coupling J. In particular, we compare the outputs for J=0.1 (black), J=0.5 (blue), and J=0.9 (bright blue). Again, the relativistic model (dashed line) and the classical model (solid line) are compared. Notice that the gap between relativistic and classical model is larger when J is relatively large.
(a)
(b)
(c)Repeating analogous calculations for the relativistic model we get(59)910=121+tanhJ2×9/10-11+2×9/10-12+12log⁡S0.9,110=121+tanhJ2×1/10-11+2×1/10-12+12log⁡S0.1,and, with some algebra,(60)log⁡S0.9=2atanh⁡45-841J=log⁡9-841J⟹S0.9=9e-8J/41,log⁡S0.1=2atanh⁡-45+841J=-log⁡9+841J⟹S0.1=19e8J/41;that is, (61)κrel=S0.9S0.1=81e-16J/41.Again, one can check that when J=0 we recover the value 81, when J>0 we get κrel<81, and when J<0 we get κrel>81. Also, κrel/κclass=e-16J/41+16J/5>1. This means that, even with this quantifier, when fixing the same coupling constant J, the emerging cooperativity is weaker for the relativistic model, as expected.Next, let us consider the cooperativity quantifier derived from the Scatchard plot. We recall that this plot is built by showing the behavior ofY/S with respect to Y. In fact, according to the simplest scenario (this corresponds to the Michaelis-Menten theory and to Clark’s theory and it requires a set of simplifying assumptions, among which the interaction is reversible; all the binding molecules are equivalent and independent; the biological response is proportional to the number of occupied binding sites; the substrate only exists in either a free (i.e., unbound) form or bound to the receptor), at equilibrium, one can write (62)S1-YY=k,where k is the proportionality constant between response and occupancy (i.e., it is the ratio between the dissociation and the association constants), and rearranging (62) we have (63)YS=-Yk+1k.The previous expression fits the equation of a line for Y/S versus Y, whose slope is -1/k. The advantages in using the Scatchard plot is that it is a very powerful tool for identifying deviations from the simple model, which, without deviations, is represented by a straight line. In particular, a concave-up curve may indicate the presence of negative cooperativity between binding sites, while a concave-down curve is indicative of positive cooperativity. Also, in the latter case, the maxima occurs at the fractional occupancy Y∗ which fulfills (64)Y∗=σ-1σ,where σ provides another quantifier for cooperativity.Starting from the classical model, we can build the functionY/S, by first getting S as a function of Y, and can be obtained by inverting formula (22); namely, (65)SY=exp⁡2atanh⁡2Y-1-2J2Y-1.By deriving Y/S with respect to Y we get (66)ddYYSY=-e2J2Y-11-4J1-Y,which is identically equal to -1 when J=0, monotonically decreasing with Y when J>0 and monotonically increasing with Y when J<0. The (possible) root therefore provides the extremal point; that is, (67)Y∗=4J-14J,and, comparing with (64), we get (68)σclass=4J.We now repeat analogous calculations for the relativistic model. First, we getS as a function of Y, by inverting formula (39), namely, (69)S=exp2atanh⁡2Y-1-2J2Y-11+2Y-12.By deriving Y/S with respect to Y we get (70)ddYYS=-eJ2Y-1/YY-1+1/21-4J1-Y2-4Y1-Y3/2,which is again identically equal to -1 when J=0, but it is no longer monotonic when J≠0. More precisely, by studying d(Y/S)/dY we can derive that when J is relatively small, Y/S does not exhibit any extremal points, but there is a flex at intermediate values of Y; for intermediate values of J there is a minimum at small values of Y and a maximum at larger values of Y; for large values of J there is a maximum. The extremal points can be found as roots of a 6th degree function of Y. We can obtain an estimate of the value Y∗ corresponding to the maximum by recalling Y≤1 and neglecting high-order terms. In this way we get (71)Y∗≈-3+2J2--9+26J22-9+J2,and, comparing with (64), we get (72)σrel=2-9+J2-15+-9+26J2.The three plots considered here (i.e., Klotz, Scatchard, and Hill) and the related estimates for the extent of cooperativity are presented in Figure2. In particular, in (d) we compare the cooperativity quantifiers for several values of J: as anticipated, in general, for a given value of J, the relativistic model gives rise to a weaker cooperativity.We proceed our analysis by deepening the role of the coupling constantJ in the binding curves related to the two models. In Figure 3 we present Klotz’s plot (a), the Scatchard plot (b), and the Hill plot (c) for the relativistic and the classic models, comparing the outcomes for different values of J. As expected, the point corresponding to S=1 and Y=1/2 is a fixed point in each plot and, in general, the gap between the two models is enhanced when J is larger (i.e., when J is closer to 1). Also, when J is not too small, the Scatchard plot for the relativistic model displays a flex at small values of Y suggesting that when the saturation is small, the cooperativity is not truly positive.In the final part of this section we want to get deeper in the comparison between the classical and the relativistic models. To this aim, we solved numerically (39), for different values of S and of J, getting a set of data Y(S,J). We can think of this set of data as the result of a set of measurements where we collect the saturation value at a given substrate concentration. Now, assuming that in this experiment we have no hints about the underlying cooperative mechanisms, we may apply the formulas for the plain positive cooperativity and infer the value of J. More practically, we calculate numerically Y from the relativistic model for different values of S and of the coupling strength, referred to as Jrel for clarity. Next, we manipulate the set of data Y(S,Jrel) by inverting the formula in (22): as the value of S is assumed to be known, we can derive the coupling strength, referred to as Jclass, expected within a classical framework. In this way, we can compare the original coupling constant Jrel with the inferred one Jclass. We can translate these procedures in formulas as follows:(73)Jclass=atanh⁡2Y-1-1/2log⁡S2Y-12Y-1=tanh⁡Jrel2Y-11+2Y-12+12log⁡S⇓Jclass=Jrel1+2Y-12≤Jrel,with equality holding only when Y=1/2.In Figure4(a) we plot Jclass versus Jrel, for different values of S. Notice that the two parameters are related by a linear law, whose slope is smaller than 1 and decreases with S. This confirms that the relativistic model yields to a weak cooperativity. The negative contributions in the relativistic model get more effective when Jrel and S are large, as further highlighted in Figure 4(b).Figure 4
Comparison between relativistic and classical model. We performed numerical experiments where we obtained Y(S,Jrel) for the relativistic model by solving the self-consistent equation (39). From this data we inferred the expected classical coupling Jclass by inverting the self-consistent equation (22). We repeated the same operations for several values of S and Jrel. In (a) we show the inferred Jclass versus the fixed Jrel: different colors represent different values of S and the identity function is also shown for reference (dashed, black curve). Notice that, in general Jclass<Jrel and the inequality is enhanced as S grows. In (b) we show a contour plot for the ratio Jclass/Jrel versus h=log(S)/2 and Jrel. Again, one can notice that, in general, Jclass/Jrel<1 and this inequality is enhanced for relatively large values of S.
(a)
(b)


## 4. Conclusions

The rewards in the overall bridge linking Chemical Kinetics and Analytical Mechanics are several, both theoretical and practical, as we briefly comment.The former lie in a deeper understanding of the mathematical scaffold for modelling real phenomena: it is far from trivial that the description of chemical/thermodynamical equilibrium is formally the same as the mechanical one. In particular, the self-consistency relation (38) that emerges from the thermodynamic principles (in fact, it stems from the requirement of simultaneous entropy maximization and energy minimization) also turns out to be, in the mechanical analogy, a direct consequence of the least action principles δA(t,x)=0. This means that the stationary point corresponds to a light perturbation of the evolution of the system in the interval [0,t]. Explicitly, we shift infinitesimally Mt,x→Mt,x+δMt,x; then(74)0=δAt,x=∂At,x∂Mt,xδMt,x=tanh⁡x-Mt,xt1+Mt,x2-tδMt,x1+Mt,x23/2+Mt,xδMt,xt1+Mt,x23/2=0,from which (38) is recovered (as usual by setting t=-Jβ and x=βh), since this holds for all variations δMt,x.Even more exciting, still by the theoretical side, is the realization of thecomplexity of systems presenting mixed reaction (i.e., where both positive and negative cooperativity are simultaneously at work) and the possible applications in information processing, as we are going to discuss.First, let us clarify that in the Literature we speak ofcomplex network orcomplex system with (mainly) two, rather distinct, meanings: in full generality, let us consider a Hamiltonian as (75)Hσ,J=∑i<jJijσiσjand let us write the two-body coupling matrix as Jij=AijWij, where A is the adjacency matrix, accounting for the bare topology of the system (its entry Aij is 1 if there is a link connecting the related nodes (i,j), which are therefore allowed to interact each other, and it is zero otherwise) and W is the weight matrix, accounting for the sign and the magnitude of the links (i.e., the type of interactions among binding sites).Dealing withA, networks where the topology is very heterogeneous (e.g., the distribution of the number of links stemming from a node has a power-law scaling) are calledcomplex networks, as it is case for the Barabasi-Albert model [54].Dealing withW, networks where the entries of the weight matrix are both positive and negative are termedcomplex systems, as the Sherrington-Kirkpatrick model [55] for the so-called spin glasses.Crucially, spin glasses spontaneously show very general information-processing skills and computational capabilities: for instance, Hopfield neural networks [56] and restricted Boltzmann machines [57], key tools in Artificial Intelligence (resp., inneural networks andmachine learning), are two types of spin glasses and it is with this last definition of complexity that we now can read the information processing capabilities of the elementary reactions we studied. For a given macromolecule under consideration, we could paste each binding site on a node and draw the links among nodes that are interacting: if two nodes are correlated (they show positive cooperativity), their relative interaction is positive, while if two nodes are anticorrelated (they show negative cooperativity), their relative interaction is negative. Dealing with mixed reactions we have to deal with spin glasses and we can thus assess how much information has been processed in a given reaction by evaluating the amount of information processed in its corresponding spin-glass representation, using our bridge. We have already started this investigation in [21, 28, 41].Finally, from a practical perspective, in the classical limit (i.e., forsimple reactions) we have an explicit expression that directly relates the Hill coefficient nH, which can be measured experimentally, to the interaction coupling J assumed in the model; namely, nH=1/(1-J). This allows designing specific models and very simple validations (at least at this coarse-grained level) and gives a new computational perspective by which analyze already developed ones (see, e.g., [58–62]). Then, regarding complex reactions, the puzzling scenario, evidenced by Koshland, finally finds out a simple descriptive framework that, crucially, also recovers to the standard cooperative scenario in the proper limit: full coherence among various, apparently antithetic, results is obtained within a unique framework.



---

*Source: 7423297-2018-01-29.xml*