A quantitative view of strategies to engineer cell-selective ligand binding

A critical property of many therapies is their selective binding to specific target populations. Exceptional specificity can arise from high-affinity binding to unique cell surface targets. In many cases, however, therapeutic targets are only expressed at subtly different levels relative to off-target cells. More complex binding strategies have been developed to overcome this limitation, including multi-specific and multi-valent molecules, but these create a combinatorial explosion of design possibilities. Therefore, guiding strategies for developing cell-specific binding are critical to employ these tools. Here, we extend a multi-valent binding model to multi-ligand and multi-receptor interactions. Using this model, we explore a series of mechanisms to engineer cell selectivity, including mixtures of molecules, affinity adjustments, and valency changes. Each of these strategies maximizes selectivity in distinct cases, leading to synergistic improvements when used in combination. Finally, we identify situations in which selectivity cannot be derived through passive binding alone to highlight areas in need of new developments. In total, this work uses a quantitative model to unify a comprehensive set of design guidelines for engineering cell-specific therapies. Summary points Affinity, valency, and other alterations to target cell binding provide enhanced selectivity in specific situations. Evidence for the effectiveness and limitations of each strategy are abundant within the drug development literature. Combining strategies can offer enhanced selectivity. A simple, multivalent ligand-receptor binding model can help to direct therapeutic engineering.


