Introduction

Signaling receptors routinely encounter a wide variety of extracellular ligands and decode their identity with remarkable precision to generate context-specific responses. This selective processing of environmental cues is essential for regulating diverse biological processes, including development, immune surveillance, and tissue homeostasis (Cantley et al., 2014). Failures in ligand discrimination underlie many diseases, including diabetes and cancer (Madsen et al., 2025; Madsen and Vanhaesebroeck, 2020).

A key determinant of ligand specificity in biochemical networks is the thermodynamic stability of molecular complexes, such as ligand–receptor or substrate–enzyme pairs. At thermal equilibrium, the abundance of complexes is determined by their equilibrium binding constants. This imposes a fundamental limit on specificity: high-affinity ligands are inevitably favored over lower-affinity competitors, with complex abundances scaling in proportion to their association constants. Notably, many biochemical networks display paradoxical behaviors that cannot be explained by equilibrium affinity alone (Clark et al., 1999; Coombs et al., 2002; Freed et al., 2017; Madsen et al., 2025; Myers et al., 2023). For example, signaling receptors such as receptor tyrosine kinases (RTKs) and T cell receptors can produce stronger signaling outputs (phosphorylation levels) in response to intermediate-affinity ligands compared to low- and high-affinity ligands (Coombs et al., 2002; Lever et al., 2014; Freed et al., 2017; Madsen et al., 2025; Myers et al., 2023). Additionally, RTKs also exhibit a non-monotonic dependence between receptor activity and kinase activity (Kiyatkin et al., 2020; Kleiman et al., 2011). These observations raise a fundamental question: how do signaling receptors overcome thermodynamic constraints to achieve robust, ligand-specific responses?

A classic scheme to bypass limitations imposed by equilibrium thermodynamics is kinetic proof-reading (KPR), a mechanism first proposed by Hopfield (Hopfield, 1974) and Ninio (Ninio, 1975). KPR enhances specificity of high affinity ligands by introducing energy-consuming, irreversible steps, such as phosphorylation/dephosphorylation cycles, that amplify differences between competing ligands. KPR has been invoked in diverse systems, including DNA replication (Hopfield, 1980), mRNA surveillance (Hilleren and Parker, 1999), protein folding (Gulukota and Wolynes, 1994), and immune receptor signaling (McKeithan, 1995; Huang et al., 2019; Lever et al., 2014). Notably, while most KPR models prefer ligands with the highest affinity, it is also known that embedding KPR schemes in larger biochemical networks may allow non-monotonic dependence between ligand affinity and network activity (Lever et al., 2014; Murugan et al., 2014). However, as we will show below, these models do not capture the non-monotonic dependence between network output and kinase activity.

In our model, high-affinity ligand-receptor complexes are sorted towards degradation-prone states, while low-affinity complexes repeatedly dissociate the ligand, resulting in maximal signaling output only from intermediate-affinity ligands. Notably, this ligand specificity can be tuned by varying easily controllable cellular parameters, e.g. enzyme abundances. This non-equilibrium kinetic sorting mechanism explains the paradoxical non-monotonic dependence of signaling activity on ligand affinity and phosphorylation rate observed in RTKs. More broadly, given the ubiquity of the signaling motifs involved, we propose that kinetic sorting provides a general strategy for achieving ligand discrimination that is likely to be broadly used across diverse signaling networks.

Results

Classic kinetic proofreading favors high-affinity ligands

Kinetic proofreading (KPR) is the standard model for non-equilibrium ligand discrimination. To set the stage, we first revisited the classic KPR model originally proposed by McKeithan to explain how T cell receptors avoid activation downstream of weak ligands (McKeithan, 1995) (Fig. 1a; see Supplementary Materials for equations).

Reaction scheme of kinetic proofreading models.

Chemical species and rate constants are shown in the figure. R denote ligand-free receptors, B denote ligand-bound inactive receptors, and Pn, n ∈ [1, N] are phosphorylated receptors. The ultimate phosphorylated species PN (marked red) is assumed to be signaling competent. (a) shows the traditional model first proposed by McKeithan (McKeithan, 1995). (b, c) show the sustained signaling model and the limited signaling model (Lever et al., 2014) which introduce additional receptor states, and I respectively, directly following receptor activation.

