## Abstract

In classical pharmacology, bioassay data are fit to general equations (e.g. the dose response equation) to determine empirical drug parameters (e.g. EC_{50} and E_{max}), which are then used to calculate chemical parameters such as affinity and efficacy. Here we used a similar approach for kinetic, time course signaling data, to allow empirical and chemical definition of signaling by G-protein-coupled receptors in kinetic terms. Experimental data are analyzed using general time course equations (model-free approach) and mechanistic model equations (mechanistic approach) in the commonly-used curve-fitting program, GraphPad Prism. A literature survey indicated signaling time course data usually conform to one of four curve shapes: the straight line, association exponential curve, rise-and-fall to zero curve, and rise-and-fall to steady-state curve. In the model-free approach, the initial rate of signaling is quantified and this is done by curve-fitting to the whole time course, avoiding the need to select the linear part of the curve. It is shown that the four shapes are consistent with a mechanistic model of signaling, based on enzyme kinetics, with the shape defined by the regulation of signaling mechanisms (e.g. receptor desensitization, signal degradation). Signaling efficacy is the initial rate of signaling by agonist-occupied receptor (*k*_{τ}), simply the rate of signal generation before it becomes affected by regulation mechanisms, measurable using the model-free analysis. Regulation of signaling parameters such as the receptor desensitization rate constant can be estimated if the mechanism is known. This study extends the empirical and mechanistic approach used in classical pharmacology to kinetic signaling data, facilitating optimization of new therapeutics in kinetic terms.

## Introduction

Pharmacological data analysis provides the drug activity parameters used to optimize the biological activity of new therapeutics and to identify mechanisms of receptor-mediated signal transduction for G-protein-coupled receptors (GPCRs). In the typical analysis process, activity is measured at various drug concentrations and the data fit to a concentration–response equation by curve fitting, for example to a sigmoid curve equation^{1,2,3}. This yields empirical drug parameters, usually the potency (EC_{50}) and a measure of the maximal signaling capacity (E_{max}). This analysis is described as “Model-free” because it makes minimal assumptions regarding the mechanism of signal transduction. These parameters can then be used to calculate mechanistic drug parameters that define the effect of the drug on the system in chemical terms (e.g. efficacy and macroscopic affinity)^{3,4,5,6,7,8}. These chemical parameters aid the translation of drug effect measurement from in vitro test systems to animal models and human disease, for example by minimizing the effect of tissue-specific parameters on the quantification of drug activity. Biased agonism measurement is an example of this process currently in use by many laboratories. Signaling assay data are fit to concentration–response equations and the empirical parameters EC_{50} and E_{max} used to calculate mechanistic parameters of biased signaling, such as transducer coefficients and log efficacy ratios (reviewed in^{9}).

The manner in which GPCR signaling changes over time, i.e. the kinetics/dynamics of response, is currently of considerable interest in quantitative pharmacology. The time course of GPCR response has been measured since the earliest studies of muscle contraction^{1}, through the discovery of G-protein-mediated signaling^{10}, the application of Ca^{2+} indicator dyes^{11}, and the most recent advances in biosensor technology that enable high-throughput signaling kinetic analysis^{12,13,14}. Measuring the time course of receptor signaling revealed fundamental mechanisms of GPCR function. The classic rise-and-fall time course of cAMP in S49 cells indicated desensitization of β_{2} adrenoceptor signaling^{15}. Sustained signaling after agonist washout by parathyroid hormone, sphingosine 1-phosphate and thyroid-stimulating hormone receptors revealed signaling by internalized receptors^{16,17,18}, a new spatiotemporal paradigm of GPCR signaling of potential therapeutic utility^{19,20,21,22,23}. Signaling kinetics can also affect quantification of drug activity. Of potential concern, drug effect can be dependent on the time point at which the response is measured. For example, biased signaling of the D_{2} dopamine receptor was shown to be time-dependent^{24}. For the CB_{1} receptor, the potency for internalization increased over time^{25}. For AT_{1} angiotensin receptor-mediated arrestin recruitment the potency increased over time, as did the E_{max} of partial agonists^{26}. Finally, differences of signaling kinetics can be a manifestation of structural differences of receptor-ligand interaction. For example, the nonpeptide GLP-1 receptor agonist TT-OAD2 protrudes outside of the receptor core, interacting with lipid, potentially explaining the unusually slow signaling kinetics of the ligand^{27}.

A data analysis framework for time course data of GPCR signaling would enable the quantification of efficacy in kinetic terms. One way to do this is to measure the initial rate of signaling, analogous to the initial rate of enzyme activity^{26,28,29,30,31,32}, an approach recently applied to receptor receptor–arrestin interaction^{26}. Kinetic analysis could also enable measurement of regulation of signaling parameters, such as the rate constant for receptor desensitization. At present, time course data for GPCR signaling are rarely analyzed by curve fitting to estimate kinetic pharmacological parameters (in contrast to the near-universal application of curve fitting to concentration-dependence data). The absence of curve fitting to time course data might be resulting in lost opportunities, particularly for biosensor assay modalities in which the time course data are collected by default. Kinetic insights into GPCR function might be being missed. Potential benefits of the kinetic dimension to ligand optimization could be going unrealized. The goal of this study was to develop a straightforward data analysis framework for quantifying the kinetics of GPCR signaling, enabling investigators to measure kinetic drug parameters useful for practical pharmacological application.

## Results

In this study, a data analysis framework for curve fitting of time course data was developed for GPCR signaling kinetics/dynamics. First, we identified the types of time course curve shape by conducting a comprehensive literature survey. It was discovered that the curve shapes, more precisely the temporal profiles, conform to a limited number of shapes (four). This survey also showed how the time course is shaped by regulation of signaling mechanisms (e.g. receptor desensitization and response degradation). The shapes were defined by simple equations that can be used to analyze time course data using familiar curve-fitting software (e.g. Prism, GraphPad Software, Inc., San Diego CA). This analysis yields a kinetic drug parameter, the initial rate of signaling by agonist-occupied receptor, that defines efficacy as the rate of signal generation before it is impacted by regulation of signaling mechanisms. This parameter, termed *k*_{τ}, provides a biologically meaningful and intuitive metric of the kinetics of signal transduction. Finally, mechanistic models are applied to the time course data using a kinetic mechanistic model of agonism developed previously^{33}, and extended here to incorporate receptor desensitization and sustained signaling by internalized receptors. This analysis demonstrated the model-free equations emerge as general forms of the mechanistic equations, providing a mechanistic foundation for the analysis approach.

### Literature survey of time course shapes

The first step in analyzing biological activity data is visual inspection of the data. In pharmacological data analysis, this led to the realization that concentration–response data are usually described by a sigmoid curve (when the x axis is the logarithm of ligand concentration)^{1,2}. Here we surveyed the GPCR signaling literature for time course data. The survey was designed to be comprehensive, spanning (1) The full range of heterotrimeric G-protein classes (G_{s}, G_{i}, G_{q/11}, G_{12/13}); (2) response times of milliseconds to hours; (3) Proximity to the receptor, from direct receptor interaction (G-protein and arrestin) to downstream signals (e.g. cell contraction/relaxation and gene expression); (4) Types of response, including chemical (second messengers), mechanical (cell structure change), electrical (ion channel currents) and genetic (gene expression). This was done for experiments when the receptor was exposed continuously to the agonist, i.e. agonist washout experiments were not considered.

Within this survey, the large majority of time course profiles conformed to one of four shapes (Fig. 1). Each shape is defined by a corresponding equation (Eqs. 1–4 in Methods). The equations are amenable to straightforward curve-fitting analysis in familiar software. The shapes are, in order of increasing complexity:

- 1.
The straight line (Fig. 1a). Response increases continuously over time at a constant rate (defined by Eq. 1).

- 2.
The association exponential curve (Fig. 1b, Eq. 2). Initially, response increases rapidly, being almost linear. The response then slows then finally approaches a plateau at which the response level remains constant over time.

- 3.
The rise-and-fall to baseline curve (Fig. 1c, Eq. 3). The first phase resembles the association exponential curve – response increases rapidly at first, then slows. The response then reaches a peak level. Subsequently, the response level decreases and finally falls back to the baseline level, i.e. the level before addition of ligand.

- 4.
The rise-and-fall to steady-state curve (Fig. 1d, Eq. 4). Response rises and falls as for the rise-and-fall to baseline curve. However, the response declines to a steady-state level, above that of the baseline level before addition of ligand.

#### Straight line time course profile

The straight line profile was evident in second messenger responses under a specific condition—when signaling was unregulated, i.e. when regulation of signaling mechanisms were blocked (Fig. 2). GPCR signaling is regulated over the short term by two primary mechanisms: receptor desensitization (involving receptor phosphorylation and subsequent arrestin binding)^{34,35,36} and degradation of the response (for example metabolism of second messengers)^{37,38,39}. The straight line profile was observed when both mechanisms were blocked. For example, AT_{1} angiotensin receptor-stimulated inositol phosphates (IP) production was linear when desensitization was blocked (using arrestin knock out cells) and response degradation was blocked (using Li^{+} to block IP breakdown) (Fig. 2a)^{40}. The same result was obtained for the proteinase-activated receptor PAR1 receptor (see Fig. 2a in^{41}). A second example is provided by the GnRH_{1} receptor, which lacks a C-terminal tail and so does not interact with arrestin^{42}. A linear time course was observed, in the presence of Li^{+} (Fig. 2b)^{43}. (See also Fig. 1a,b in^{44} and Fig. 3 in^{45}.) A third example is provided by the glucagon receptor. A linear time course of cAMP accumulation was observed in hepatic membranes^{46}. The membranes lacked cAMP phosphodiesterase activity^{47}, and likely lacked receptor desensitization components owing to extensive washing of the preparation (Fig. 2c).

Linear time course data were also observed for long duration responses, far downstream, at the level of DNA, specifically gene expression and DNA synthesis (Supplementary Fig. S1)^{48,49}. In these responses there was a delay before the onset of the response. This has been rationalized for gene expression responses by a need for a build-up of signal transduction intermediates to a threshold level that then initiates the process^{48}.

#### Association exponential time course profile

The association exponential profile was the most commonly-observed time course shape, especially for second messenger molecules such as cAMP and inositol phosphates, but also for a variety of other signals, from upstream events such as arrestin recruitment to downstream cellular functions such as muscle cell relaxation (Fig. 3, Supplementary Fig. S2). Regarding regulation of signaling, the profile was evident when one of the two mechanisms was operative and one was blocked (the mechanisms being receptor desensitization and response degradation). Examples include AT_{1} receptor-stimulated IP production, with receptor desensitization in operation (arrestin wild-type cells) but without response degradation (IP metabolism blocked by Li^{+}) (Fig. 3a)^{40}; GnRH_{1} receptor-simulated IP production with response degradation in operation (no Li^{+} present) but without receptor desensitization (the GnRH receptor lacking a C-terminal tail) (Fig. 3b)^{45}; and β_{2} adrenoceptor-stimulated cAMP generation with response degradation (cAMP metabolism) but without receptor desensitization (arrestin knockout cells) (Fig. 3c)^{53}.