Introduction
Many drugs both derive their therapeutic bene!t and avoid toxicity through selective binding to speci!c cells within the body. Often, target cells can di#er from o#-target populations only subtly in surface receptor expression, making selective activation of target cells di"cult to achieve. Even with a drug of very speci!c molecular binding, genetic and non-genetic heterogeneity can create a wide distribution of cell responses. This can result in reduced e#ectiveness and increased toxicity. Speci!c cell targeting is a universal challenge in protein-based therapies. For example, in cancer, resistance to anti-tumor antibodies, 1 targeted inhibitors, 2 chemotherapies, 3 and chimeric antigen receptor T cells 4,5 all can arise through the selection of target cells among heterogeneous cell populations. While the immune system takes advantage of heterogeneity at the single-cell level to translate noisy in$ammatory signals into robust yet sensitive responses, 6 this heterogeneity impedes our e#ort of creating a highly speci!c drug. The intricacies of both inter-population receptor expression di#erences and intra-population receptor expression heterogeneity present signi!cant challenges that limit the selectivity of therapies within the body.
Further improving cell-speci!c targeting will require new strategies of engineering speci!city. Noncellular therapies such as protein therapies have most extensively been engineered to target speci!c cell types through mutations that provide high-a"nity binding to unique surface antigens. 7 This Figure 1: A model system for exploring the factors contributing to cell selectivity. a) A simpli!ed schematic of the binding model. There are two types of receptors and two types of ligand monomers that form a tetravalent complex. b) A cartoon for four cell populations expressing two di#erent receptors at low or high amounts. c) A sample heat/contour map for the modelpredicted log ligand bound given the expression of two types of receptors. d) Eight arbitrary theoretical cell populations with various receptor expression pro!les.
Here we investigate cell-speci!c targeting quantitatively by extending a multi-valent, multi-ligand equilibrium binding model. Virtually any therapy, including monoclonal antibodies, small-molecule inhibitors, and cytokines, can be thought of as ligands for respective receptors expressed on target cells. Ligand binding to target cells is a necessary !rst step, and essential for a drug's intended action. In contrast, binding to o#-target cells can result in unintended e#ects or toxicity. Some cell populations can be distinguished by their expression of a unique receptor not expressed by other populations, but more commonly, target and o#-target cells express the same collection of receptors and di#er only in their magnitudes of receptor expression. In such situations, engineering drugs to optimize target cell binding while minimizing their binding to o#-target cells is an area of ongoing research, and has inspired a myriad of drug design strategies. [12][13][14] In this work, we de!ne cell population selectivity as the ratio of the number of ligands bound to target cell populations divided by the number of ligands bound to o#-target cell populations. We will use a quantitative binding estimation for each cell population to examine these strategies.
While ligand-receptor binding events in biology are diverse, they are governed by thermodynamic properties and the law of mass action. For the binding between monomer ligands and receptors, their a"nity can be described by the association constant , or its reciprocal, the dissociation constant . Binding behavior is more complicated when the ligands are multivalent complexes consisting of multiple units, each of which can bind to a receptor (Fig. 1a). During initial association, we assume that the !rst subunit on a ligand complex binds according to the same dynamics that govern monovalent binding. Subsequent binding events exhibit di#erent behavior, however, due to the increased local concentration of the complex and steric e#ects. Here, we assume that the e#ective association constant for the subsequent bindings is proportional to that of the free binding, but scaled by a crosslinking constant, , which describes how easily a multivalent ligand bound to a cell monovalently can attain secondary binding. Another consideration that must be made when modeling multivalent binding processes is whether ligand complexes are formed via random assortment of monomers, or whether the monomer composition is uniform across complexes by engineering. We developed a multivalent binding model that calculates the amount of ligand bound at equilibrium taking each of these factors into account (see methods).
As a simpli!cation, we will consider theoretical cell populations that express only two receptors capable of binding ligand (Fig. 1b), ranging in abundance from 100 to 10,000 per cell. Figure 1c shows the logscaled predicted amount of binding of a monovalent ligand given the abundance of two receptors, with the concentration of ligand, , to be , and its dissociation constants to the two receptors to be and , respectively. Because all axes are log-scaled, the number of contour lines between two points indicates the ratio of ligand binding between populations. For instance, in !gure 1c, cell populations at points 1 and 2 are on the same contour line and thus have the same amount of ligand bound; the cell populations at points 1 and 3 are separated by multiple contour lines, indicating that cells at point 3 bind more ligands (In fact, the ratio can be read as the exponent of the contour line di#erence. For point 3 to point 1, the ratio is one point to another as a change of expression pro!le for a cell population. This situation might correspond to some cues inducing expression of a receptor, such as interferon-induced upregulation of MHC and regulatory T cell upregulation of IL-2R⍺. 15 When the amount of receptor 1 ( ) increases on a cell (a point moves rightward, e.g. from 1 to 2), the amount of binding doesn't increase signi!cantly. On the contrary, a cell with increased expression of (moving upward, e.g. from 1 to 3) will bind signi!cantly more ligands. Here, the ligand's high a"nity for and relatively low a"nity for lead to binding varying more strongly with changes to expression than .
To analyze more general cases, we arbitrarily de!ned eight theoretical cell populations according to their di#erential expression of two receptor types ( and plotted on x and y axes). As shown in Fig  1d, they either have high ( ), medium ( ), or low ( ) expression of and . The receptor expression pro!le within each cell population can also vary widely. To demonstrate cell-to-cell heterogeneity, we also de!ned intrapopulation variability of each population arbitrarily. For instance, the expression pro!le of has a wider range. The intrapopulation variance will be accounted for using sigma point !lters. 16 We will use this binding model to examine how engineering a ligand using various strategies can improve cell-speci!c targeting. Although we will only consider two receptor and ligand subunit types respectively, the principles we present can generalize to more complex cases. populations exposed to a monovalent ligand with dissociation constants to receptor 1 and 2 ranging from , at a concentration . Ligand bound ratio of (b) to , (c) to , (d) to , and (e) to .

A!nity Provides Selectivity Toward Cell Populations with Divergent Receptor Expression
We !rst explored altering the ligand binding a"nity as an engineering strategy to enhance its cell population speci!city. Here, we showed the binding pattern of monovalent ligands with various a"nities ranging from to to both receptors and (Fig. 2a). We found that when a target cell population expresses a receptor not expressed by o#-target cell populations, enhancing the a"nity of the drug to this receptor is a clear and e#ective strategy to increase selective binding to this population. For example, only signi!cantly expresses , while only signi!cantly expresses . When the a"nity to is enhanced and the a"nity to is reduced, the binding selectivity towards increases (Fig. 2b). The contour plots in Fig. 2a shows this trend more intuitively: when ligand a"nity to increases, which corresponds to shifting from subplots on the left to the ones on the right, binding is shown to vary more strongly according to expression, indicating an increase in the amount of ligand bound for the populations with high expression, such as , , and .
However, situations where both on-and o#-target cell populations express the same set of receptors and di#er only in their magnitude of expression are just as common. In these cases, we found out that it is bene!cial for the drug to bind tightly to the most comparatively highly expressed receptor on the target population. Some examples of this pattern are the selective binding to over , and over (Fig. 2a). For these two pairs, the bene!t of a"nity changes is limited by the relative discrepancy in receptor abundances (Fig. 2c,d). Since the ratios of receptor expressions between and are signi!cantly lower than those between and , the greatest binding selectivity that can be achieved for these two populations is also signi!cantly lower. When both receptors are uniformally more abundant in a target population than the o# target population, such as when comparing to , a"nity tuning fails to modulate binding selectivity (Fig. 2d). Therefore, to use a"nity changes for selectivity enhancement, it is critical to identify which receptors a cell population of interest expresses the most uniquely as compared to o#target populations.
Our binding contour plots are also useful for considering the interplay between a"nity modulation and intrapopulation receptor expression heterogeneity. For example, has relatively high variance in both and expression (Fig. 2a). When the binding a"nities to and are divergent, as shown in the subplots in the upper left corner and bottom right corner of Fig. 2a, the ellipse representing its range of receptor abundances rides through multiple contour lines, indicating wide variation in the quantity of bound ligand. This intrapopulation variation in the amount of ligand bound, however, is not observed when the a"nities to and are more balanced. Thus, a population's receptor expression heterogeneity and variation are important factors to consider when modulating ligand a"nity.

Valency Enables Selectivity Based on Quantitative Di"erences in Receptor Abundance
Given the limitations of deriving selectivity from a"nity changes, we next explored the e#ect of valency (Fig. 3). Our previous work has shown that this binding model can accurately predict the multivalent binding between IgG antibodies and Fcγ receptors. 17 To further demonstrate our model's capacity to predict multivalent binding activity, we !t our model to a set of published experimental measurements wherein $uorescently labeled nanorings were assembled with speci!c numbers of binding units. 18 After !tting to determine the crosslinking coe"cient of multivalent nanorings ( ) and receptor a"nity values ( ), our model was able to accurately match the binding of nanorings with four or eight subunits to cells expressing a known abundance of receptor partners (Fig. S1a).
Multivalent ligand binding di#ers from monovalent binding in its nonlinear relationship with receptor density, allowing multivalent ligands to potentially selectively target cells based on receptor abundance. 19 The e#ects of varying ligand valency are visualized in Fig. 3a. Here, the ligand only binds to receptor , not receptor , and the valency of the ligand is varied across subplots. Varying ligand valency results in a nonlinear relationship between receptor abundance and binding magnitude. For example, in the monovalent case, there are roughly the same amount of contour lines between and as there are between and . However, in the tetravalent case, there are comparatively more contour lines between and , indicating that ligand bound is a nonlinear function of receptor expression when considering multivalent binding (Fig. 3a).
Selectivity derived by changing valency requires coordinate changes in a"nity. As many previous studies have suggested, multi-valent ligands can demonstrate particular selectivity to target populations with high receptor abundances when their subunits have low a"nity to the receptors. 13,18 Comparing binding between cell populations exposed to ligands of variable a"nities over a range of ligand valency demonstrates the interplay between these factors ( Fig. 3b-d). Here, the ligand binding ratios between cell populations are shown for ligands of high, medium, and low a"nities for ( of , , and ), and the shaded areas indicate the variance caused by intrapopulation heterogeneity, determined by sigma point !lters. The ligand binding ratio between and is maximized by low a"nity ligands, but requires greater valency to achieve peak binding selectivity when compared to ligands of greater a"nity ( Fig. 3b). A similar valency optimum for a given a"nity is seen for binding selectivity between and , and and ( Fig. 3c-d). Ligands with lower a"nities achieve optimal binding with higher valencies and exhibit higher selectivity for cells expressing greater amounts of receptor.
For any speci!c multivalent ligand at equilibrium, it may bind to any number of receptors up to its valency, which we will refer to as its binding degree. Our model allowed us to identify the mechanism of valency-mediated receptor abundance selectivity by examining the distribution of ligands bound at each degree to di#erent cells for octavalent ligands (valency = 8). A cell expressing receptors displays similar amounts of binding at each degree for ligands with dissociation constants of , , and (Fig. 3e). However, a cell expressing 10% as many receptors exhibits extremely low amounts of higher-degree binding for ligand complexes of low binding a"nity ( Fig. 3f). This e#ect arises due to the ligand's rate of dissociation being larger than its multivalent binding association rate at low receptor abundances (Fig. 3g). Multivalent ligands undergo initial binding events at rates una#ected by steric e#ects and changes in local concentration, but low a"nity and receptor density severely limit the possibility of secondary binding events. In contrast, cells with a higher abundance of receptors accumulate ligand bound by having them bind multivalently, as the forward rate of 10 4 1000 100 10 nM secondary binding events is greater than that of receptor-ligand disassociation (Fig. 3f). This e#ect allows multivalent ligands to achieve selective binding to cells based on their receptor abundances, which monovalent ligands, even with engineered a"nities, are unable to do. to . e-f) Number of ligand bound to each possible number of receptors for cells exposed to octavalent ligand complexes composed of subunits with dissociation constants of , , or for receptor 1. e) Number of octavalent complex bound at each degree for a cell with receptors. f) Number of octavalent complex bound at each degree for a cell with receptors. g) Ratio of forward to reverse binding rate for secondary binding events for multivalent ligands to cells expressing variable amounts of receptors.

Non-Overlapping Ligand Targeting Drives a Limited Selectivity Enhancement of Mixtures
While most therapies rely on the action of a single molecular species, mixtures may enhance selectivity through combinations of actions. 20 Furthermore, some biologics inevitably act as mixtures of species through heterogeneity in glycosylation and the presence of endogenous ligands. 21 Therefore, it is important to understand how mixtures of complexes in$uence overall response.
To evaluate the contribution of mixtures, we evaluated model-predicted binding while varying the composition between two distinct monovalent ligands, each exhibiting preferential binding to either or (Fig. 4a). The trend that arises is very similar to an additive combination of the single ligand cases. This pattern highlights a key limitation of using mixtures for selectivity: selectivity between two populations varies monotonically with the composition, so any mixture combination is no better than using the more speci!c ligand entirely (Fig. 4b).
While mixture engineering fails to enhance binding selectivity between two cell populations, it is potentially bene!cial when considering two or more o#-target cell populations. More speci!cally, when a target population expresses two target receptors, but o#-target populations express each receptor individually in high amounts, drug mixtures can o#er enhanced selectivity. For example, when maximizing targeting to over and , which individually express high levels of the receptors found on the , a uniform mixture of ligands with high a"nity for receptor 1 and 2 provides a modest improvement in targeting selectivity (Fig. 4c). However, even in these cases, the magnitude of selectivity enhancement is modest. Finally, although we only consider the amount of binding, ligands can have non-overlapping signaling e#ects even with identical amounts of binding. In these cases, the e#ect of combinations can be distinct from either individual ligand. 17

Heterovalent Bispeci#c Ligands Exhibit Unique Charateristics When Activated Fully Bound
Constructing multispeci!c drugs has become a promising new strategy for !ner target cell speci!city with the advancement of engineering techniques. 23 However, the number of possible con!gurations of multispeci!c drugs is combinatorially large and impossible to enumerate. Here, we use heterovalent bispeci!c ligands as examples to explore the unique bene!t of multispeci!city distinct from any strategy analyzed before. Here, we compare bispeci!city with a 50%-50% mixture of two monovalent ligands, and a 50%-50% mixture of two di#erent homogeneous bivalent ligands (Fig. 5i). These two strategies both have some similarities to bispeci!c therapeutics. A bispeci!c ligand contains two di#erent ligand monomers with equal proportion, similar to a 50%-50% monovalent mixture. By comparing bispeci!c with a 50%-50% monovalent mixture, we can elucidate the extra bene!t of tethering these two monomers into one complex. A bispeci!c ligand is also by nature bivalent, and by comparing it with homogeneous bivalent drugs, we can understand how having two di#erent subunits in the same complex modify the behavior of a drug. Together, we sought to identify the unique advantages a#orded by heterovalent ligands.
We !rst applied the binding model to predict the amount of ligand bound in bispeci!c drugs (Fig. 5a), a 50%-50% mixture of two monomers (Fig. 5b), and a 50%-50% mixture of two di#erent homogeneous bivalents (Fig. 5c), with the same set of parameters in ligand concentration and a"nities. Surprisingly, the patterns of ligand binding in these three cases are almost exactly the same, and bispeci!city appeared to o#er no unique properties. However, many bispeci!cs only impart their therapeutic action when both of their subunits bind to a target population. For example, in the design of bispeci!c antibodies, it is common to require both subunits bound to be e#ective. 24,25 When the !rst antigenbinding region of a bispeci!c antibody docks at an antigen, it will only make a loose connection between the target and the antibody, unable to induce a strong immune response, and only when both antigen-binding regions bind to their target will this bispeci!c antibody be fully e#ective. We investigated whether bispeci!c ligands that have to be activated with both of their subunits binding to a target display any special cell population selectivity characteristics.
To scrutinize the case for bispeci!c ligands where the target binding of both subunits is required, we extended our model to calculate the amount of ligand fully bound (Fig. 5i). For the heterovalent bispeci!c case, ligand fully bound will only account for the amount of ligand with both of their ligand monomers bound to a target. With the same set of parameters, the predictions made for bispeci!c fully bound show a very distinct pattern from general ligand bound (Fig. 5e). The contour plot of fully bound bispeci!c ligands has more convex contour lines. For example, has about the same level of general ligand bound as (Fig. 5a), but it has signi!cantly more ligand fully bound than (Fig. 5e). This convexity of contour lines indicates that for bispeci!c complexes, doublepositive cells bind more ligands fully than cells only strongly expressing one type of receptors. This should be obvious since these two subunits of the bispeci!c complex prefer di#erent receptors. The results indicate that bispeci!c ligands will only exhibit special characteristics when it requires both receptor binding subunits to bind a receptor.
The speci!c amount of fully bound ligands is dependent on the ligand's propensity for crosslinking captured by the constant (Fig. 5i). Typically in biochemistry, when there is less steric hindrance among the subunits of a multivalent drug molecule (e.g. longer tether and smaller subunit size), or when there is local receptor clustering on the target cell, secondary binding will be easier to achieve, corresponding to larger . 26 We plotted the pattern of bispeci!c fully bound with , , and ( Fig.5d-f). In general, when is larger, the ligands are more capable of multimerization, and there will be more fully bound units. To demonstrate how this characteristic of fully bound bispeci!c ligand imparts cell population speci!city, we compared the selectivities between some chosen target-to-o#-target cell population pairs under bispeci!c drug versus another drug given a range of (Fig. 5g,h). These plotted numbers are the selectivities imparted by a bispeci!c drug divided by the selectivities from a drug mixture of either monovalent (Fig. 5g) or homogeneous bivalent (Fig. 5h). When these quotients are larger, it implies that fully bound bispeci!c ligand has greater advantages than its counterpart to impart selective binding for the chosen target cell Figure 5g compares the selectivities under bispeci!c ligands versus a 50%-50% mixture of monovalent ligands. The results show that fully bound bispeci!c can grant better selective binding when is small enough. This is logical, since when is small and crosslinking is rarer, most ligands will bind monovalently, and fully bound bispeci!c will especially favor the cell populations with higher receptor expression. However, when we compare fully bound bispeci!c to fully bound homogeneous bivalent mixtures (Fig. 5h), the advantage of bispeci!c drugs is not obvious except for to selectivity. Given that we only account for fully bound ligand for both therapeutics, the e#ect demonstrated in Figure 5g no longer holds. Together, we showed that bispeci!c ligands only exhibit unique advantages in inducing selective binding when they are only e#ective when both of their subunits bind and crosslinking is more di"cult.  . g,h) Comparing bispeci!c selectivities with mixture selectivities, varies with , the crosslinking constant. When the ratios are larger, bispeci!c ligands bind to target populations more speci!cally. g) bispeci!c selectivities divided by monovalent 50%-50% mixture selectivities, h) bispeci!c selectivities divided by a 50%-50% homo-bivalent mixture selectivities. i) Cartoons of bispeci!c ligands and ligand fully bound binding model.

Combining Strategies For Superior Selectivity
Each strategy described above provided selectivity bene!ts in distinct situations, suggesting that they might synergistically improve selectivity when combined. We explored this potential synergy using optimization to determine the ligand speci!cations which provided optimal selectivity for our theoretical cell populations. More speci!cally, we optimized the selectivity of a ligand for a particular population, while considering all other populations to be o#-target. Our optimization allowed ligand characteristics to vary within biologically plausible bounds; which characteristics we allowed to vary de!ned our examination of the e"cacy of our various ligand engineering strategies. We examined optimizing a"nity alone, mixture along with a"nity, and valency along with a"nity, and !nally combined all three strategies (Fig. 6).
Optimizing a ligand for selectivity to highlights a situation in which a"nity imparts greater speci!city, and optimal selectivity is achieved by combining a"nity and valency modulations (Fig. 6a-f). Here selectivity is optimized by ligands with selective binding to receptor 2, and those with valencies allowing for selective binding to cells with higher abundances of receptor expression. One case contradictory to this trend is shown during the optimization for selectivity towards the ( Fig. 6gl). While a"nity engineering is unable to impart some small contributions to enhanced selectivity, signi!cant improvement is only achieved when utilizing valency modulation techniques. A more di"cult design problem is featured in the optimization of ( Fig. 4m-r). Since it lies in the midst of the other populations in receptor expression space, any modulation of a"nity, valency or combining it with mixture-based strategies seems ine#ective. It should be noted that in all described cases, engineering the composition of a mixture of ligands is generaly ine#ective for imparting selectivity when the ligand's design speci!cations are $exible, and is likely only e"cacious when using ligands with static properties and considering multiple o#-target populations.
Our results highlight that both in singular and combined strategies for therapeutic manipulation, the target and o#-target populations dictate the optimal approach. It is also clear that combined approaches do o#er synergies which can be harnessed, but that those are only emergent in particular therapeutic situations. 25 Ligand concentration . Xo ligands are monovalent ligands with a"nities of for both receptor 1 and 2. The dissociation constant was allowed to vary between and for both receptors using the "a"nity" approach. Valency was allowed to vary from 1 to 16 for the "valency" approach in addition to a"nities varying. Mixtures were assumed to be composed of two monovalent ligands, and a"nities were allowed to vary in the "mixture" approach. The combined "all" approach allowed all of these quantities to vary simultaneously. The crosslinking constant was allowed to vary between and for all approaches. b-f,h-l,n-f) Heatmap of magnitude of ligand bound for ligand with optimized characteristics according to various ligand engineering strategies. Target population is shown in red. a-f) Pertains to optimal targeting of , g-l) pertains to optimal targeting of , and m-r) pertains to optimal targeting of .

Using Binding Competition to Invert Receptor Targeting
While the strategies above provided selectivity in many cases, we recognized that they are all limited to a positive relationship between receptor abundance and binding. Therefore, we wondered if binding competition with a receptor antagonist, or "dead ligand", could invert this relationship.
To investigate the e#ect of ligand competition with an antagonist, we modeled mixtures of ligands, but only quanti!ed the amount of binding for the active ligand. Here we chose to only consider a monovalent agonist and tetravalent antagonist (Fig. 7). We found that combinations of monovalent agonistic ligands and multivalent antagonistic ligands were able to uniquely target populations expressing small or intermediate amounts of receptors, which is demonstrated when comparing ligand binding ratios between to (Fig. 7a). Here, a nearly tenfold increase in selectivity can be granted to monovalent agonists when combined with a tetravalent antagonist. In this case, there are greater quantities of agonist bound to than (Fig. 7b). This is striking, as expresses either as many or fewer abundances of receptors one and two when compared to . This phenomenon, which could not be achieved without multivalent antagonists, occurs due to the preferential binding of multivalent antagonists to populations expressing higher abundances of (Fig. 7c, 3e). Thus, in cases where previously discussed ligand engineering strategies and approaches fail to achieve selective binding to cells expressing smaller or similar amounts of receptors to o#-target populations, combinations of agonistic and antagonistic ligands may provide unique bene!ts.

Discussion
Here, we developed and employed a multivalent, multi-ligand, multi-receptor binding model and used it to explore the e#ectiveness of various ligand engineering strategies for population-selective binding (Fig. 1). Using a representative set of theoretical cell populations de!ned by their distinct expression of two receptors, we examined the e"cacy of several potential ligand engineering strategies, including changes to a"nity, ligand valency, mixtures of species, multi-speci!city, antagonist competition, and these in combination.
Each strategy's contribution can be summarized by general patterns. We found that a"nity changes were most e#ective when the target and o#-target populations expressed distinct combinations of receptors (Fig. 2). Binding selectivity was enhanced by increasing the a"nity of the ligand for those receptors that are more abundantly expressed on the target cells and decreasing its a"nity for those that are not. Selectivity depended upon large di#erences in receptor expression levels between onand o#-target populations. When target and o#-target populations expressed the same pattern of receptors and only di#ered in receptor abundances, modulations in valency, but not a"nity, were e#ective (Fig. 3). A key determinant of valency's e#ectiveness was the secondary binding and unbinding rate which is dependent on both the kinetics of the receptor-ligand interaction and receptor abundance. Ligand mixtures were mostly ine#ective for imparting binding selectivity, and only had modest bene!ts when considering two or more o#-target populations (Fig. 4). Heterovalent bispeci!c ligands only showed unique advantages over mixtures of monovalent ligands or bivalent ligands when solely considering fully bound ligand (Fig. 5). These ligands exhibit preferential binding to target populations that have high expression of both receptors over those with high expression of only a single receptor, with the propensity for secondary binding acting as a key determinant for selectivity. We found that, while a single ligand engineering strategy dominated in its contributions to cell type selectivity, synergies between these strategies existed in some cases to derive even greater speci!city (Fig. 6). Finally, we found that combinations of monovalent therapeutic ligands with multivalent antagonistic ligands allow for the unique selection of target populations expressing relatively fewer receptors than o#-target populations (Fig. 7).
While our multivalent binding model identi!ed strategies for selective targeting in many cases, it also identi!ed situations for which selective binding is challenging. For example, selectively targeting populations based on their absence of receptor expression remains challenging. While we computationally show the potential of using multivalent antagonists with monovalent agonists to selectively target such populations, implementing this may be challenging in practice. In cases where a target population expresses fewer receptors of any kind than an o#-target population, our analysis suggests that targeting other receptors should be considered. However, in cases where target populations express more of any type of receptor than an o#-target population, we show that one or more of our formulated ligand engineering strategies can be employed to improve binding selectivity. While we expect the same patterns to apply with greater than two receptors, certain emergent behaviors may exist with tri-speci!c and more complex ligand binding.
A few of the strategies that we explored have been utilized in existing engineered therapies. For example, a"nity changes to the cytokine IL-2 have been used to bias its e#ects towards either e#ector or regulatory immune populations. 27,28 Varying the valency of tumor-targeted antibodies leads to selective cell clearance based upon the levels of expressed antigen. 19 Manipulating of the a"nities of the !bronectin domains on octovalent nanorings was shown to enhance the selectivity of binding to cancerous cells displaying relatively higher densities of !bronectin receptors compared to native tissue. 18 The tendency of low-a"nity, multivalent interactions to target cells expressing high receptor abundances was also described in a study describing the selectivity of multivalent antibody binding to tumor cells bound by a bispeci!c therapeutic ligand. 13 These examples lend support to the accuracy and translational value of our model. At the same time, recognizing these previously described ligand engineering approaches as separable strategies provides clearer guidance for future engineering.
Some of the optimization strategies described here have not been exploited in part due to the complexity of real biological applications. First, some strategies may not be biochemically practical. For example, the manipulating ligand a"nity requires intricate protein engineering. Potential dynamic changes in the receptor expression pro!le of a target population also complicate the matter. It is well documented that cancer cells can escape therapeutic targeting by upregulating 29,30 or downregulating 31 the expression of certain receptors. In this case, both the current and potential abundance of each receptor must be considered. While this work does not address many of these issues, we propose that using a computational binding model can analyze these strategies quantitatively and collectively from a mechanistic perspective. Even when the absolute mathematically optimal ligand characteristics cannot be achieved biochemically, our analyses provide guidance within the bounds of what is attainable and how to approach the optimum, accounting for implementation feasibility and facilitating the implementation of strategy combinations.
In many therapeutic applications where selective engagement of target cell populations is an important performance metric, such as the treatment of cancer, cellular therapies are becoming increasingly popular. 32 Human engineered chimeric antigen receptor (CAR) T cells have enhanced the potential to selectively recognize and attack malignant tissues. 33 These technologies bypass ligand-receptor binding restrictions by allowing recognition in signaling response. However, we have shown here that selectivity can often be attained with relatively simple therapeutic ligands. This study lays the framework for how ligand engineering can be directed using computational modeling. It should be noted that the application of this logic is reliant on knowledge of the target and o#-target cell population receptor expression levels. Future application of the ligand binding logic described in this study could be guided using high-throughput single-cell pro!ling techniques, such as RNA-seq or highparameter $ow cytometry. A computational tool that could directly translate such datasets into ligand design criteria, selecting among potential receptor targets, may represent a potential avenue for the translation of our analyses into a more broadly applicable ligand engineering tool.

Data and Software Availability
All analysis was implemented in Python v3.8, and can be found at https://github.com/meyer-lab/cellselective-ligands.

Generalized multi-ligand, multi-receptor multivalent binding model
To model multivalent ligand complexes, we extended our previous binding model to the multi-ligand case. 17 We de!ne as the number of distinct monomer ligands, the number of distinct receptors, and the association constant of monovalent binding between ligand and receptor as . Multivalent binding interactions after the initial interaction have an association constant of , proportional to their corresponding monovalent a"nity. The concentration of complexes is , and the complexes consist of random ligand monomer assortments according to their relative proportion. The number of ligand complexes in the solution is usually much greater than that of the receptors, so we assume binding does not deplete the ligand concentration. The proportion of ligand in all monomers is . By this setup, we know . is the total number of receptor expressed on the cell surface, and the number of unbound receptors on a cell at the equilibrium state during the ligand complex-receptor interaction.
The binding con!guration at the equilibrium state between an individual complex and a cell expressing various receptors can be described as a vector of length , where is the number of ligand bound to receptor , and is the number of unbound ligand on that complex in this con!guration. The sum of elements in is equal to , the e#ective avidity. For all in , let when is in , and . The relative amount of complexes in the con!guration described by at equilibrium is being the multinomial coe"cient. Then the total relative amount of bound receptor type at equilibrium is By conservation of mass, we know that for each receptor type , while is a function of . Therefore, each can be solved numerically using . Similarly, the total relative amount of complexes bound to at least one receptor on the cell is

Generalized multivalent binding model with de#ned complexes
When complexes are engineered and ligands are not randomly sorted into multivalent complexes, such as with the Fabs of bispeci!c antibodies, the proportions of each kind of complex become exogenous variables and are no longer decided by the monomer composition 's. The monomer composition of a ligand complex can be represented by a vector , where each is the number of monomer ligand type on that complex. Let be the proportion of the complexes in all ligand complexes, and be the set of all possible 's. We have .
The binding between a ligand complex and a cell expressing several types of receptors can still be represented by a series of . The relationship between 's and is given by . Let the vector , and the corresponding of a binding con!guration be . For all in , we de!ne where and . The relative amount of complexes bound to a cell with con!guration at equilibrium is v q,eq = ( ) x R tot,n = R eq,n + R bound,n n R bound,n R eq,n R eq,n R tot,n i {1, 2, . . . , N L } ψ ij = R eq,j K a,ij K * x j = {1, 2, . . . , N R } ψ i0 = 1 q Then we can calculate the relative amount of bound receptor as By , we can solve numerically for each type of receptor. The total relative amount of ligand binding at equilibrium is

Mathematical optimization
We used the SciPy function scipy.optimize.minimize to combine several strategies and achieve optimal selectivity. 34 Unless speci!ed otherwise, the initial values for optimization were for crosslinking coe"cient , 1 for valency , 100% ligand 1 for mixture composition, and for the a"nity dissociation constants. The boundaries were to for , 1-16 for , 0-100% ligand 1 for mixture composition, and to for the a"nity dissociation constants.

Sigma point #lter
To consider the intrapopulation variance of a cell population in the optimization, we implemented the sigma point !lters, 16 a computationally e"cient method to approximate the variance proprogated through an ordinary di#erential equation-based model.
Binding activity of nanorings displaying both high (C5) and low (B22) a"nity EpCAM binding domains was measured. Binding to an EpCAM high ( antigens/cell) population was measured using $ow cytometry. We used nonlinear least squares optimization, as described above, to !t our multivalent binding model to the binding data, using a unique scaling factor for each !bronectin clone to convert between measured $uorescent intensity and magnitude of ligand bound. We allowed a"nity of the !bronectin clones to vary during optimization.
R tot,n = R eq,n + R bound,n R eq,n 1 µM 10 −15 10 −9 K * four or eight !bronectin binding domains at concentrations of . C5 and B22 are high and low a"nity !bronectin binding domains respectively. Using nonlinear least-squares regression, the crosslinking coe"cient was found to be ; the dissociation constants for C5 and B22 were found to be and respectively, and ligand to $uorescent conversion factor for C5 and B22 were found to be and respectively. b) Predicted $uorescent values for MDA cells expressing , SK cells expressing , LNCaP cells expressing , and MCF-7 cells expressing receptors per cell as described in Csizmar et al. 18 . Fluorescence was predicted for interactions with octa-and tetravalent ligands expressing either C5 or B22 !bronecting binding domains. c) Predicted binding ratios of MCF-7 cells to MDA cells when bound with octa-and tetravalent ligands expressing either C5 or B22 !bronecting binding domains.