In this model, ligand-bound receptors undergo a series of phosphorylation steps, with the final state PN representing the active, signaling-competent form. Importantly, ligand unbinding at any phosphorylation stage returns the receptor to the unbound state R. We parameterized the model using dimensionless quantities: the ligand dissociation rate δ = kdτ, phosphorylation rate ω = kpτ, and ligand concentration u = L/KD, where KD = kd/kon. Assuming saturating ligand (u → ∞), the steady-state abundance of the active state is:

As expected, increasing the phosphorylation cascade length N amplifies the preference for low-dissociation (high-affinity) ligands (Fig. 2a), reflecting the classical KPR outcome.

Ligand discrimination in kinetic proofreading models.

(a) Activity PN plotted as a function of non-dimensional ligand dissociation rate δ for the traditional KPR scheme (Fig. 1a). (b) Activity PN plotted as a function of non-dimensional ligand dissociation rate δ for the limited signaling model (Fig. 1b). (c) The dependence of the activity on the dimensionless phosphorylation rate ω for the limited signaling model. All figures plotted for a sequence of N = 1, 5, and 10 phosphorylation sites.

Modified KPR schemes do not explain paradoxical RTK behavior

Before introducing our model, we briefly review two previously proposed extensions of receptor-level KPR that exhibit non monotonic ligand discrimination: the sustained signaling model and the limited signaling model (Lever et al., 2014) (Fig. 1b,c). Both models introduce an additional state to Mckeithan’s KPR scheme. The sustained signaling model adds an active but ligand-free state , while the limited signaling model introduces an inactivated state I downstream of PN.

While both models show non-monotonic dependence of signaling activity on ligand affinity (Lever et al., 2014), only the limited signaling model retains this non-monotonic dependence at saturating ligand concentrations (Lever et al., 2014) (Fig. 2b), consistent with some paradoxical features observed in RTKs (Freed et al., 2017; Madsen et al., 2025; Myers et al., 2023). However, the limited signaling model fails to reproduce a second key observation in RTKs: receptor activity in this model increases monotonically with kinase activity, whereas RTK experiments show that partial kinase inhibition can paradoxically increase receptor activity (Kiyatkin et al., 2020; Kleiman et al., 2011) (Fig. 2c). Thus, these models are insufficient to explain RTK signaling dynamics.

Notably, these models neglect a key feature of many receptor signaling pathways; preferential degradation of activated receptors (Sorkin and Goh, 2009; Koenig and Edwardson, 1997; Charpentier and King, 2021; Cendrowski et al., 2016). Below, we incorporate preferential degradation in our model to investigate how it governs receptor activity.

A kinetic sorting model integrates active receptor degradation

We build a model to study the effect of two widespread signaling motifs: sequential multi-site phosphorylation and ligand-induced receptor degradation (Fig. 3) on ligand discrimination. In our model, receptors are delivered to the surface at a constant rate, internalized at a basal rate kint, and degraded more rapidly when highly phosphorylated . Ligand-bound receptors undergo irreversible phosphorylation and dephosphorylation through distinct irreversible mechanisms. We note that both kinase and phosphatase are irreversible reactions carried out by separate enzymes. While their effect on the coarse-grained model of the receptor may appear reversible, it is important to note that receptor phosphorylation via ATP hydrolysis and removal of the phosphate group from the receptor corresponds to a futile cycle that does not recharge the ADP molecule to an ATP molecule. In addition to the previously defined dimensionless parameters, we define the dimensionless active receptor degradation rate, , and the relative rate of dephosphorylation, ρ = kdp/kp. A key feature of our model is that all phosphorylated species are signaling competent. Indeed, in many signaling pathways all phosphorylation sites on the receptor (Schulze et al., 2005; Tobin, 2008; Kollewe et al., 2004; Lemmon and Schlessinger, 2010; Latorraca et al., 2020) have downstream effects. Therefore, we define the net activity An of phosphorylation site n as all receptor states where the site n is phosphorylated: An = ∑mn Pm.