The association exponential profile was observed in direct receptor-transducer interaction assays, including receptor-G-protein interaction and receptor-arrestin interaction (Supplementary Fig. S2a,b)^{26,54}. This was expected; mechanistically, the response is likely resultant from a bimolecular interaction and the association exponential equation is a general form of the familiar bimolecular association kinetics equation. G-protein activation was also described by the profile, measured by [^{35}S]GTPγS binding, or by a G-protein activation biosensor (Supplementary Fig. S2c,d)^{54,55}. cAMP production responses also conform to the profile, both stimulation (Supplementary Fig. S2e) and inhibition (Supplementary Fig. S2f)^{56,57}. A specialized signaling pathway, β-catenin stabilization via the Wnt-Frizzled system, was also described by the profile^{58} (Supplementary Fig. S2g). Electrical signaling was also described by the association exponential profile, as shown by GIRK channel gating in *Xenopus* oocytes (Supplementary Fig. S2h)^{59}. In the final example, a downstream cellular response was found to conform to the profile, relaxation of human airway smooth muscle cells via the β_{2} adrenoceptor (Supplementary Fig. S2i)^{60}.

#### Rise-and-fall to baseline time course profile

The rise-and-fall to baseline profile is a classic curve shape in pharmacology, leading to the definition of “Fade” (decline in the response to a continuous application of agonist^{62}) and the discovery of the underlying mechanism (regulation of signaling processes, especially receptor desensitization)^{15,36}. An equation was identified that defines this shape, here termed the rise-and-fall to baseline equation (Eq. 3). This equation is familiar in pharmacokinetics, being the equation defining drug concentration in the oral dosing absorption and elimination model^{63,64}. We have loaded it into a custom Prism template, “Signaling kinetic model-free equations” available in the supplementary files.

Two examples for second messenger molecules are shown in Fig. 4, diacylglycerol production via the AT_{1} receptor (Fig. 4a)^{65} and cAMP production via the β_{2} adrenoceptor^{53} (Fig. 4b). In both studies the mechanisms underlying the curve shape were investigated. The rise-and-fall to baseline mechanism was evident when both regulation of signaling mechanisms were operative, i.e. when there were no experimental manipulations of receptor desensitization or degradation of the response. When the regulation mechanisms were manipulated, the shape of the curve changed. When receptor desensitization was blocked or attenuated, the resulting curve shape approached the association exponential curve, as shown in Fig. 2a of the original study for the AT_{1} receptor^{65} and Fig. 3c for the β_{2} adrenoceptor. Response degradation was in operation, demonstrated for the β_{2} adrenoceptor by the effect of phosphodiesterase inhibition^{53} and assumed for the AT_{1} receptor because diacylglycerol is typically cleared rapidly by diacylglycerol kinases^{38}.

The rise-and-fall to baseline response is also well known in calcium signaling, representing the change of cytoplasmic Ca^{2+} concentration on stimulation by GPCR agonists (usually when there is no extracellular Ca^{2+} in the assay). An unusual regulation mechanism is involved in Ca^{2+} signaling via intracellular stores. The amount of cytoplasmic Ca^{2+} that can be obtained is limited by the amount of Ca^{2+} in the intracellular stores^{66,67}. Depletion of Ca^{2+} in the store as it is being released into the cytoplasm limits the amount and rate of cytoplasmic Ca^{2+} release. This regulation mechanism can be described as depletion of response precursor^{33}. The response is also regulated by clearance of cytoplasmic Ca^{2+} out of the cell^{68,69}. This can be described as response degradation. These two regulation processes result in the rise-and-fall to baseline profile as described previously^{33}. An example is provided by GnRH-stimulated Ca^{2+} mobilization in pituitary gonadotrophs, as shown in Fig. 4c^{70}. Other receptors giving this profile include the BB1 receptor^{71} and AT_{1} receptor^{26}. It is important to note that other time course shapes are observed for cytoplasmic Ca^{2+} mobilization, including the rise-and-fall to steady-state curve (considered separately below), and calcium oscillations and waves^{72,73}. In addition, rise-and-fall profiles are observed in Ca^{2+} responses which are not fit well by the rise-and-fall equations used in this study (e.g. in Supplementary Fig. S3).

#### Rise-and-fall to steady-state time course profile

A second rise-and-fall curve shape is encountered in GPCR signaling in which the response, after rising to the peak, falls to a plateau level of response which is above baseline, i.e. response declines to a steady-state after peaking. Examples of this shape are shown in Fig. 5. In this study, an equation was identified that describes these data, termed the rise-and-fall to steady-state equation (Eq. 4). This equation was identified as a general form of explicit GPCR signaling mechanism equations, as described in Appendix 2. These mechanisms are more complex than the simplest mechanisms considered previously^{33} and include recently-discovered mechanisms, such as persistent signaling by internalized receptors^{16,17}. These mechanisms are considered below (“More complex models”). We have loaded the equation into a custom Prism template, “Signaling kinetic model-free equations” available in the supplementary files.

The most commonly-encountered example was calcium mobilization, the increase of cytoplasmic Ca^{2+} upon application of the GPCR agonist. A representative example is shown in Fig. 5a, cytoplasmic Ca^{2+} concentration stimulated by GnRH via the GnRH_{1} receptor^{70}. It is well known that the cytoplasmic Ca^{2+} level usually reaches a plateau after sustained application of the agonist, particularly when Ca^{2+} is included in the extracellular medium. Under this condition, Ca^{2+} can re-enter the cell and this can result in a steady-state being reached between the numerous processes controlling cytoplasmic Ca^{2+} concentration^{68,72,73,75}. This mechanism can be represented in the context of the regulation of signaling mechanisms as a reformation of the response from the response decay product, in this case reappearance of cytoplasmic Ca^{2+} from the extracellular space. This mechanism was formulated in Appendix 2.3.3 (see supplementary files). This mechanism is described in general form by the rise-and-fall to steady-state equation (Eq. 3) and in explicit form by Eq. (15). Data are analyzed by the general form (see “Measuring the initial rate from curve fit parameters” below) and explicit form (“Measuring model parameters”).

The rise-and-fall to steady-state profile is also commonly observed in numerous downstream signaling responses, examples of which are shown in Fig. 5: ERK1/2 phosphorylation by lysophosphatidic acid^{76}, and protein kinase C and Rho activation via the AT_{1} receptor^{77}. In some cases, the data fit the rise-and-fall to steady-state equation well (high R^{2} value) but the fitted parameters were ambiguous, a scenario encountered when the two rate constant values were almost equal. An example is provided in Supplementary Fig. S4. This observation requires further investigation, for example using structural identifiability analysis.

### Initial rate measurement using model-free analysis

In order for time course analysis to be useful for pharmacological investigations, pharmacological parameters need to be extracted from the data that capture the temporal dimension of activity. One kinetic parameter used routinely for another target class (enzymes) and occasionally for GPCRs is the initial rate of activity. This is the rate at the earliest part of the time course, when the response increases linearly over time. This parameter offers several benefits. (1) It is a biologically-meaningful descriptor of signaling activity. (2) In principle, the initial rate is effectively a pure efficacy parameter, being the rate of signal generation before it is affected by regulation of signaling mechanisms, which affect the shape of the later part of the time course. This is formally demonstrated using the kinetic pharmacological model below. (3) From a practical perspective, the initial rate provides a single parameter, as opposed to multiple parameters, of drug efficacy, simplifying application to ligand optimization in medicinal chemistry projects. (4) The initial rate can be determined regardless of the shape of the time course, enabling comparison of a ligand’s efficacy between responses with different temporal profiles (see below). We recently applied the initial rate approach to arrestin–receptor interaction, a special case in which there was no signal transduction^{26}. Here it is applied universally to GPCR signaling responses. It is important to note the initial rate is system-dependent because it is dependent on system parameters such as the receptor and effector concentrations (see “Comparison of models with model-free analysis” below).

#### Measuring the initial rate from curve fit parameters

For curved time course data, the initial rate is often measured by assessing which data points lie on the linear portion at the start of the curve, then fitting a straight line equation to those points. Here a more efficient method is developed. The entire time course is fit to the equation that defines the curve, then the initial rate calculated from the fitted parameters. This requires an additional equation that defines the initial rate (\(\mathrm{IR}\)) in terms of the fitted parameters. This equation was obtained here for each of the time course shapes, as the limit of the time course equation as time approaches zero (the formal definition of the initial rate condition) (Appendix 1, see supplementary files):

Straight line, from Eq. (1):

Association exponential, from Eq. (2):

Rise-and-fall to baseline, from Eq. (3):

Rise-and-fall to steady-state, from Eq. (4):

This method was used to quantify the initial rate of signaling for the responses in Figs. 2, 3, 4 and 5 and Supplementary Figs. S1, S2. Data were fit to the appropriate time course equation using Prism 8. The fitted parameter values and the fitted standard error are shown in the panels of the figures. For linear and association exponential profiles data were analyzed using built-in equations provided with the software. For the rise-and-fall profiles, user-entered custom equations were used. (A ready-to-use Prism template containing these equations, “Signaling kinetic model-free equations,” is provided in the supplementary files.) In some cases, the signal initiation time was a variable in the equation (see Eqs. 5–9), allowing for delay between addition of ligand and initiation of signal (e.g. in gene expression assays, Supplementary Fig. S1), or accommodating the scenario where the ligand addition time is not precisely defined. The data presented in this study were fit well by the equations; in all but one case, the correlation coefficient was > 0.95, and the standard error of the fit parameters was typically < 10% of the fitted value (see “Curve fit results” Excel file in Supporting Material). With two exceptions, data were fit well using the default fitting method of Prism 8 (least squares regression with medium convergence criteria^{78}). The exceptions were the rise-and-fall fit in Fig. 4b (cAMP production via the β_{2} adrenoceptor without phosphodiesterase inhibition) and Fig. 5b (LPA-stimulated ERK phosphorylation). In these cases, the default method yielded an ambiguous fit^{79} (see “Curve fit results” Excel file in Supporting Material), which was resolved using the “Robust regression” option which minimizes the contribution of outliers to the fit^{74}.

The initial rate value was then calculated from the fitted parameter estimates using the equations above. Note the simplicity of obtaining the initial rate in most cases. For example, for the association exponential curve, the initial rate is simply the rate constant multiplied by the steady-state response. For the rise-and-fall to baseline curve, it is simply the value of the fitted parameter \(C\). The calculated initial rate values are shown in the figures (Figs. 2, 3, 4 and 5, Supplementary Figs. S1, S2). A maximally-stimulating concentration of agonist was employed in all the examples provided. Under this condition, the initial rate is the efficacy parameter of the agonist in kinetic terms, here termed \({\mathrm{IR}}_{\mathrm{max}}\).

