Critical Re-evaluation of the Slope Factor of the Operational Model of Agonism Short Title: Slope Factor of the Operational Model of Agonism

Although being a relative term, agonist efficacy is a cornerstone in the proper assessment of agonist selectivity and signalling bias. The operational model of agonism (OMA) has become successful in the determination of agonist efficacies and ranking them. In 1985, Black et al. introduced the slope factor to the OMA to make it more flexible and allow for fitting steep as well as flat concentration-response curves. Functional analysis of OMA demonstrates that the slope factor implemented by Black et al. affects relationships among parameters of the OMA. Fitting of the OMA with Black et al. slope factor to concentration-response curves of experimental as well as theoretical data (homotropic allosteric modulation, substrate inhibition and non-competitive auto-inhibition) resulted in wrong estimates of operational efficacy and affinity. In contrast, fitting of the OMA modified by the Hill coefficient to the same data resulted in correct estimates of operational efficacy and affinity. Therefore OMA modified by the Hill coefficient should be preferred over Black et al. equation for ranking of agonism and subsequent analysis, like quantification of signalling bias, when concentration response curves differ in the slope factor.


Introduction
Since its introduction in 1983, the operational model of agonism (OMA) has become a golden standard in an evaluation of agonism and subsequently also of signalling bias 1 . The OMA describes the response of the system as a function of ligand concentration using three parameters: 1) The equilibrium dissociation constant of agonist (KA) to the receptor initiating functional response; 2) The maximal possible response of the system (EMAX); 3) The operational efficacy of agonist (τ). Original OMA describes the dependence of functional response on the concentration of ligand by rectangular hyperbola (Eq. 2). In practice, positive cooperativity or positive feedback lead to steep and negative cooperativity or negative feedback lead to flat concentration-response curves that the OMA does not fit to. Therefore, an equation intended for the description of non-hyperbolic concentration curves (Eq. 7), particularly for systems where the number of receptors varies, was introduced by Black et al. 2 . Since then, Eq. 7 is commonly used. The presented mathematical analysis of Eq. 7 shows that the slope factor affects the relationship between observed maximal response to agonist and operational efficacy and the relationship between the concentration of agonist for half-maximal response and its equilibrium dissociation constant. The presented case study shows that fitting Eq. 7 results in wrong estimates of model parameters.
The presented case study shows that the use of the Hill coefficient as a slope factor to the original rectangular hyperbola of OMA instead of the slope factor by Black et al. is a viable alternative. The Hill-Langmuir equation was originally formulated to describe the cooperative binding of oxygen molecules to haemoglobin 3 . The Hill coefficient reflects the sense and magnitude of cooperativity among concurrently binding ligands. A value of the Hill coefficient lower than 1 indicates negative cooperativity that is manifested as a flat concentration-response curve. The smaller value of the Hill coefficient the flatter curve is. A value of the Hill coefficient greater than 1 indicates positive cooperativity that is manifested as a steep concentration-response curve. The greater value of the Hill coefficient the steeper curve is. In contrast to Black et al. slope factor, the proposed implementation of the Hill coefficient (Eq. 10) does not affect the centre or asymptotes of the hyperbola describing concentration-response curves. Therefore, the Hill coefficient does not affect relationships among three parameters of the OMA, KA, EMAX and τ. This makes OMA (Eq. 11) more practical in many ways.

The general concept of the operational model of agonism
In general, OMA consist of two functions. One function describes the binding of an agonist to a receptor as the dependence of the concentration of agonist-receptor complexes [AR] on the concentration of an agonist [A]. The second function describes the dependence of functional response (E) to an agonist on [AR]. OMA expresses the dependence of E on [A].

Definition of OMA
In the simplest case, when both binding and response function are described by rectangular hyperbola (Supplementary Information Figure S1), the resulting function is also rectangular hyperbola. For example, in bi-molecular reaction, the dependence of ligand binding to the re-