Reaction scheme of kinetic sorting model.

Chemical species and rate constants are shown in the figure. R denote ligand-free receptors, B denote ligand-bound inactive receptors, and Pn, n ∈ [1, N] are phosphorylated receptors. ϕ represents an implicit source and sink, corresponding to receptor delivery and internalization, respectively. It does not denote a physical chemical species.

Parameter ranges

To ensure that the phenomena captured by our model are relevant to real signaling networks, we selected ranges for the dimensionless parameters based on direct experimental measurements and model fits. Importantly, many of these kinetic processes have comparable rates across diverse receptor systems (Koenig and Edwardson, 1997; Subtil et al., 1994; Liu et al., 2000). Specifically, basal receptor internalization occurs at rates of kint ≈ 10−4–10−3, s−1(Wiley, 2003), while active receptor internalization is typically faster, at (Wiley, 2003; Lyashenko et al., 2020). Ligand dissociation rates typically fall in the range kd ≈ 10−2–10−1, s−1(Chen et al., 2009; Lyashenko et al., 2020), and receptor phosphorylation (kp) and dephosphorylation (kdp) occur at ~ 10−1–100, s−1(Kleiman et al., 2011; Chen et al., 2009; Lyashenko et al., 2020). For EGFR, equilibrium dissociation constants range from ~ 0.1, nM for the high-affinity ligand Beta-cellulin to ~ 25, nM for the low-affinity ligand AREG (Hu et al., 2022; Macdonald-Obermann and Pike, 2014). Based on these values, we set the following ranges for dimensionless parameters: , and δ = kd/kint ∈ [1, 1000]. Finally, the number of phosphorylation sites with known functional roles typically ranges from 5 to 25 (Schulze et al., 2005). These broad ranges comfortably encompass experimentally measured estimates. Unless otherwise specified, our default parameter values are δ = 20, ω = 200, ρ = 0.01, β = 50, and N = 10.

Before examining how phosphorylation levels depend on model parameters, we illustrate the mechanism of kinetic sorting of receptor states, which tunes ligand specificity beyond pure thermo-dynamic preference, using a simple example. To that end, we consider a signaling network with N = 5 phosphorylation sites interacting with three ligands of distinct affinities—high, medium, and low. We assume the dissociation rates for these ligands are δH = 20, δM = 200, and δL = 1000, respectively. In order to compare our model with aforementioned paradoxical experimental observations which have been performed at saturating ligand concentration, we take the limit u → ∞.

Fig. 4 shows that low affinity ligands (δL = 1000) predominantly sort receptors towards the inactive state B and early phosphorylation states Pn, n ~ 1 as frequent ligand unbinding prevents progression to later phosphorylation states. This behavior resembles the traditional KPR mechanism described by McKeithan (McKeithan, 1995). In contrast, receptors bound to high affinity ligands are sorted toward later phosphorylation states, which mark them for enhanced degradation. Here, similar to traditional KPR, the fraction of receptors reaching the final phosphorylation state is highest for high affinity ligands. Yet, the overall receptor pool is reduced due to ligand-induced degradation, lowering net phosphorylation activity. Strikingly, receptors bound to intermediate affinity ligands (δM = 200) are sorted towards intermediate phosphorylation states, resulting in maximal phosphorylation output. Below, we show how kinetic parameters govern the ability of the network to overcome thermodynamic preference and acquire ligand specificity.

Kinetic sorting of receptor species.

Abundances of networks species B (ligand bound inactive receptor) and Pn, n ∈ [1, 5] for a signaling receptor with N = 5 phosphorylation sites. Abundances are shown for ligands of three different affinities. The inset shows the activity of the first phosphorylation site A1. Species abundances below 10−3 are not shown.

Early phosphorylation sites show ligand-specificity

Figure 5a illustrates how total phospho-rylation activity at each site, An, n ∈ [1, N] varies with ligand dissociation rate δ. We note that the activity of the nth site is given by the total concentration of all species that have the nth site phos-phorylated;. We find that early phosphorylation sites (n ~ 1) exhibit maximal activity at intermediate values of δ while both high- and low-affinity ligands suppress net receptor phosphorylation. Our model predicts that this ligand-specificity diminishes for later sites, where outputs increasingly resemble traditional KPR which favors high-affinity ligands.