It is also possible to perform global fitting of the signaling kinetic data with time and ligand concentrations as independent variables. This can require knowledge of which parameters in the model are dependent on the ligand concentration, which in turn can require knowledge of the underlying mechanism. The methods used for this approach are described in ref^{33}.

#### Comparison with single time point analysis

Historically, pharmacological responses have been measured for multiple concentrations of ligand at a single time point and ligand activity described by the potency (e.g. EC_{50}) and maximal response (E_{max}). Here we compared this single time point concentration–response analysis with a concentration–response analysis of the initial rate of signaling. For this purpose, we measured cAMP accumulation and arrestin recruitment via the V_{2} vasopressin receptor. This receptor is activated by two cognate endogenous ligands, vasopressin and oxytocin^{80}. The receptor is well known to couple stably to arrestin when bound by vasopressin (a so-called class B arrestin binding profile^{34}). The receptor is activated by two endogenous ligands, vasopressin and oxytocin, that differ in their spatiotemporal mechanisms of signaling^{81}. Genetically-encoded biosensors were used to quantify these responses in real time in live cells expressing the receptor^{26,82}, as described in Methods.

Figure 6 shows the time courses for a range of concentrations of the two ligands for cAMP accumulation and arrestin recruitment. In all cases, the time course data were fit well by the association exponential equation (Eq. 6, Fig. 6, see also “Curve fit results” Excel file in Supporting Material). From the fitted parameters, the initial rate was calculated. This was done by multiplying the steady-state response by the rate constant as described above (see “Curve fit results” Excel file for values). The initial rate for maximally-stimulating concentrations of ligand is shown graphically by the dashed straight lines in Fig. 6c,f. The concentration–response of the initial rate was then evaluated, as shown in Fig. 7a. The concentration response was fit well by a standard sigmoid curve equation, from which the efficacy and potency of the ligand could be determined (Table 1). The efficacy is the maximal initial rate (\({\mathrm{IR}}_{\mathrm{max}}\)), i.e. the upper plateau of the curve. The potency is the midpoint concentration of the ligand, specifically the concentration of ligand giving an initial rate half of the maximum value. Here this is termed L_{50} (Table 1).

The concentration–response analysis indicated oxytocin was less active for recruiting arrestin than vasopressin (Fig. 7a, Table 1), consistent with a previous report^{81}. This was manifest as a 6.2-fold lower potency for oxytocin compared with vasopressin (580 nM vs 93 nM respectively, − log L_{50} values significantly different, Table 1). The maximal initial rate of arrestin recruitment (\({\mathrm{IR}}_{\mathrm{max}}\)) for oxytocin was value 70% of the vasopressin value but the \({\mathrm{IR}}_{\mathrm{max}}\) values were not significantly different (Table 1). These results indicate the weaker arrestin recruitment by oxytocin can be accounted for by a lower affinity of oxytocin compared with vasopressin for the arrestin-bound receptor. By contrast the two ligands were effectively equivalent for stimulation of cAMP production (Fig. 6c, Table 1). These findings imply biased agonism of oxytocin relative to vasopressin (G-protein bias) in terms of the initial rate.

Next, data were analyzed using a single time point. A time late in the time course was selected when the responses for high concentrations was approaching steady-state (18 min after ligand addition, 22 min after the read start, Fig. 6). This reflects a modality in which the assay is left long enough for the maximal signal to accumulate, often referred to as an “Endpoint” modality. The concentration–response curves are shown in Fig. 7b and results of the sigmoid curve analysis in Table 2. Using this metric, the parameters were slightly different from the initial rate analysis (Fig. 7a, Table 1). Specifically, the maximal response of oxytocin was similar to that of vasopressin for both arrestin and cAMP responses. In addition, rather than being uquipotent in the cAMP assay, vasopressin was more potent than oxytocin (L_{50} of 0.089 and 0.36 nM, − log L_{50} values significantly different, Table 2). This means that the biased agonism in terms of potency observed in the initial rate analysis for oxytocin was not observed in the endpoint modality (since the potency difference for arrestin, 4.0-fold, was similar to the potency difference for cAMP, 6.9-fold, Table 2). One possible explanation for this difference is that the initial rate for vasopressin in generating cAMP was limited by slow receptor-ligand equilibration, a scenario discussed in theoretical terms previously^{33} and in the Discussion below. Another difference was that the potency for arrestin recruitment was higher using the endpoint approach compared with the initial rate approach (by 4.4-fold for vasopressin and 3.2-fold for oxytocin, compare Tables 1 and 2). This difference was also observed for the AT_{1} angiotensin receptor and reflects the kinetic mechanism of arrestin recruitment, as described in ref.^{26}.

### Mechanistic pharmacological model of GPCR signaling and regulation kinetics

Recently a mechanistic kinetic pharmacological model was developed to quantify GPCR signaling activity in kinetic terms^{33} that is being applied to understanding signaling efficacy in kinetic terms^{8}. The model comprises a minimal number of parameters to enable parameter estimation by curve fitting. The approach is macroscopic—ligand signaling activity is defined at the level of the whole system, like the operational model of agonism^{7}. (Numerous microscopic systems-biology type models have been developed to simulate GPCR signaling dynamics that quantify each interaction in the signaling cascade^{8,53,84,85,86,87,88} but, while they enable close examination of mechanism, they contain too many parameters to be useful in medicinal chemistry campaigns.) The model is illustrated in Fig. 8. Signaling activity is quantified using a single, readily-measurable parameter, *k*_{τ}, which is the initial rate of signaling by the agonist-bound receptor in enzymatic terms^{33}. The model is extended to incorporate known regulation of signaling mechanisms. GPCR signaling is regulated to limit signaling, preventing over-stimulation of the cell. We consider the properties of the model, how it relates to the model-free analysis, and how to estimate the efficacy and regulation parameters by curve fitting. The model gives rise to the four curve shapes observed experimentally. It emerged that the specific shape is simply dependent on the number of regulatory mechanisms (0, straight line; 1, association exponential; 2, rise-and-fall to baseline). It is shown the model-free analysis emerges from the mechanistic model, the equations of the former being general forms of the equations for the latter. Estimating efficacy is straightforward and does not require knowledge of the mechanism—it is shown *k*_{τ} is equal to \({\mathrm{IR}}_{\mathrm{max}}\). The regulation parameters can be estimated if the mechanism is known, for example the receptor desensitization rate constant.

#### Model description

The kinetic mechanistic model describes GPCR signaling in terms of enzyme activity. A signal results from conversion of a signal precursor (analogous to the substrate) to the signal (the product) by the agonist-bound GPCR (the enzyme) (Fig. 8), as described previously^{33}. The rate of ligand efficacy for generating the response is quantified as the initial rate of signaling, termed *k*_{τ}, analogous to the initial rate of enzyme activity. *k*_{τ} is the rate of response generation by the agonist-occupied receptor, defined as the product of the total receptor concentration (\({[R]}_{TOT}\)), total signal precursor (\({E}_{P})\) and the response-generation rate constant (\({k}_{E})\):

(By comparison, the initial rate of enzyme activity is the product of the enzyme concentration, substrate concentration and the catalytic rate constant.) The model as formulated in this study assumes receptor-ligand occupancy does not change over time, that the equilibrium level of occupancy is rapidly achieved and is not rate-liming. This condition is likely met for maximally-stimulating ligand concentrations (used to quantify efficacy), and for all concentrations of lower potency ligands, scenarios which result in rapid association of ligand with receptor (see Discussion). Note the model framework is amenable to extension to incorporate the kinetics of receptor-ligand binding, as described previously^{33}.

The model can be extended to incorporate regulation of signaling. The canonical mechanisms are response degradation and receptor desensitization. Response degradation is the process by which the signal, once generated, is cleared over time. Examples include breakdown of second messenger molecules, such as cAMP by phosphodiesterases^{37}, and de-activation of G-protein by hydrolysis of bound GTP to GDP by the intrinsic GTPase activity of the G-protein^{89}. Some signals are decreased by clearance of the signaling species from the relevant compartment, for example efflux of cytosolic Ca^{2+} ions^{69}. In the model, response degradation (pink region of Fig. 8) is represented simply by exponential decay, governed by \({k}_{D}\), the response degradation rate constant. This component of the model was described previously^{33}.

#### Receptor desensitization model

In the present study, the model is extended to incorporate the second canonical regulation process, receptor desensitization (Fig. 8 yellow region). Receptor desensitization typically results from receptor phosphorylation and subsequent arrestin binding, which inhibits G-protein interaction^{34,35,36}. This process is represented simply as an exponential decay of the agonist-bound receptor concentration that can couple to the signaling machinery. This is governed by the desensitization rate constant, \({k}_{DES}\). The basic model is formulated in Appendix 2.1 and is represented by Scheme 1 below:

The time course shape predicted by this model is an association exponential curve (Fig. 9b). This is consistent with experimental results; systems in which the sole regulation mechanism is receptor desensitization yield an association exponential time course. Examples include IP production by the AT_{1} (Fig. 3a)^{40}, PAR1 (see Fig. 2b in^{41}) and C-terminally-extended GnRH receptor (see Fig. 3c in^{43}) (when response degradation is blocked by Li^{+}). This curve shape makes sense intuitively. At early times response is generated rapidly. The response then slows, because response generation is attenuated by the loss of active receptor. Ultimately the response approaches a limit (the plateau). At this limit the response level does not change because no new response is being generated and the existing response is not degraded.

This model can now be extended to incorporate both regulation mechanisms, receptor desensitization and response degradation (Appendix 2.2), represented by Scheme 2 below:

The time course shape for this model is a rise-and-fall to baseline curve (Fig. 9e). This makes sense intuitively. The response rises rapidly then slows as receptor desensitization starts to slow the rate of response generation. The response becomes further limited owing to response degradation. The response reaches an upper limit (the peak) when the rate of response generation equals the rate of degradation. After this, response degradation predominates over response generation. Less and less new response is generated because the active receptor concentration is declining, and ultimately no new response will be generated because the active receptor concentration will decline to zero. Ultimately the response level falls to zero once all existing response has been degraded. This profile is evident in systems where both mechanisms have been shown to be in operation. Examples include diacylglycerol production via the AT_{1} receptor (Fig. 4a)^{65}, cAMP accumulation via the β_{2} adrenoceptor in the absence of phosphodiesterase inhibition (Fig. 4b)^{53}, and IP production via the CCK1 receptor in the absence of Li^{+} (Fig. 10a of^{90}).

#### More complex models

More complex regulation and signaling mechanisms discovered experimentally can be represented using the model formulation, as described in Appendix 2.3. Receptor can resensitize and then contribute to signaling. This model is formulated in combination with response degradation in Appendix 2.3.1. The receptor, once desensitized and internalized, can continue to signal from internal compartments^{16,17,18}. This model is formulated in Appendix 2.3.2. A third complex model describes calcium signaling; response is degraded (cleared from the cytoplasm) but the response can reform (from calcium that re-enters the cell). This process, in combination with depletion of the response precursor, is formulated in Appendix 2.3.3. In all three models the time course is described by a rise-and-fall to steady-state equation (see Appendix 2.3).