Analysis of OMA
Eq.4 is the equation of OMA 1 . It has three parameters: The equilibrium dissociation constant of agonist (KA) at the response-inducing state of the receptor that is specific to a combination of ligand and receptor. The maximal possible response of the system (EMAX) is specific to the system. And the "transducer ratio" (τ) that is specific to a combination of ligand and system.
Eq. 4 is a rectangular hyperbola with the horizontal asymptote, the observed maximal response to agonist A (E'MAX), given by Eq. 5.

′ = +1
Eq. 5 A more efficacious agonist (having a high value of parameter τ) elicits higher E'MAX than less efficacious agonists (having a low value of parameter τ). Thus, τ is actually operational efficacy.
The relationship between parameter τ and E'MAX is hyperbolic meaning that two highly efficacious agonists (e.g., τ values 10 and 20) differ in E'MAX values less than two weak agonists (e.g., τ values 0.1 and 0.2).
In Eq. 4, the concentration of agonist A for half-maximal response (EC50), is given by Eq. 6.

= +1
Eq. 6 According to Eq. 6, for τ > 0, the EC50 value is always lower than the KA value. The KA to EC50 ratio is greater for efficacious agonists than for weak agonists. Similarly to E'MAX, the relationship between parameter τ and EC50 is hyperbolic. In contrast to E'MAX values, the ratio KA to EC50 ratio is more profound for two highly efficacious agonists (e.g., τ values 10 and 20) than two weak agonists (e.g., τ values 0.1 and 0.2).

Limitations of OMA
The OMA has several weak points. The major drawback of OMA is the lack of physical basis of the agonist equilibrium dissociation constant KA. In the Eq. 2 considers agonist binding [AR] refers to agonist binding to the conformation from which the functional response is initiated.
The agonist binding to the receptor in an inactive conformation is not observed in the response (KE→∞; τ = 0) In the radioligand binding experiments agonists bind to all receptor conformations including the inactive ones. For various reasons the receptors in the conformation activating functional response may be scarce or absent from radioligand binding experiments.
Then it may be impossible to determine the KA value in the radioligand binding experiments.
All three parameters of OMA (EMAX, KA and τ) are interdependent 4 . To fit Eq. 4 to the experimental data one of the parameters must be fixed. Thus, the maximal response of the system EMAX has to be determined before fitting Eq. 4. It can be archived by comparing functional response to a given agonist at system with reduced population of receptors by irreversible alkylation 5 . Another limitation of the OMA is that the shape of the functional response is a rectangular hyperbola.

Definition of non-hyperbolic OMA
In practice, concentration-response curves steeper or flatter than the ones described by Eq. 4 are observed. In such cases, Eq. 4 does not fit experimental data. As stated by the authors, Eq. 7 was devised for non-hyperbolic dependence of functional response on concentration of ago- Eq. 7

Analysis of non-hyperbolic OMA
Introduced power factor n changes slope and shape of functional response curve (Supplementary information Figure S3). Nevertheless, Eq. 7 as a mathematical function has rectangular Eq. 9 Evidently, the introduced slope factor n affects both the observed maximal response E'MAX and half-efficient concentration of agonist EC50 ( Figure 1A and C). The influence of the slope factor on E'MAX is bidirectional (Supplementary Information Table S1, Figure S2). For operational efficacies τ > 1, an increase in the value of slope factor increases E'MAX. ( Figure 1A and 2A blue lines). For operational efficacies τ < 1, an increase in slope factor decreases E'MAX ( Figure   1C and 2A yellow lines). The effect of slope factor on E'MAX is the most eminent for low values of operational efficacy τ, making the estimation of model parameters of weak partial agonists impractical. Imagine full agonist τ=10 and very weak agonist τ=0.1. For n=1: Full agonist E'MAX is 90 % and weak agonist E'MAX is 10 % of system EMAX. For n=2: Full agonist E'MAX is 99 % (one-tenth more) and weak agonist E'MAX is just 1 % (ten times less). An increase in the value of the slope factor increases the EC50 value ( Figure 2B). Again, the effect of the slope factor on the EC50 value is more eminent at low values of operational efficacy τ (red lines). Paradoxically, any combination of operational efficacy τ and slope factor fulfilling the inequality in Figure 2C (blue area) results in EC50 values greater than KA (e.g., Figure 1C, yellow lines). For example, EC50 > KA applies if τ = 0.5 and n > 1.6, or if τ = 1 and n > 1.6, or if τ = 1.5 and n > 2. 15. The upper asymptote of inequality is 2. Thus, the possibility of EC50 > KA applies to τ < 2.

Figure 2 Analysis of Black & Leff equation (Eq. 7)
A, Dependency of observed E'MAX to system EMAX ratio (ordinate) on slope factor n (abscissa) and operational efficacy τ (legend). B, Dependency of EC50 to KA ratio (ordinate) on slope factor n (abscissa) and operational efficacy τ (legend). C, Inequality plot of slope factor n (abscissa) and operational efficacy τ (ordinate) yielding half-efficient concentration EC50 greater than equilibrium dissociation constant KA.
The operational efficacy τ may be also considered as a measure of "receptor reserve". In the system with a relatively small capacity of a functional response output, the full agonist reaches its maximal response before reaching full receptor occupancy. Thus, the agonist EC50 value is lower than its affinity for the receptor. The smaller occupancy fraction is needed for the full response to a given agonist the greater is difference between agonist EC50 and KA values. According to OMA (Eq. 2), the relation between EC50 and KA is described by Eq. 6. The greater value of operational efficacy τ, the smaller EC50 value and the greater the difference from KA.
Thus, the value of operational efficacy τ is a measure of receptor reserve of a given agonist at a given system. In a system with a large capacity of functional output, agonists do not have a receptor reserve and must reach full receptor occupancy to elicit a full signal. In such a system, even full agonists have small operational efficacies. Thus, the parameter τ is specific to a combination of ligand and system.
Nevertheless, for agonists that elicit at least some response in a given system, the parameter τ must be greater than 0. Then according to Eq.6 of the operational model of agonism, the EC50 value must be smaller than the KA value. In principle, the EC50 value greater than the KA can be achieved only by some parallel mechanism that increases the apparent K'A, provided that a difference K'A to KA is greater that EC50 to KA. For example, such mechanism may be negative allosteric modulation of agonist binding or non-competitive inhibiton of functional response.

Limitations of the non-hyperbolic OMA
Besides all limitations of the hyperbolic OMA, the non-hyperbolic version of OMA has additional drawbacks. The most important is the lack of mechanistic background for factor n. Exponentiation of agonist concentration [A] to power factor n results in S-shaped curves. As it will be shown later, a combination of common binding mechanisms and functional responses does not result in S-shaped curves. Importantly, as shown above, exponentiation of operational efficacy τ to power factor n breaks the logical relationship between observed maximal response E'MAX and operational efficacy τ. That, as it will be shown later, impedes correct estimation of τ and KA values.

Definition of OMA with Hill coefficient
Hill coefficient may serve as an alternative slope factor in the OMA. Hill equation incorporates the Hill coefficient as a slope factor to rectangular hyperbola 3 . The major advantage of the Hill coefficient as a slope factor is that it allows for a change in the eccentricity (vertices) of the hyperbola-like curves without changing centre (EC50) and asymptotes (E'MAX) (Supplementary Information Figure S3). [ ] +( +1 ) Eq. 11 As expected, the Hill coefficient does not influence maximal observed response E'MAX either half-efficient concentration of agonist EC50 ( Figure 1B and D).

Implications of OMA with Hill coefficient
Analysis of the OMA with slope factor by Black et al. (Eq. 7) have shown that the slope factor n has a bidirectional effect on the relationship between the parameters E'MAX and τ. and that the slope factor n affects the relationship between the parameters EC50 and KA. In contrast, in Eq. 11 neither value of E'MAX nor the value of EC50 is affected by the Hill coefficient ( Figure  1B and D). The parameters E'MAX and EC'50 can be readily obtained by fitting Eq. 10 to the single concentration-response data.

Limitations of OMA with Hill coefficient
The major criticism of the Hill equation is its parsimonious character. It is relatively simple and its parameters are easy to estimate. However, as a model, it is just an approximation. In an experiment, the slope of the concentration-response curve different from unity may be a result of the parallel signalling mechanism providing feedback or allosteric cooperativity. E.g., positive cooperativity results in steep concentration-response curves and negative cooperativity results in flat concentration-response curves.

The case study
How fitting the Black & Leff equation to experimental data can affect estimates of the operational efficacy and subsequent analysis is demonstrated on the following example of measurement of the GTPγS binding as a functional response of M2 receptor to muscarinic agonists carbachol and oxotremorine (Figure 3). At Gi1 G-protein (Figure 3, left) carbachol stimulated the GTPγS binding more than oxotremorine. In contrast, at GoA G-protein oxotremorine stimulated GTPγS binding slightly more than carbachol (Figure 3, right). While the functional response to carbachol has normal slope (nH ≈ 1), the response to oxotremorine is shallow (nH = 0.84 and 0.64, respectively) ( Table 1) Table  1. Figure 3. Black&Leff (Eq. 7) and Hill equation (Eq. 10) were fitted to concentration-response data of GTPγS binding to G-proteins upon stimulation of M2 receptors by carbachol or oxotremorine in Figure 3. The EMAX was fixed to 1.

OMA of homotropic allosteric of binding
The simplest mode of interaction that leads to variation in the slope of concentration-response curves is allosteric modulation. Let's take for example homotropic allosteric modulation and suppose that ligand A binds to the orthosteric (AR) and the allosteric binding site (RA) with the same affinity. For the sake of simplicity suppose that a ligand A activates the receptor only from the orthosteric binding site. And that A bound to the orthosteric binding site (AR) has the same operational efficacy as A in the ternary complex (ARA) (Figure 4). In such a case, the amount of active complexes AR and ARA is given by Eq. 12.
As can be seen in Figure 3

Table 2 Results of fitting Black & Leff and Hill equations to the model of homotropic allosteric modulation of binding
Balck & Leff (Eq. 7) and Hill (Eq. 10) were fitted to model data (Eq. 13) with EMAX fixed to 1. Apparent pK'A values were calculated from model Eq. 13 according to Eq. 14 or from EC50 values according to Eq. 11.
Fitting Eq. 13 with fixed system EMAX to the model of functional response homotropic allosteric modulation of binding yields correct parameter estimates that are associated with the low level of uncertainty (Supplementary Information Figure S6).

OMA of homotropic allosteric modulation of efficacy
Homotropic allosteric modulation of agonist binding may result in the change in operational efficacy of agonist τ ( Figure 5). If the operational efficacy of AR is τ and the operational efficacy of ARA is βτ (KE/β) then the apparent activation constant K'E changes with the concentration of an agonist according to Eq. 16.

Eq. 16
Where β is the cooperativity factor of efficacy. Values of β greater than 1 denote positive cooperativity, a decrease in apparent K'E. Values of β smaller than 1 denote negative cooperativity, an increase in apparent K'E.
For the sake of simplicity suppose neutral cooperativity (α = 1) in the binding of A to the orthosteric and allosteric sites. The EC50 value in homotropic allosteric modulation of efficacy is given by Eq. 20 (Supplementary information Eq. A42).
This should not represent a problem as for neutral cooperativity (β=1) Hill coefficient equals one and thus OMA in the form of Eq. 4 can be used.

OMA of substrate inhibition
Another common mechanism affecting the functional response to an agonist is substrate inhibition. In substrate inhibition, substrate concentration-dependently inhibits the reaction. In the case of functional response to receptor activation, agonist-receptor complexes [AR] are the substrate ( Figure 6). Functional response is given by Eq. 21. Eq. 23
Eq. 24 The EC50 value for the model of substrate inhibition is given by Eq. 25 (Supplementary information Eq. A53).
Eq. 25 The shape for functional response is similar to homotropic allosteric modulation of efficacy   Table 4.  Figure S8).

OMA of non-competitive inhibition
As shown in Figure 2, OMA with slope factor n allows for EC50 values higher than KA. The simplest mode of interaction that increases observed EC50 above KA is non-competitive autoin- Eq. 28 Thus, apparent operational efficacy τ' is given by Eq. 29.
Eq. 29 The EC50 value for the model of non-competitive auto-inhibition is given by Eq. 30 (Supplementary information Eq. A61).   Table 5.  Figure S9 and S10). In the case of KI=5 (Supplementary Information Figure S9), estimates of operational efficacy τ and inhibition factor σ are swapped pointing to the symmetry of Eq. 27. This symmetry makes calculation of τ and σ impossible as any τ and σ combination resulting in an appropriate apparent efficacy τ' (Eq. 28) fits well a given functional-response data (Supplementary Information Figure S10).