Kinetic sorting model predicts ligand specificity.

(a) The activity An of the nth phosphorylation site as a function of dimensionless dissociation rate δ. The activity is normalized to the maximum activity. The maximum An as a function of n is shown in the inset. (b) Activity of the first phosphorylation site A1 plotted as a function of the dissociation rate δ for different values of the phosphorylation rate ω. (c, d) Activity of the first phosphorylation site A1 plotted as a function of phosphorylation rate ω (dephosphorylation rate ρ in panel d) for different values of the dissociation rate δ.

To examine how model parameters shape ligand specificity, we focused on the activity at the first phosphorylation site, A1, which exhibits the strongest discriminatory behavior (Fig.5a). As shown in Fig.5b, achieving ligand specificity at high dissociation rates δ requires sufficiently high phosphorylation rates ω. Notably, our model captures a puzzling observation from EGFR signaling: the high-affinity ligand EGF produces lower/comparable steady-state phosphorylation compared to lower-affinity ligands such as Epigen and Epiregulin (Freed et al., 2017; Myers et al., 2023; Madsen et al., 2025). Experimental estimates place the basal EGFR internalization rate at kint ≈ 1.3 × 10−3, s−1 (Chen et al., 2009), the EGF dissociation rate at kd ≈ 3 × 10−2, s−1 (Chen et al., 2009), and the phosphorylation rate at kp ≈ 10−1 −100, s−1, yielding δEGF ≈ 10–20 and ωEGFR ≈ 100–1000. Low-affinity ligands such as Epigen (EPGN) and Epiregulin (EREG) have equilibrium dissociation constants about 10-fold higher than EGF (Hu et al., 2022), corresponding to δEPGNδEREG ≈ 100–200. The effective degradation rate of fully activated receptors is estimated to be 10–50 times higher than that of inactive receptors (Lyashenko et al., 2020), implying β = 50. Under these conditions, our model predicts a switch in phosphorylation levels: as δ increases from δEGF to δEPGN, receptor phosphorylation increases—reversing the expectation based purely on thermodynamic affinity. This effect arises because EGF-bound receptors are efficiently sorted towards degradation-prone states compared to those bound to lower-affinity ligands.

Our model also explains another paradox in EGFR signaling. Experimental studies have shown that EGF-stimulated receptors exhibit higher steady-state phosphorylation when kinase activity is partially inhibited (Kiyatkin et al., 2020; Kleiman et al., 2011). As shown in Fig.5c, at low δ values (e.g., δ = 16), decreasing the phosphorylation rate ω from levels typical of EGFR (ωEGFR ≈ 100–1000) paradoxically increases overall receptor phosphorylation. A similar effect is observed when receptor dephosphorylation is enhanced (Fig.5d). Importantly, our model makes a testable prediction: the reversal of thermodynamic preference observed between EGF and EPGN/EREG will disappear when kinase activity is mildly suppressed (see, e.g., the curves for ω = 256 and ω = 16 over δ ∈ [10, 100]), such as by treatment with low doses of the kinase inhibitor gefitinib (Herbst et al., 2004). This non-monotonic trend may help prevent cells with abnormally high kinase activity from becoming constitutively active, thereby preserving their sensitivity to extracellular cues.

Multi-site phosphorylation and ligand-induced degradation are both essential for ligand-specificity

To assess the importance of sequential multi-site phosphorylation on ligand specificity, we analyzed , the phosphorylation of the first site for signaling networks with N phosphorylation sites. Fig. 6a shows that multi-site phosphorylation is essential to endow signaling networks with ligand specificity and ligand-induced receptor degradation alone is not sufficient. This is because the non-monotonic preference for intermediate affinity ligands arises only when the receptors can be sorted among multiple phosphorylation sites: earlier ones for low affinity ligands and later ones for high affinity ligands.

Multiple phosphorylation sites and receptor degradation dictate ligand specificity.