#### Comparison of models with experimental data

The receptor desensitization model completes the set of kinetic mechanistic signaling models developed for GPCRs, initiated in ref.^{33}. This allows comparison of a complete kinetic model of GPCR signaling and regulation with experimental data, and with the model-free analysis. The model mechanisms are shown schematically in Fig. 9 and are: no regulation; receptor desensitization; response degradation; response recycling; receptor desensitization and response degradation; precursor depletion and response degradation; and the three more complicated models described above and in Appendix 2.3. Simulators of the models are provided in Supporting Material, enabling investigators to evaluate changes of the parameter values.

The time course curve shapes predicted by the model are those observed experimentally – the straight line, the association exponential, and the two rise-and-fall curves (Fig. 9). It emerges that the shape is dependent on the number of regulatory mechanisms. With no regulatory mechanisms, the time course is a straight line (Fig. 9a). With one mechanism, an association exponential curve results (receptor desensitization, response degradation, or response recycling, Fig. 9b–d). With two mechanisms, a rise-and-fall to baseline curve results (Fig. 9e,f). More precisely, when an input regulation mechanism (desensitization or precursor depletion) is coupled with an output mechanism (response degradation). With the more complicated mechanisms, a rise-and-fall to baseline profile results (Fig. 9g–i).

#### Comparison of models with model-free analysis

The four curve shapes of the mechanistic model are the same as those of the model-free analysis. This is shown graphically by comparing Fig. 9 with Fig. 1. It is also evident mathematically by inspection of the equations; the model equations listed in Supplementary Tables S1–S4 conform to the forms of the general equations (Eqs. 1–4). This finding indicates the model-free analysis emerges from known mechanisms of GPCR signaling and regulation.

The maximal initial rate from the model-free analysis (\({\mathrm{IR}}_{\mathrm{max}}\)) is equal to *k*_{τ} from the mechanistic model. This is evident by taking the initial rate of the equations, i.e. the limit as time approaches zero. By definition this limit is \({\mathrm{IR}}_{\mathrm{max}}\) for the model-free analysis. For the mechanistic equations, this limit is *k*_{τ} (for a maximally-stimulating concentration of agonist, Appendix 1.2). This supports the hypothesis that the initial rate is purely an efficacy term, being the rate of response generation before it is impacted by regulation of signaling mechanisms, because *k*_{τ} contains efficacy terms but no regulation terms. It also provides a mechanistic foundation of \({\mathrm{IR}}_{\mathrm{max}}\), indicating it is analogous to the initial of an enzyme’s activity (formally, the product of the receptor and precursor concentrations and a rate constant). Finally, this formulation demonstrates the maximal initial rate is system-dependent, being dependent on the receptor and precursor levels, which can differ between systems.

#### Measuring model parameters

The models contain parameters for signal generation and for regulation of signaling. These parameters can be estimated by curve fitting by two methods. The time course data can be fit directly to equations that explicitly describe the model. These equations are listed in Supplementary Information and are provided in a custom Prism template in the supplementary files (“Signaling kinetic mechanism equations”). Alternatively, the data can be fit to general time course equations (Eqs. 1–4) and the fitted parameter values used to calculate the model parameters. These approaches are used here to estimate efficacy and regulation parameters, which can be done using just a maximally-stimulating concentration of agonist. Note in all the literature examples in this study, a maximally-stimulating concentration was used.

The general equation fitting approach is used here for the linear, association exponential, and rise-and-fall to baseline examples. For the linear examples, the only parameter to be estimated is *k*_{τ} and this is equal to the slope of the line for a maximally-stimulating concentration of agonist. The resulting *k*_{τ} values are the Slope values in Fig. 2 and Supplementary Fig. 1.

For the association exponential time course, there are two parameters, the efficacy parameter *k*_{τ}, and a regulation parameter specified by the model (for example, the desensitization rate constant \({k}_{DES}\)). *k*_{τ} can be estimated using the same method as that used in the model-free analysis to calculate \({\mathrm{IR}}_{\mathrm{max}}\), since these parameters are equivalent (see “Measuring the initial rate from curve fit parameters” above). The rate constant \(k\) is multiplied by the steady-state level of response, for a maximally-stimulating concentration of agonist. This resulting *k*_{τ} values are equal to the initial rate values given in Fig. 3 and Supplementary Fig. 2, and the \({\mathrm{IR}}_{\mathrm{max}}\) values in Table 1. Note in all the model variants that give rise to the association exponential profile, *k*_{τ} is the product of \(k\) and the steady-state response, i.e. *k*_{τ} is calculated in the same way. This means *k*_{τ} can be estimated without knowing the specific mechanism underlying the curve shape. The regulation parameter of the association exponential time course is dependent on the mechanism. In all cases, the regulation parameter is given by the rate constant parameter \(k\) of the general equation (Supplementary Table S2). If the mechanism is known, the value of \(k\) can be ascribed to a regulatory process. This is shown using the studies cited here where the regulation of signaling mechanism was evaluated. In Fig. 3a, IP generation by the AT_{1} receptor^{40}, the mechanism was most likely receptor desensitization since deletion of arrestin changes the shape to a straight line (Fig. 2a) and response degradation was blocked by Li^{+}. The rate constant \(k\) then represents the desensitization rate constant \({k}_{DES}\). The fitted value of 0.040 min^{−1} indicated a receptor desensitization half time of 17 min. In Fig. 3b,c, the mechanism is most likely response degradation since receptor desensitization was minimized, the GnRH receptor lacking a C-terminal tail^{45} (Fig. 3b), the β_{2} adrenoceptor-expressing cells devoid of arrestin^{40} (Fig. 3c), and in both cases inhibitors of degradation excluded. The \(k\) value then represents the response degradation rate constant \({k}_{D}\). The fitted values, 0.12 min^{−1} for the GnRH_{1} receptor and 1.1 min^{−1} for the β_{2} adrenoceptor, indicated response degradation half-times of 5.8 min and 38 s, respectively.

For the rise-and-fall to baseline mechanisms, there are three parameters – the efficacy parameter *k*_{τ}, and two regulation parameters. *k*_{τ} can again be estimated using the general equation without knowing the specific mechanism underlying the time course profile. *k*_{τ} is equal to the parameter \(C\) of the general equation (Eq. 3) when a maximally-stimulating agonist concentration is employed, and these values are shown for the literature examples in Fig. 4. The rate constant values can be ascribed to the mechanism rate constants if the mechanism is known. In Fig. 4a, diacyglycerol production by the AT_{1} receptor, the mechanism is most likely receptor desensitization and response degradation, since blocking arrestin recruitment changes the shape to an association exponential (see Fig. 2a of ref.^{65}) and diacyclglycerol is rapidly degraded by diacylglycerol kinases^{38}. It is then necessary to determine which of the general rate constants (\({k}_{1}\) and \({k}_{2}\)) corresponds to which of the mechanism rate constants (\({k}_{DES}\) and \({k}_{D}\)). In this case, it is likely the faster of the rates (\({k}_{1}\)) corresponds to \({k}_{D}\), since blocking arrestin recruitment gives an association exponential profile with rate similar to \({k}_{1}\) (see Fig. 2a of ref.^{65}). This approach gives a \({k}_{D}\) (\({k}_{1}\)) value of 5.5 min^{−1} and a \({k}_{DES}\) (\({k}_{2}\)) value of 0.80 min^{−1}, corresponding to degradation and desensitization half times of 7.6 s and 52 s, respectively. In Fig. 4c, cytoplasmic Ca^{2+} elevation via the GnRH receptor, the mechanism is most likely precursor depletion and response degradation, based on the known mechanism of calcium responses (see “Rise-and-fall to baseline time course profile” above). The more rapid of the two rate constants is most likely precursor depletion since this is the first event in the signaling cascade. This logic gives a \({k}_{DEP}\) (\({k}_{1}\)) value of 0.46 s^{−1} and a \({k}_{D}\) (\({k}_{2}\)) value of 0.072 s^{−1}, corresponding to depletion and degradation half times of 1.5 s and 9.6 s, respectively. Figure 4b shows cAMP generation by the β_{2} adrenoceptor. In this example, likely resulting from receptor desensitization and response degradation based on the experiments in ref.^{53}, the two rate constant values are close to one another (1.5 and 0.96 min^{−1}), preventing the ascribing of the constants to the regulatory processes but demonstrating the rate of the processes is similar in this system.

For the rise-and-fall to steady-state models, *k*_{τ} can be estimated using the general equation and a maximally-stimulating concentration of agonist. *k*_{τ} is equal to \({\mathrm{IR}}_{\mathrm{max}}\), calculated from the parameters as described above. For estimating the regulation parameters, we recommend using the mechanistic model equation rather than the general equation. For this curve shape, the general equation parameters correspond to complicated combinations of the mechanism parameters, which introduces constraints of the parameter values that are not incorporated into the general model equations (Appendix 2.3). This approach was used to analyze the Ca^{2+} mobilization data in Fig. 5a. In this experiment, Ca^{2+} was included in the extracellular medium, which results in a steady-state being reached between cytoplasmic Ca^{2+} entry and export. This is the “Response degradation to steady-state with precursor depletion” model (Appendix 2.3.3). The equation (Eq. (17), see Methods) fitted the data well (R^{2} value 0.995, see curve in Fig. 5a, standard error of the fitted model parameters less than 10%—see “Curve fit results” Excel file in Supporting Material). The regulation rate constant values were 0.37 s^{−1} for \({k}_{DEP}\), 0.052 s^{−1} for \({k}_{D}\) and 0.0033 s^{−1} for \({k}_{R}\), corresponding to half times for precursor depletion, response degradation and response reformation of 1.9 s, 13 s and 3.5 min, respectively. As expected, the *k*_{τ} value from the mechanism equation fit was the same as the value from the model-free analysis (240 nM s^{−1}).

## Discussion

New technologies have enabled efficient measurement of the kinetics of GPCR signaling using continuous-read, real-time bioassay modalities^{12,13,14}. Methods are now required to extract pharmacological parameters from the time course data. Such methods have been developed for specific responses that are immediately proximal to the receptor (internalization^{25} and arrestin-receptor interaction^{26}). Here methods were designed for universal application to GPCR signaling and regulation responses. Curve-fitting methods are described for analyzing time course signaling data, designed for routine use by pharmacologists, enabling application to medicinal chemistry and receptor research. These models could account for a broad array of historical time course signaling data. The equations are simple and built into commonly-used curve-fitting programs, or are provided in ready-to-use templates. The efficacy for signaling is quantified kinetically by an intuitive, familiar and biologically meaningful parameter, the initial rate of signaling, which is obtained from the curve fit parameters. The underlying theoretical framework is familiar, based on the concepts of enzyme kinetics and the known mechanisms of regulation of signaling. Resources are provided in Supporting Material to facilitate implementation and application of the models, including Prism templates containing the equations, and a time course data simulator for the mechanistic models. With these new analysis tools, the model can be tested prospectively using new data, testing the effects of regulation of signaling, evaluating how the initial rate characterization compares with historical single time point pharmacological methods, and investigating the effect of receptor reserve (see below).