Discussion
The operational model of agonism (OMA) 1 is widely used in the evaluation of agonism. The OMA characterizes a functional response to an agonist by the equilibrium dissociation constant of agonist (KA), the maximal possible response of the system (EMAX) and the operational efficacy of agonist (τ) (Eq. 4). To fit non-hyperbolic functional responses slope factor n was introduced to the OMA (Eq. 7) 2 . Analysis of the Black & Leff equation (Eq. 7) has shown that the slope factor n has a bidirectional effect on the relationship between the parameters E'MAX and τ ( Figure 1A versus 1C) and also affects the relationship between the parameters EC50 and KA.
In practice, as exemplified in Figure 3, Despite the dire effects of slope factor n, the Black & Leff equation is widely accepted 7,10-16 .
It even entered textbooks 17 . Very little concern on factor n has risen. For example, Kenakin et al. 9 analysed in detail the effects of slope factor n on EC50 and τ to KA ratio but did not deal with the bi-directional effect of n on τ neither proposed an alternative approach to avoid potential pitfalls. To force a proper shape on functional-response curves, Gregory et al. 18 Figure 1B and D). Therefore, biased signalling may be inferred from the comparison of the ratio of intrinsic activity (E'MAX/EC50) of tested agonist to the intrinsic activity of reference agonist at two signalling pathways as in the case of Hill equation the E'MAX/EC50 ratio is equivalent to τ / KA ratio 6,20 . As shown for experimental data in Figure 3 and Table 1 shown that the Hill coefficient does not affect the relationship between the parameters E'MAX and τ neither between the parameters EC50 and KA. Fitting Hill equation to the concentrationresponse data gives good estimates of EC50 and E'MAX values, which are suitable for further analysis of OMA, OMAMA and signalling bias, e.g. using relative intrinsic activities.
USA). Plates were dried in microwave oven 800W 3 min and then 40µl of ROTISZINT® Eco Plus (ROTH) was added. The plates were counted in the Wallac Microbeta scintillation counter.