(a) Activity of the first phosphorylation site, A1, as a function of the dissociation rate δ for signaling networks with different number of phosphorylation sites. (b) The optimal dissociation rate δopt that leads to maximum phosphorylation activity as a function of dimensionless degradation rate β for different values of ω. δopt is shown only if δopt ∈ [1, 1000]. (c) The relative activity of a ligand with dissociation rate that differs by kBT compared to δopt plotted as a function of β for different values of ω (see inset). Of the two ligands that differ in stability by kBT, the ligand exhibiting maximum activity is considered.

To assess how receptor degradation shapes ligand specificity for a multi-site phosphorylation network, we examined how altering receptor turnover influences model behavior. As shown in Fig. 6b, the optimal dissociation rate δopt, which maximizes receptor phosphorylation levels, increases with ligand-induced degradation rate β. Crucially, this optimal δopt emerges only when receptor degradation is strong (β ≫ 1). These predictions can be tested by blocking receptor degradation, e.g., via mutation of ubiquitination sites (Gerritsen et al., 2023).

To quantify ligand specificity, we computed receptor phosphorylation in response to ligands differing by at least one kBT in binding free energy from the optimal ligand. Figure 6c shows that as β increases, phosphorylation downstream of suboptimal ligands (red line in inset) declines relative to the optimal ligand. This enhanced specificity is further amplified by increasing kinase activity ω.

These results show that both multi-site phosphorylation and ligand induced degradation are key features controlling ligand specificity in our kinetic sorting mechanism.

Discussion

Cells face the formidable task of decoding multiple chemically distinct extracellular signals to generate appropriate, context-specific responses. This challenge is especially acute for cell surface receptors like RTKs, GPCRs, and interleukin receptors, which bind multiple cognate ligands and yet elicit distinct downstream outcomes. While equilibrium affinity provides a baseline expectation for ligand specificity, it cannot fully explain the rich and often counterintuitive behaviors observed in many signaling systems.

Here, we show that a non-equilibrium mechanism of kinetic sorting which operates through multi-site phosphorylation and active receptor degradation can explain how signaling networks achieve ligand specificity beyond equilibrium limits. In kinetic sorting, high-affinity ligand-receptor complexes are sorted toward degradation-prone states, low-affinity complexes are sorted towards inactivated states, and intermediate-affinity ligands strike the optimal balance between progression and degradation to maximize signaling. This framework explains paradoxical features observed in RTK systems, including the non-monotonic dependence of phosphorylation on ligand affinity and kinase activity. Importantly, our model predicts that early phosphorylation sites show the strongest ligand discrimination, consistent with recent experimental observations. It also makes the testable prediction that impairing receptor degradation should reduce specificity by eliminating the kinetic sorting effect. Given the ubiquity of the essential motifs of our mechanism, i.e., multi-site phosphorylation and receptor degradation, we believe that kinetic sorting may be a common mechanism to modulate ligand specificity at the receptor level, potentially in addition to other mechanisms that endow signaling networks with ligand-specificity, both at the receptor level (Lever et al., 2014) as well as in downstream signaling pathways (Singh and Nemenman, 2017).

In contrast to what has been shown previously for KPR models (Coombs et al., 2002; Lever et al., 2014), the kinetic sorting model also captures the non-monotonic relationship between signaling output and kinase/phosphatase activity observed in RTK systems such as EGFR (Kleiman et al., 2011; Kiyatkin et al., 2020). In these systems, partial inhibition of kinase activity paradoxically increases steady-state receptor phosphorylation, a behavior not accounted for by equilibrium models (see Supplementary Materials) or by prior non-equilibrium schemes such as the limited signaling model (Lever et al., 2014). This type of protective filtering can ensure that downstream signaling remains contingent on extracellular cues and is not constitutively active, thereby preventing persistent, cue-independent activation. Such regulation could help maintain control in pathways such as those governing growth, where deregulated activity can have severe consequences. The potential benefit of this regulatory pattern suggests it could be advantageous in other signaling contexts. Consistent with this idea, non-monotonic regulation by kinase or phosphatase activity is found in other systems through distinct mechanisms (e.g., the non-monotonic effects of the phosphatase CD45 in T-cell receptor signaling (Courtney et al., 2019)). This indicates that selective filtering based on enzymatic activity is a strategy employed in diverse biological settings. While direct evidence for the kinetic sorting mechanism remains limited to RTKs, similar filtering behavior emerges in theoretical analyses of phosphorylation–dephosphorylation cycles in more general settings (Martins and Swain, 2013), suggesting it may represent a broader principle of enzymatic signaling networks.