This approach has practical benefits for optimization of new molecules and for GPCR signaling research. In project workflow, it allows representation of the time course data set by a minimal number of informative parameters (e.g. *k*_{τ}, \({k}_{DES}\)). This enables tabulation of kinetic signaling data for a series of ligands, necessary for medicinal chemists to establish kinetic structure–activity relationships. If the response at a specific time point is required, for example when comparing in vivo and in vitro activity, this can be calculated from the curve fit parameters, avoiding the necessity of measuring precisely the same time points in every experiment. The approach solves the time-dependency problem, in which pharmacological parameter values can be dependent on the time point at which they are measured, complicating quantification of efficacy and biased agonism^{24,25,26}. This is because: (1) Activity is quantified as a rate constant value, which is constant over time; (2) Differences of curve shape can be accommodated because the same efficacy parameter, the initial rate, can be extracted from all four of the commonly-encountered time course curve shapes. Finally, the method separates efficacy and regulation, allowing these parameters to be quantified separately (e.g. *k*_{τ} and \({k}_{DES}\)), enabling independent structure–activity optimization of these processes.

The mechanistic model yields insights that explain the time course curve shapes and provide mechanistic meaning to the parameters from the model-free analysis. The model and experimental data show that regulation of signaling defines the time course curve shape. Regulation changes the shape from a straight line (where response is generated indefinitely), to an association exponential curve or a rise-and-fall curve, with the type of curve dependent on the number and nature of the regulation mechanisms (Fig. 9). The underlying mechanism can be interrogated by blocking receptor desensitization and/or response degradation. The model supports the hypothesis that the initial rate of activity is purely an efficacy term, unaffected by regulation of signaling: It is shown that *k*_{τ}, devoid of regulation terms, is the limit of the equations as time approaches zero, the formal definition of the initial rate (Appendix 1). The model translates the empirical parameters of the model-free analysis into biologically-meaningful parameters of GPCR signaling mechanisms. \({\mathrm{IR}}_{\mathrm{max}}\) is equivalent to *k*_{τ}, the initial rate of signaling generation by the agonist occupied receptor, analogous to the initial rate of an enzyme’s activity. The rate constants \(k\), \({k}_{1}\) and \({k}_{2}\) correspond simply to the rate constants of the regulation processes, e.g. receptor desensitization and response degradation. If the mechanism is known, the rate constant value can be assigned to a regulatory activity, enabling regulation of signaling to be quantified in simple, kinetic terms.

The relationship between the kinetic analysis and existing pharmacological analysis is now considered. Established pharmacological analysis quantifies signaling in terms of efficacy and affinity, usually at a single time point at which equilibrium between receptor and ligand is assumed^{3,4,5,6,7,8}. The efficacy term in the kinetic approach (\({\mathrm{IR}}_{\mathrm{max}}\) or *k*_{τ}) is a direct analogue of established efficacy parameters (E_{max}, and τ of the operational model^{7}). Consequently, the rank order of efficacy should be the same using both approaches. Affinity is more difficult to measure using the kinetic approach because of the equilibration issue. Properly measuring affinity from a concentration–response requires the assay to be incubated long enough for equilibrium between receptor and ligand to be closely approached^{91}. By contrast, the initial rate in the kinetic method is quantified using the earliest time points after the addition of ligand and under these conditions it cannot be generally assumed the ligand and receptor are at equilibrium^{5,33,92,93}. The impact can be evaluated by considering typical values of receptor-ligand binding rate constants, using the mass-action receptor-ligand association kinetics equation^{94}. For low affinity ligands, the equilibration time is short, potentially enabling accurate estimation of affinity when using responses that proceed over several minutes. A compound with 10 μM affinity, binding with an association rate constant (\({k}_{\mathrm{on}}\)) of 10^{8} M^{−1} min^{−1}, reaches 97% of its equilibrium occupancy within 0.1 s at its \({K}_{\mathrm{d}}\) concentration (i.e. 10 μM). However, for high affinity ligands the equilibration time can be long relative to the duration of the response. A compound with 1 nM affinity and a \({k}_{\mathrm{on}}\) of 10^{8} M^{−1} min^{−1} takes 18 min to reach 97% of its equilibrium receptor occupancy when applied at its \({K}_{\mathrm{d}}\) concentration. The problem is magnified when the response is very rapid, as previously described^{92}, particularly in Ca^{2+} signaling, which proceeds in seconds. This issue can be circumvented when receptor and ligand can be pre-incubated prior to the initiation of the response, for example in GTPγS binding assays, but this option is rarely available in whole cell assays, the most commonly used assay modality. Importantly, this issue does not affect estimation of efficacy. \({\mathrm{IR}}_{\mathrm{max}}\) and *k*_{τ} are estimated using maximally-effective concentrations which, due to mass action, bind more rapidly. For example, when applied at 10 μM concentration, the 1 nM \({K}_{\mathrm{d}}\), 10^{8} M^{−1} min^{−1} \({k}_{\mathrm{on}}\) compound reaches 97% of equilibrium occupancy within 0.2 s. Finally, the effect of signal amplification, manifest as receptor reserve, in the kinetic mechanistic model requires further investigation. While the model can incorporate receptor reserve intrinsically (as a depletion of response precursor)^{33}, the rate of signaling in a multistep pathway is ultimately a function of the rates of the individual steps and the stoichiometric relationships between them. For example, the initial rate could be susceptible to receptor reserve effects if there is a rate-limiting step upstream of the response being measured. This issue requires experimental investigation.

In conclusion, this study introduces straightforward data analysis methods to quantify the kinetics of GPCR signaling. The simplicity of the analysis procedures, the intuitive nature of the parameters, and the mechanistic foundation in known and emerging GPCR signaling paradigms will facilitate pharmacological discovery and optimization in kinetic terms.

## Methods

### Extraction of time course data

Data were extracted from figures in published articles using a plot digitizer (Graph Grabber v2, Quintessa Limited, Henley-on-Thames, UK).

### Equations for model-free analysis

For model-free analysis, time course data were fit to empirical time course equations. The following equations were used when the response started at the moment of ligand addition. The straight line equation is Eq. (1),

where *y* is response, Slope is the gradient of the line and *t* is the time of response measurement. Baseline is the response in the absence of ligand. Note Baseline is assumed to remain constant over time. The corresponding equation in Prism is named “Straight line” and the parameter Yintercept in this equation corresponds to Baseline^{50}.

The association exponential equation is Eq. (2):

where SSR is the steady-state response, specifically the ligand-specific response as time approaches infinity. Note the *y* asymptote value as time approaches infinity is SSR + Baseline. \(k\) is the observed rate constant in units of *t*^{−1}. The corresponding equation in Prism 8 is named “One-phase association,”^{61} in which Span corresponds to SSR, Y0 corresponds to Baseline, and Plateau corresponds to SSR + Baseline.

The rise-and-fall to baseline equation is Eq. (3):

where \(C\) is a fitting constant in units of y-units.*t*^{−1} (which is the initial rate of signaling, see Appendix 1.1) and \({k}_{1}\) and \({k}_{2}\) are observed rate constants in units of *t*^{−1}. In the analysis, \({k}_{1}\) is the larger of the two rate constant values, constrained to be greater than \({k}_{2}\) (see “Fitting procedures” below). This equation has been loaded into a Prism template as a user-defined equation named “Rise-and-fall to baseline time course” (see “Signaling kinetic model-free equations” file available in the supplementary files).

The rise-and-fall to steady-state equation is Eq. (4):

where, as above, SSR is the steady-state response, specifically the ligand-specific response as time approaches infinity. Note the *y* asymptote value as time approaches infinity is SSR + Baseline. \(D\) is a unitless fitting constant, and \({k}_{1}\) and \({k}_{2}\) are observed rate constants in units of *t*^{−1}. In the analysis, \({k}_{1}\) is the larger of the two rate constant values, constrained to be greater than \({k}_{2}\) (see “Fitting procedures” below). This equation has been loaded into a Prism template as a user-defined equation named “Rise-and-fall to steady state time course” in which the parameter SteadyState corresponds to SSR (see “Signaling kinetic model-free equations” file available in the supplementary files).

In certain circumstances it is desirable to float the initiation of signaling time in the analysis. This enables the analysis to accommodate slight uncertainty as to the precise time point of ligand addition, especially for rapidly-generated signals. It also allows for signaling mechanisms where there is a delay between ligand-receptor binding and initiation of detectable signaling, often observed in gene expression assays, and thought to be due to a necessary build-up of signal transduction intermediates to a threshold level^{48}. When floating the start time in the analysis it is necessary to establish the baseline response level prior to the start time. Equations 1–4 above can be adapted to incorporate a fitted start time, the corresponding equations being Eqs. (5–8) below:

where \({t}_{0}\) is the signal initiation time and *t* is time uncorrected for the signal initiation time. Equation (6) is equivalent to the equation named “Plateau followed by one phase association” in Prism 8, in which Span corresponds to SSR and Y0 corresponds to Baseline^{95}. Equations 5,7 and 8 have been loaded into a Prism template as User-defined equations named “Baseline then straight line time course”, “Baseline then rise-and-fall to baseline time course” and “Baseline then rise-and-fall to steady state time course” respectively (see “Signaling kinetic model-free equations” Prism file available in the supplementary files). In Prism the equations are coded with an IF statement to fit the baseline before \({t}_{0}\). For example, Eq. (7) is written in Prism as,

In some cases, stimulation of signaling results in a decrease of response signal. In other words, the detected response signal decreases over time. This can result from technical considerations, for example, certain biosensors decrease in signal upon binding the response analyte (so called “Downward” sensors). It can also result from the signaling mechanism, for example stimulation of G_{i} is often detected by the resulting inhibition of cAMP accumulation. This can be handled by taking the inverse of the detected signal and plotting this versus time (creating an “Upward” time course, see for example Supplementary Fig. S2d). Alternatively, downward data can be fit to the downward analogue of the equations. Examples of this approach are used in this study (inhibition of cAMP accumulation, see Supplementary Fig. S2f, and GIRK channel currents, see Supplementary Fig. S2h); data were analyzed with the downward analogue of the association exponential equation, which is the exponential decay equation shown below (Eq. 9):

Note the response at \(t=0\) is equal to SSR + Baseline. The corresponding equation in Prism 8 is named “Plateau followed by one phase decay,” in which Span corresponds to SSR, Plateau corresponds to Baseline, and Y0 corresponds to SSR + Baseline^{96}.

### Calculation of initial rate

The initial rate of signaling was determined by first fitting the data to the time course equations, then using the fitted parameter values to calculate the initial rate (\(\mathrm{IR}\)). This calculation utilized an equation defining the initial rate in terms of the fitted parameters, the equation obtained by taking the limit of the time course equation as time approaches zero (the formal definition of the initial rate condition). These equations were derived in Appendix 1:

Straight line:

Association exponential:

Rise-and-fall to baseline:

Rise-and-fall to steady-state:

### Equations for kinetic mechanistic model analysis and simulation

The kinetic mechanistic model was used to simulate and analyze time course data using the model equations, derived in Ref.^{33} and Appendix 2. The equations employed are given in Supplementary Tables S1–S4 and are listed in Supplementary information. A Prism template containing the equations called “Signaling kinetic mechanism equations” is available in the supplementary files. The terms in the models are defined in Supplementary Table S5. Data in Fig. 5a (Ca^{2+} mobilization via the GnRH_{1} receptor) were analyzed using the equation for the precursor depletion & response degradation to steady-state model (Eq. (16), Appendix 2.3.3) that incorporated a variable signal initiation time (\({t}_{0}\)):

where,

### Fitting procedures

Data were fit by nonlinear regression to the equations using Prism 8. In almost all cases, the default fitting settings were used, specifically least-squares nonlinear regression with medium convergence criteria, no special handling of outliers and no weighting^{78}. In two cases, Figs. 4b and 5b, the default settings yielded an ambiguous fit, in that some of the parameter values were not resolved (see “Curve fit results” Excel file in Supporting Material)^{79}. In these cases, unambiguous parameter values were obtained using the “Robust regression” option which minimizes the contribution of outliers to the fit^{74}. In the rise-and-fall fits (Eqs. 3, 4, 7, 8), which contain two rate constants \({k}_{1}\) and \({k}_{2}\), \({k}_{1}\) was assumed to be the larger of the rate constant values and this was handled by constraining \({k}_{1}\) to be greater than \({k}_{2}\). In all cases, rate constant values were constrained to be greater than zero. The fitted values, the standard error and the correlation coefficient R^{2} are given in the “Curve fit results” Excel file in Supporting Material. Values were given to four significant figures. When the signal initiation time was floated (Eqs. 5–9), the initial value (X0) was entered manually.

In order to obtain an estimate of the error of the fitted parameters, the standard error was computed as described^{51}, assuming a symmetrical confidence interval. This required changing the default “Calculate CI or parameters” from asymmetrical to symmetrical in the “Confidence” dialogue^{97}.

The baseline response, that in the absence of ligand, was handled by incorporating a baseline term into the equations and fitting this parameter as part of the analysis. A commonly-used alternative, subtracting baseline from the data set and excluding the baseline term from the equation, was found to distort the fit significantly for the rise-and-fall to baseline equation, unless the subtracted baseline value was in very close agreement with that determined by incorporating the baseline into the curve-fitting procedure.

For the V_{2} vasopressin receptor experiments, the biosensor fluorescence intensity data were normalized to baseline. Specifically, response was quantified as the fluorescence after ligand addition divided by the mean baseline signal before addition, a metric termed F/ΔF. (The mean baseline signal is the mean of the fluorescence intensity measurements before ligand addition.) The data was handled in this way because it provides a convenient intra-well control, minimizing error due to any slight well-to-well differences in the amount of sensor or number of cells. The arrestin sensor data was inverted by calculating the value 1 − F/ΔF. This was done because the arrestin sensor is a downward sensor (a decrease in fluorescence resulted from receptor interaction)^{26}. In the fitting procedure, each technical replicate was considered an individual point.

Concentration–response data for the initial rate were fit to a sigmoid curve equation, the “Log(agonist) vs. response − Variable slope” equation in Prism^{83}:

The “Bottom” parameter was constrained to zero. (Note in the Prism formulation, L_{50} is written as EC_{50}. The EC_{50} term is not used here because it has an explicit pharmacological definition^{62}.)

### Measurement of cAMP generation and arrestin recruitment using biosensors