Our findings complement prior studies on mechanisms of ligand specificity that operate at thermal equilibrium, such as those described in the Bone Morphogenetic Protein (BMP) pathway (Antebi et al., 2017; Su et al., 2022; Parres-Gold et al., 2025). BMP signaling relies on promiscuous ligand–receptor interactions, with specificity emerging from differences in receptor abundance, binding affinity, and complex activity. In contrast, our work shows that non-equilibrium mechanisms—such as phosphorylation cycles and ligand-induced receptor degradation—can achieve ligand discrimination even for a single receptor type. Given that ligand–receptor promiscuity, multisite phosphorylation, and receptor turnover are common features across signaling systems (e.g., in the EGFR/ErbB family (Linggi and Carpenter, 2006)), it is likely that biological networks integrate both equilibrium and non-equilibrium strategies to achieve robust and tunable ligand specificity.

In recent years, there has been growing interest in engineering synthetic physical and chemical circuits capable of carrying out complex computational tasks, including input discrimination, classification, prediction, and the generation of multiple stable cell states (Shakiba et al., 2021; Ma et al., 2022; Benzinger et al., 2022; Zhu et al., 2022; Floyd et al., 2024; Parres-Gold et al., 2025; Aoki et al., 2019). Some of these synthetic strategies rely on equilibrium thermodynamics (Parres-Gold et al., 2025), while others exploit non-equilibrium steady states (Floyd et al., 2024). We propose that non-equilibrium kinetic sorting, which harnesses receptor synthesis and degradation, could provide synthetic biologists with a powerful framework for achieving precise control over molecular abundances and dynamic system behavior.

Finally, we address a major concern in non-equilibrium signaling circuits: the energetic cost of operation. Previous theoretical work has shown that free energy dissipation places fundamental constraints on the performance of signaling networks (Bryant and Machta, 2023; Govern and ten Wolde, 2014; Lan et al., 2012; Mehta and Schwab, 2012; Qian and Reluga, 2005; Cao et al., 2015; Azeloglu and Iyengar, 2015; Floyd et al., 2024; Mahdavi et al., 2024). These studies typically focus on futile cycles of reversible modifications such as phosphorylation or methylation. In contrast, ligand-induced receptor degradation—a central feature of many signaling networks—is a far more energy-intensive process. For example, MCF10A cells maintain approximately 105 EGFR molecules on the surface (each 1,210 amino acids in length) (Shi et al., 2016), with a synthesis rate of about 15 receptors per second(Lyashenko et al., 2020), corresponding to an energetic cost of roughly ~ 8 × 104 ATP/sec (assuming 4.5 ATP per peptide bond (Milo et al., 2010)). By comparison, EGFR dephosphorylation occurs over ~ 15 seconds (Kleiman et al., 2011), and only 5 − 10% of receptors are phosphorylated at steady state (Shi et al., 2016; Feng et al., 2023), resulting in a much lower energetic cost of ~ 6 × 102 ATP/sec for dephosphorylation. Thus, the energetic burden of receptor turnover can exceed that of reversible modification cycles by up to two orders of magnitude. These estimates suggest that, at least in eukaryotic cells where signaling proteins may turnover multiple times within cellular lifetime (Milo et al., 2010), non-equilibrium modification cycles are unlikely to pose a fundamental energetic limitation on the functionality of signaling networks. Here, the energetic demands of signaling networks must account for protein turnover in addition to non-equilibrium modification cycles.

Supplementary materials

Equations for proofreading models

The equations describing species abundances in the traditional KPR model similar to that of McK-eithan (McKeithan, 1995) are as follows:

For the limited signaling model, the dynamics of B, and Pi, i ∈ [1, N − 1] are identical to the traditional KPR model. The dynamics of R and PN are modified as follows:

Equations for the model with receptor degradation

Signaling receptors participate in a variety of complex regulatory processes, including non-linear ligand binding dynamics (Limbird et al., 1975; Macdonald and Pike, 2008), receptor oligomerization (Mudumbi et al., 2024; Huang et al., 2016), context-specific interactions with adapter proteins (Madsen and Vanhaesebroeck, 2020; Feng et al., 2023), and trafficking between cellular compartments leading to degradation (Sorkin and Goh, 2009; Wiley, 2003; Irannejad and Von Zastrow, 2014).

While computational models that incorporate these mechanistic details are powerful tools for hypothesis generation (Chen et al., 2010; Qiao et al., 2025), they often require large-scale datasets for accurate parameterization (Feng et al., 2023). As an alternative, simplified models that intentionally omit certain mechanistic details can still yield deep qualitative insights, even if they cannot quantitatively reproduce experimental data.

In this study, we present such a simplified model aimed at explaining two paradoxical features of receptor tyrosine kinase (RTK) signaling: (1) the non-monotonic relationship between ligandreceptor affinity and steady-state receptor phosphorylation (Freed et al., 2017; Madsen et al., 2025; Myers et al., 2023), and (2) the counterintuitive increase in receptor phosphorylation following mild kinase inhibition (Kleiman et al., 2011; Kiyatkin et al., 2020).

To keep the model simple and tractable, we neglect receptor recycling and oligomerization. Previously, we showed that the combined effects of endocytosis, recycling, and degradation can be captured by a single effective dimensionless parameter, β in this study, which reflects the degradation bias of fully phosphorylated receptors compared to partially phosphorylated receptors(Lyashenko et al., 2020). Similarly, receptor dimerization and negative cooperativity can be abstracted into a Hill coefficient η < 1 (Lyashenko et al., 2020). For the phenomena explored here, including oligomerization would modify the shape of the response curves but not their qualitative behavior.

Under these assumptions, the governing equations for the model are given by:

All equations are solved at steady state and in the limit u → ∞. All codes required to generate the figures in the manuscript can be found at https://github.com/BarriosJer0/KineticSorting.

Equations for a model at thermal equilibrium

To confirm the role of non-equilibrium thermodynamics on ligand specificity, we consider the closest equivalent equilibrium model. The strongest requirement of an equilibrium model is that all reactions must be bidirectional. Another requirement is that microscopic reversibility or detailed balance. Specifically, ratios of rate constants around loops must equal to unity for all loops. The first requirement implies that unidirectional reactions: synthesis and degradation of receptors and the irreversible loss of activity due to ligand dissociation cannot exist in a reaction network that operates at equilibrium. The simplest equilibrium model closest to the kinetic sorting scheme is governed by the following equations:

In the above equations, we use the notation P0B.

We note that phosphorylation/dephosphorylation reactions are unidirectional non-equilibrium reactions carried out by different enzymes: Phosphorylation hydrolyses ATP to ADP and attaches a phosphate group to the receptor. In contrast, while dephosphorylation removes phosphate group from the receptor, it does not recharge an ADP molecule back to ATP. Notably, however, this non-equilibrium nature of the phosphorylation/dephosphorylation cycle is not apparent in our coarse-grained kinetic scheme where ATP and ADP are not explicitly considered. We retain this part of the non-equilibrium model since a corresponding equilibrium model can be imagined where different sites on the receptor change conformation between an inactive and an active state and that these changes occur in a sequential manner.

Solving these equations at steady state and taking the limit u = Lkon/kd → ∞, we have

where RT is the total number of receptors and ρ = kdp/kp. Note that as expected, this equilibrium model has no dependence on ligand dissociation rate kd at saturation, further confirming that non-equilibrium reactions are needed to endow cells with ligand specificity.

Additional information

Funding

National Institutes of Health (R35GM142547)

Sir Henry Wellcome Fellowship (220464/Z/20/Z)