Genetically-encoded biosensors were used to measure cAMP generation and arrestin recruitment via the V_{2} vasopressin receptor. The sensors have been described previously (red cADDis for cAMP^{82} and see ref.^{26}. for arrestin recruitment), and are packaged in the BacMam vector for delivery to cells. The cDNA for the V_{2} vasopressin receptor was obtained from the cDNA Resource Center (Bloomsburg University, Bloomsburg, PA). The experiments were conducted in HEK293T cells transfected with the receptor and each sensor individually (i.e. one batch of cells with receptor and the cAMP sensor and a separate batch with receptor and arrestin sensor). HEK 293T cells were cultured in Eagle’s minimum essential media (EMEM) supplemented with 10% fetal bovine serum (FBS) and penicillin–streptomycin at 37 °C in 5% CO_{2}. One day before the transfection, HEK293T cells were seeded at a density of 27,000 cells/well in a 96-well plate (Greiner Cellcoat #655946). Approximately 16–20 h later, cells were transfected with 50 ng of the V_{2} receptor plasmid per well, using Lipofectamine 2000 from Invitrogen (Waltham, MA USA). After an approximately 4 h incubation period at 37 °C and 5% CO_{2}, the transfection mix was replaced with 100uL complete cell culture media and the cells were allowed to recover for approximately 1.5 h under normal growth conditions. The cells were then transduced with either the red cADDis or green arrestin sensor BacMam stocks. To prepare the transduction mixture, the BacMam containing cADDis or the arrestin sensor, 2 mM sodium butyrate, and EMEM were combined in a final volume of 50 μL. For each experiment, 6.18 × 10^{8} viral genes of cADDis virus or 4.24 × 10^{8} viral genes of arrestin sensor virus were added to each well. The 50 uL transduction mixture was then added to the 96-well plate (50 uL/well) and incubated for approximately 24 h at 37 °C in 5% CO_{2}.

Thirty minutes prior to fluorescence plate reader experiments, the media in each well was replaced with 150 μL of Dulbecco’s phosphate buffered saline (DPBS) supplemented with Ca^{2+} (0.9 mM) and Mg^{2+} (0.5 mM). Fluorescence plate reader experiments were performed on the Synergy Mx reader (BioTek, Winooski, VT). The green fluorescence detection was recorded using 488/20 nm excitation and 525/20 nm fluorescence emission, while red fluorescence detection was recorded using 565/20 nm excitation wavelength and 603/20 nm fluorescence emission. After acquiring the baseline fluorescence for several minutes, drug was added manually with a multichannel pipette in a volume of 50 µL at the indicated time points.

Oxytocin and vasopressin were obtained from Cayman Chemical (Ann Arbor, MI USA). All working concentration of drugs were dissolved in DPBS and added manually to the HEK 293T cells at the indicated concentrations and time points.

### Statistical analysis

Concentration–response data for V_{2} vasopressin receptor responses to vasopressin and oxytocin were compared statistically by paired, two-tailed t-test using Prism 8. Comparisons were performed for maximal effect and L_{50} (Tables 1, 2).

## References

- 1.
Clark, A. J. The reaction between acetyl choline and muscle cells.

*J. Physiol.***61**, 530–546. https://doi.org/10.1113/jphysiol.1926.sp002314 (1926). - 2.
Gaddum, J. H. The action of adrenalin and ergotamine on the uterus of the rabbit.

*J. Physiol.***61**, 141–150. https://doi.org/10.1113/jphysiol.1926.sp002280 (1926). - 3.
Rang, H. P. The receptor concept: pharmacology’s big idea.

*Br. J. Pharmacol.***147**(Suppl 1), S9–S16. https://doi.org/10.1038/sj.bjp.0706457 (2006). - 4.
Kenakin, T. Quantifying biological activity in chemical terms: a pharmacology primer to describe drug effect.

*ACS Chem. Biol.***4**, 249–260. https://doi.org/10.1021/cb800299s (2009). - 5.
Slack, R. J. & Hall, D. A. Development of operational models of receptor activation including constitutive receptor activity and their use to determine the efficacy of the chemokine CCL17 at the CC chemokine receptor CCR4.

*Br. J. Pharmacol.***166**, 1774–1792. https://doi.org/10.1111/j.1476-5381.2012.01901.x (2012). - 6.
Finlay, D. B., Duffull, S. B. & Glass, M. 100 years of modelling ligand-receptor binding and response: a focus on GPCRs.

*Br. J. Pharmacol.*https://doi.org/10.1111/bph.14988 (2020). - 7.
Black, J. W. & Leff, P. Operational models of pharmacological agonism.

*Proc. R. Soc. Lond. B Biol. Sci.***220**, 141–162 (1983). - 8.
Zhao, P. & Furness, S. G. B. The nature of efficacy at G protein-coupled receptors.

*Biochem. Pharmacol.***170**, 113647. https://doi.org/10.1016/j.bcp.2019.113647 (2019). - 9.
Kenakin, T. Biased receptor signaling in drug discovery.

*Pharmacol. Rev.***71**, 267–315. https://doi.org/10.1124/pr.118.016790 (2019). - 10.
Rodbell, M., Birnbaumer, L., Pohl, S. L. & Krans, H. M. The glucagon-sensitive adenyl cyclase system in plasma membranes of rat liver. V. An obligatory role of guanylnucleotides in glucagon action.

*J. Biol. Chem.***246**, 1877–1882 (1971). - 11.
Tsien, R. Y., Pozzan, T. & Rink, T. J. Calcium homeostasis in intact lymphocytes: cytoplasmic free calcium monitored with a new, intracellularly trapped fluorescent indicator.

*J. Cell. Biol.***94**, 325–334. https://doi.org/10.1083/jcb.94.2.325 (1982). - 12.
Lohse, M. J., Nuber, S. & Hoffmann, C. Fluorescence/bioluminescence resonance energy transfer techniques to study G-protein-coupled receptor activation and signaling.

*Pharmacol. Rev.***64**, 299–336. https://doi.org/10.1124/pr.110.004309 (2012). - 13.
Marullo, S. & Bouvier, M. Resonance energy transfer approaches in molecular pharmacology and beyond.

*Trends. Pharmacol. Sci.***28**, 362–365. https://doi.org/10.1016/j.tips.2007.06.007 (2007). - 14.
Tewson, P.

*et al.*Simultaneous detection of Ca^{2+}and diacylglycerol signaling in living cells.*PLoS ONE***7**, e42791. https://doi.org/10.1371/journal.pone.0042791 (2012). - 15.
Shear, M., Insel, P. A., Melmon, K. L. & Coffino, P. Agonist-specific refractoriness induced by isoproterenol. Studies with mutant cells.

*J. Biol. Chem.***251**, 7572–7576 (1976). - 16.
Ferrandon, S.

*et al.*Sustained cyclic AMP production by parathyroid hormone receptor endocytosis.*Nat. Chem. Biol.***5**, 734–742. https://doi.org/10.1038/nchembio.206 (2009). - 17.
Mullershausen, F.

*et al.*Persistent signaling induced by FTY720-phosphate is mediated by internalized S1P1 receptors.*Nat. Chem. Biol.***5**, 428–434. https://doi.org/10.1038/nchembio.173 (2009). - 18.
Calebiro, D.

*et al.*Persistent cAMP-signals triggered by internalized G-protein-coupled receptors.*PLoS Biol.***7**, e1000172. https://doi.org/10.1371/journal.pbio.1000172 (2009). - 19.
Jensen, D. D.

*et al.*Neurokinin 1 receptor signaling in endosomes mediates sustained nociception and is a viable therapeutic target for prolonged pain relief.*Sci. Transl. Med.*https://doi.org/10.1126/scitranslmed.aal3447 (2017). - 20.
Eichel, K. & von Zastrow, M. Subcellular organization of GPCR signaling.

*Trends Pharmacol. Sci.***39**, 200–208. https://doi.org/10.1016/j.tips.2017.11.009 (2018). - 21.
Hothersall, J. D., Brown, A. J., Dale, I. & Rawlins, P. Can residence time offer a useful strategy to target agonist drugs for sustained GPCR responses?.

*Drug Discov. Today***21**, 90–96. https://doi.org/10.1016/j.drudis.2015.07.015 (2016). - 22.
Lobingier, B. T. & von Zastrow, M. When trafficking and signaling mix: how subcellular location shapes G protein-coupled receptor activation of heterotrimeric G proteins.

*Traffic***20**, 130–136. https://doi.org/10.1111/tra.12634 (2019). - 23.
Stoeber, M.

*et al.*A genetically encoded biosensor reveals location bias of opioid drug action.*Neuron***98**, 963–976. https://doi.org/10.1016/j.neuron.2018.04.021 (2018). - 24.
Klein Herenbrink, C.

*et al.*The role of kinetic context in apparent biased agonism at GPCRs.*Nat. Commun.***7**, 10842. https://doi.org/10.1038/ncomms10842 (2016). - 25.
Zhu, X., Finlay, D. B., Glass, M. & Duffull, S. B. Model-free and kinetic modelling approaches for characterising non-equilibrium pharmacological pathway activity: internalisation of cannabinoid CB1 receptors.

*Br. J. Pharmacol.***176**, 2593–2607. https://doi.org/10.1111/bph.14684 (2019). - 26.
Hoare, S. R. J., Tewson, P. H., Quinn, A. M. & Hughes, T. E. A kinetic method for measuring agonist efficacy and ligand bias using high resolution biosensors and a kinetic data analysis framework.

*Sci. Rep.***10**, 1766. https://doi.org/10.1038/s41598-020-58421-9 (2020). - 27.
Zhao, P.

*et al.*Activation of the GLP-1 receptor by a non-peptidic agonist.*Nature***577**, 432–436. https://doi.org/10.1038/s41586-019-1902-z (2020). - 28.
Willars, G. B., McArdle, C. A. & Nahorski, S. R. Acute desensitization of phospholipase C-coupled muscarinic M3 receptors but not gonadotropin-releasing hormone receptors co-expressed in alphaT3-1 cells: implications for mechanisms of rapid desensitization.

*Biochem. J.***333**(Pt 2), 301–308. https://doi.org/10.1042/bj3330301 (1998). - 29.
Wessling-Resnick, M. & Johnson, G. L. Allosteric behavior in transducin activation mediated by rhodopsin. Initial rate analysis of guanine nucleotide exchange.

*J. Biol. Chem.***262**, 3697–3705 (1987). - 30.
Rubenstein, R. C., Linder, M. E. & Ross, E. M. Selectivity of the beta-adrenergic receptor among Gs, Gi’s, and Go: assay using recombinant alpha subunits in reconstituted phospholipid vesicles.

*Biochemistry***30**, 10769–10777. https://doi.org/10.1021/bi00108a023 (1991). - 31.
Gimenez, L. E., Baameur, F., Vayttaden, S. J. & Clark, R. B. Salmeterol efficacy and bias in the activation and kinase-mediated desensitization of beta2-adrenergic receptors.

*Mol. Pharmacol.***87**, 954–964. https://doi.org/10.1124/mol.114.096800 (2015). - 32.
Pedrosa, R., Gomes, P. & Soares-da-Silva, P. Distinct signalling cascades downstream to Gsalpha coupled dopamine D1-like NHE3 inhibition in rat and opossum renal epithelial cells.

*Cell Physiol. Biochem.***14**, 91–100. https://doi.org/10.1159/000076930 (2004). - 33.
Hoare, S. R. J., Pierre, N., Moya, A. G. & Larson, B. Kinetic operational models of agonism for G-protein-coupled receptors.

*J. Theor. Biol.***446**, 168–204. https://doi.org/10.1016/j.jtbi.2018.02.014 (2018). - 34.
Moore, C. A., Milano, S. K. & Benovic, J. L. Regulation of receptor trafficking by GRKs and arrestins.

*Annu Rev Physiol***69**, 451–482. https://doi.org/10.1146/annurev.physiol.69.022405.154712 (2007). - 35.
Hausdorff, W. P., Caron, M. G. & Lefkowitz, R. J. Turning off the signal: desensitization of beta-adrenergic receptor function.

*FASEB J.***4**, 2881–2889 (1990). - 36.
Krupnick, J. G. & Benovic, J. L. The role of receptor kinases and arrestins in G protein-coupled receptor regulation.

*Annu. Rev. Pharmacol. Toxicol.***38**, 289–319. https://doi.org/10.1146/annurev.pharmtox.38.1.289 (1998). - 37.
Chang, Y. Y. Cyclic 3′,5′-adenosine monophosphate phosphodiesterase produced by the slime mold

*Dictyostelium discoideum*.*Science***161**, 57–59 (1968). - 38.
Merida, I., Avila-Flores, A. & Merino, E. Diacylglycerol kinases: at the hub of cell signalling.

*Biochem. J.***409**, 1–18. https://doi.org/10.1042/BJ20071040 (2008). - 39.
Naccarato, W. F., Ray, R. E. & Wells, W. W. Biosynthesis of myo-inositol in rat mammary gland. Isolation and properties of the enzymes.

*Arch. Biochem. Biophys.***164**, 194–201 (1974). - 40.
Kohout, T. A., Lin, F. S., Perry, S. J., Conner, D. A. & Lefkowitz, R. J. beta-Arrestin 1 and 2 differentially regulate heptahelical receptor signaling and trafficking.

*Proc. Natl. Acad. Sci. U. S. A.***98**, 1601–1606. https://doi.org/10.1073/pnas.041608198 (2001). - 41.
Paing, M. M., Stutts, A. B., Kohout, T. A., Lefkowitz, R. J. & Trejo, J. beta -Arrestins regulate protease-activated receptor-1 desensitization but not internalization or down-regulation.

*J. Biol. Chem.***277**, 1292–1300. https://doi.org/10.1074/jbc.M109160200 (2002). - 42.
Heding, A.

*et al.*The rat gonadotropin-releasing hormone receptor internalizes via a beta-arrestin-independent, but dynamin-dependent, pathway: addition of a carboxyl-terminal tail confers beta-arrestin dependency.*Endocrinology***141**, 299–306. https://doi.org/10.1210/endo.141.1.7269 (2000). - 43.
Heding, A.

*et al.*Gonadotropin-releasing hormone receptors with intracellular carboxyl-terminal tails undergo acute desensitization of total inositol phosphate production and exhibit accelerated internalization kinetics.*J. Biol. Chem.***273**, 11472–11477 (1998). - 44.
Davidson, J. S., Wakefield, I. K. & Millar, R. P. Absence of rapid desensitization of the mouse gonadotropin-releasing hormone receptor.

*Biochem. J.***300**(Pt 2), 299–302 (1994). - 45.
Davis, J. S., West, L. A. & Farese, R. V. Gonadotropin-releasing hormone (GnRH) rapidly stimulates the formation of inositol phosphates and diacyglycerol in rat granulosa cells: further evidence for the involvement of Ca

^{2+}and protein kinase C in the action of GnRH.*Endocrinology***118**, 2561–2571. https://doi.org/10.1210/endo-118-6-2561 (1986). - 46.
Rodbell, M., Lin, M. C. & Salomon, Y. Evidence for interdependent action of glucagon and nucleotides on the hepatic adenylate cyclase system.

*J. Biol. Chem.***249**, 59–65 (1974). - 47.
Pohl, S. L., Birnbaumer, L. & Rodbell, M. Glucagon-sensitive adenyl cylase in plasma membrane of hepatic parenchymal cells.

*Science***164**, 566–567 (1969). - 48.
Baker, J. G., Hall, I. P. & Hill, S. J. Temporal characteristics of cAMP response element-mediated gene transcription: requirement for sustained cAMP production.

*Mol. Pharmacol.***65**, 986–998. https://doi.org/10.1124/mol.65.4.986 (2004). - 49.
Withers, D. J.

*et al.*Rapamycin dissociates p70(S6K) activation from DNA synthesis stimulated by bombesin and insulin in Swiss 3T3 cells.*J. Biol. Chem.***272**, 2509–2514. https://doi.org/10.1074/jbc.272.4.2509 (1997). - 50.
Motulsky, H. J. Equation: Fitting a straight line with nonlinear regression. https://www.graphpad.com/guides/prism/8/curve-fitting/reg_linear_with_nonlinear.htm?q=equation+straight+line (2019)

- 51.
Motulsky, H. J. How standard errors are computed. https://www.graphpad.com/guides/prism/8/curve-fitting/reg_how_standard_errors_are_comput.htm?q=standard+error (2019).

- 52.
Motulsky, H. J. Standarad error of parameters. https://www.graphpad.com/guides/prism/8/curve-fitting/reg_standarad_error_of_parameters.htm?q=standard+error (2019).

- 53.
Violin, J. D.

*et al.*beta2-adrenergic receptor signaling and desensitization elucidated by quantitative modeling of real time cAMP dynamics.*J. Biol. Chem.***283**, 2949–2961. https://doi.org/10.1074/jbc.M707009200 (2008). - 54.
Hein, P.

*et al.*Gs activation is time-limiting in initiating receptor-mediated signaling.*J. Biol. Chem.***281**, 33345–33351. https://doi.org/10.1074/jbc.M606713200 (2006). - 55.
Traynor, J. R., Clark, M. J. & Remmers, A. E. Relationship between rate and extent of G protein activation: comparison between full and partial opioid agonists.

*J. Pharmacol. Exp. Ther.***300**, 157–161 (2002). - 56.
Fernandez, N.

*et al.*Reduction of G protein-coupled receptor kinase 2 expression in U-937 cells attenuates H2 histamine receptor desensitization and induces cell maturation.*Mol. Pharmacol.***62**, 1506–1514. https://doi.org/10.1124/mol.62.6.1506 (2002). - 57.
Tewson, P.

*et al.*Assay for detecting Galphai-mediated decreases in cAMP in living cells.*SLAS Discov.***23**, 898–906. https://doi.org/10.1177/2472555218786238 (2018). - 58.
Hannoush, R. N. Kinetics of Wnt-driven beta-catenin stabilization revealed by quantitative and temporal imaging.

*PLoS ONE***3**, e3498. https://doi.org/10.1371/journal.pone.0003498 (2008). - 59.
Zhang, Q., Pacheco, M. A. & Doupnik, C. A. Gating properties of GIRK channels activated by Galpha(o)- and Galpha(i)-coupled muscarinic m2 receptors in Xenopus oocytes: the role of receptor precoupling in RGS modulation.

*J. Physiol.***545**, 355–373. https://doi.org/10.1113/jphysiol.2002.032151 (2002). - 60.
Ojiaku, C. A.

*et al.*Transforming growth factor-beta1 decreases beta2-agonist-induced relaxation in human airway smooth muscle.*Am. J. Respir. Cell Mol. Biol.***61**, 209–218. https://doi.org/10.1165/rcmb.2018-0301OC (2019). - 61.
Motulsky, H. J. Equation: One phase association. https://www.graphpad.com/guides/prism/8/curve-fitting/index.htm?reg_exponential_association.htm (2019).

- 62.
Neubig, R. R., Spedding, M., Kenakin, T. & Christopoulos, A. International Union of Pharmacology Committee on Receptor Nomenclature and Drug Classification. XXXVIII. Update on terms and symbols in quantitative pharmacology.

*Pharmacol. Rev.***55**, 597–606. https://doi.org/10.1124/pr.55.4.4 (2003). - 63.
Mayersohn, M. & Gibaldi, M.

*Am. J. Pharm. Educ.***34**, 608–614 (1970). - 64.
Bourne, D. W. A. Oral dose adminsitration. https://www.boomer.org/c/p1/Ch08/Ch0803.html (2017).

- 65.
Violin, J. D., Dewire, S. M., Barnes, W. G. & Lefkowitz, R. J. G protein-coupled receptor kinase and beta-arrestin-mediated desensitization of the angiotensin II type 1A receptor elucidated by diacylglycerol dynamics.

*J. Biol. Chem.***281**, 36411–36419. https://doi.org/10.1074/jbc.M607956200 (2006). - 66.
Hofer, A. M., Landolfi, B., Debellis, L., Pozzan, T. & Curci, S. Free [Ca2+] dynamics measured in agonist-sensitive stores of single living intact cells: a new look at the refilling process.

*EMBO J.***17**, 1986–1995. https://doi.org/10.1093/emboj/17.7.1986 (1998). - 67.
Yu, R. & Hinkle, P. M. Desensitization of thyrotropin-releasing hormone receptor-mediated responses involves multiple steps.

*J. Biol. Chem.***272**, 28301–28307 (1997). - 68.
Berridge, M. J., Bootman, M. D. & Roderick, H. L. Calcium signalling: dynamics, homeostasis and remodelling.

*Nat. Rev. Mol. Cell Biol.***4**, 517–529. https://doi.org/10.1038/nrm1155 (2003). - 69.
Carafoli, E. Calcium pump of the plasma membrane.

*Physiol. Rev.***71**, 129–153. https://doi.org/10.1152/physrev.1991.71.1.129 (1991). - 70.
Merelli, F.

*et al.*Gonadotropin-releasing hormone-induced calcium signaling in clonal pituitary gonadotrophs.*Endocrinology***131**, 925–932. https://doi.org/10.1210/endo.131.2.1379169 (1992). - 71.
Liu, A. M.

*et al.*Galpha(16/z) chimeras efficiently link a wide range of G protein-coupled receptors to calcium mobilization.*J. Biomol. Screen***8**, 39–49. https://doi.org/10.1177/1087057102239665 (2003). - 72.
Tsien, R. W. & Tsien, R. Y. Calcium channels, stores, and oscillations.

*Annu. Rev. Cell Biol.***6**, 715–760. https://doi.org/10.1146/annurev.cb.06.110190.003435 (1990). - 73.
Clapham, D. E. Calcium signaling.

*Cell***80**, 259–268. https://doi.org/10.1016/0092-8674(95)90408-5 (1995). - 74.
Motulsky, H. J. Robust nonlinear regression. https://www.graphpad.com/guides/prism/8/curve-fitting/reg_fit_tab_2.htm?q=robust (2019).

- 75.
Putney, J. W. Jr. Capacitative calcium entry revisited.

*Cell Calcium***11**, 611–624 (1990). - 76.
Gesty-Palmer, D., El Shewy, H., Kohout, T. A. & Luttrell, L. M. beta-Arrestin 2 expression determines the transcriptional response to lysophosphatidic acid stimulation in murine embryo fibroblasts.

*J. Biol. Chem.***280**, 32157–32167. https://doi.org/10.1074/jbc.M507460200 (2005). - 77.
Namkung, Y.

*et al.*Functional selectivity profiling of the angiotensin II type 1 receptor using pathway-wide BRET signaling sensors.*Sci. Signal*https://doi.org/10.1126/scisignal.aat1631 (2018). - 78.
Motulsky, H. J. Method tab. https://www.graphpad.com/guides/prism/8/curve-fitting/reg_method-tab.htm (2019).

- 79.
Motulsky, H. J. Ambiguous. https://www.graphpad.com/guides/prism/8/curve-fitting/reg_analysischeck_nonlin_ambiguous.htm?q=ambiguous+fit (2019).

- 80.
Bichet, D.

*et al.*Vasopressin and oxytocin receptors (version 2019.4) in the IUPHAR/BPS guide to pharmacology database.*IUPHAR/BPS Guide Pharmacol.*https://doi.org/10.2218/gtopdb/F66/2019.4 (2019). - 81.
Feinstein, T. N.

*et al.*Noncanonical control of vasopressin receptor type 2 signaling by retromer and arrestin.*J. Biol. Chem.***288**, 27849–27860. https://doi.org/10.1074/jbc.M112.445098 (2013). - 82.
Tewson, P. H., Martinka, S., Shaner, N. C., Hughes, T. E. & Quinn, A. M. New DAG and cAMP sensors optimized for live-cell assays in automated laboratories.

*J. Biomol. Screen***21**, 298–305. https://doi.org/10.1177/1087057115618608 (2016). - 83.
Motulsky, H. J. Equation: log(agonist) vs. response—variable slope. www.graphpad.com/guides/prism/8/curve-fitting/index.htm?reg_dr_stim_variable.htm (2019).

- 84.
Falkenburger, B. H., Jensen, J. B. & Hille, B. Kinetics of M1 muscarinic receptor and G protein signaling to phospholipase C in living cells.

*J. Gen. Physiol.***135**, 81–97. https://doi.org/10.1085/jgp.200910344 (2010). - 85.
Riccobene, T. A., Omann, G. M. & Linderman, J. J. Modeling activation and desensitization of G-protein coupled receptors provides insight into ligand efficacy.

*J. Theor. Biol.***200**, 207–222. https://doi.org/10.1006/jtbi.1999.0988 (1999). - 86.
Waelbroeck, M. Activation of guanosine 5’-[gamma-(35)S]thio-triphosphate binding through M(1) muscarinic receptors in transfected Chinese hamster ovary cell membranes; 1. Mathematical analysis of catalytic G protein activation.

*Mol. Pharmacol.***59**, 875–885 (2001). - 87.
Bridge, L. J., Mead, J., Frattini, E., Winfield, I. & Ladds, G. Modelling and simulation of biased agonism dynamics at a G protein-coupled receptor.

*J. Theor. Biol.***442**, 44–65. https://doi.org/10.1016/j.jtbi.2018.01.010 (2018). - 88.
Xin, W., Tran, T. M., Richter, W., Clark, R. B. & Rich, T. C. Roles of GRK and PDE4 activities in the regulation of beta2 adrenergic signaling.

*J. Gen. Physiol.***131**, 349–364. https://doi.org/10.1085/jgp.200709881 (2008). - 89.
Gilman, A. G. G proteins: transducers of receptor-generated signals.

*Annu. Rev. Biochem.***56**, 615–649. https://doi.org/10.1146/annurev.bi.56.070187.003151 (1987). - 90.
Streb, H., Heslop, J. P., Irvine, R. F., Schulz, I. & Berridge, M. J. Relationship between secretagogue-induced Ca

^{2+}release and inositol polyphosphate production in permeabilized pancreatic acinar cells.*J. Biol. Chem.***260**, 7309–7315 (1985). - 91.
Motulsky, H. J. & Mahan, L. C. The kinetics of competitive radioligand binding predicted by the law of mass action.

*Mol. Pharmacol.***25**, 1–9 (1984). - 92.
Charlton, S. J. & Vauquelin, G. Elusive equilibrium: the challenge of interpreting receptor pharmacology using calcium assays.

*Br. J. Pharmacol.***161**, 1250–1265. https://doi.org/10.1111/j.1476-5381.2010.00863.x (2010). - 93.
Sum, C. S.

*et al.*Pharmacological characterization of GPCR agonists, antagonists, allosteric modulators and biased ligands from HTS hits to lead optimization. In*Assay Guidance Manual*(eds G. S. Sittampalam*et al.*) https://pubmed.ncbi.nlm.nih.gov/31693331/ (2020). - 94.
Motulsky, H. J. Equation: association kinetics (one ligand concentration). https://www.graphpad.com/guides/prism/8/curve-fitting/reg_association1.htm?q=association+kinetics (2019)

- 95.
Motulsky, H. J. How to: Plateau followed by one phase association. www.graphpad.com/guides/prism/8/curve-fitting/reg_exponential_plateau_then_association.htm?q=association (2019).

- 96.
Motulsky, H. J. Equation: Plateau followed by one phase decay. www.graphpad.com/guides/prism/8/curve-fitting/reg_exponential_decay_plateau.htm?q=decay+plateau (2019).

- 97.
Motulsky, H. J. Confidence tab. https://www.graphpad.com/guides/prism/8/curve-fitting/reg_confidence_tab.htm (2020).

## Acknowledgements

Research reported in this publication was supported by National Institute of General Medical Sciences, National Institutes of Health under Award Number R44GM125390 and National Institute of Neurologic Disorders and Stroke, National Institutes of Health R44NS082222. This content is solely the responsibility of the authors and does not necessarily represent the official views of the National Institutes of Health.

## Author information

### Affiliations

### Contributions

S.R.J.H. and L.J.B. designed the data analysis framework and analyzed data. P.H.T, A.M.Q. and T.E.H. conceived, designed and executed experiments and analyzed data. S.R.J.H wrote the manuscript. All authors reviewed the manuscript.

### Corresponding author

## Ethics declarations

### Competing interests

The authors declare no competing interests.

## Additional information

### Publisher's note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

## Rights and permissions

**Open Access** This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The images or other third party material in this article are included in the article’s Creative Commons license, unless indicated otherwise in a credit line to the material. If material is not included in the article’s Creative Commons license and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this license, visit http://creativecommons.org/licenses/by/4.0/.

## About this article

### Cite this article

Hoare, S.R.J., Tewson, P.H., Quinn, A.M. *et al.* Analyzing kinetic signaling data for G-protein-coupled receptors.
*Sci Rep* **10, **12263 (2020). https://doi.org/10.1038/s41598-020-67844-3

Received:

Accepted:

Published:

## Comments

By submitting a comment you agree to abide by our Terms and Community Guidelines. If you find something abusive or that does not comply with our terms or guidelines please flag it as inappropriate.