Density-dependent cooperative non-specific binding in solid-phase SELEX affinity selection

The non-specific binding of undesired ligands to a target is the primary factor limiting the enrichment of tight-binding ligands in affinity selection. Solution-phase non-specific affinity is determined by the free-energy of ligand binding to a single target. However, the solid-phase affinity might be higher if a ligand bound concurrently to multiple adjacent immobilized targets in a cooperative manner. Cooperativity could emerge in this case as a simple consequence of the relationship between the free energy of binding, localization entropy and the spatial distribution of the immobilized targets. We tested this hypothesis using a SELEX experimental design and found that non-specific RNA aptamer ligands can concurrently bind up to four bead-immobilized peptide targets, and that this can increase their effective binding affinity by two orders-of-magnitude. Binding curves were quantitatively explained by a new statistical mechanical model of density-dependent cooperative binding, which relates cooperative binding to both the target concentration and the target surface density on the immobilizing substrate. Target immobilization plays a key role in SELEX and other ligand enrichment methods, particularly in new multiplexed microfluidic purification devices, and these results have strong implications for optimizing their performance.

1 Nonlinear regression is equivalent to linear regression using linear approximations to the model function at the data points. Accuracy is improved by using parameters that improve the linear approximation, i.e., parameters that reduce the contributions from the second derivatives. When K d is used, this contribution is [(∆K d ) 2 /2]|d 2 B K d /dK 2 d | = (∆K d ) 2 /(K d + T ) 3 . When logK d is used, this is [(∆logK d ) 2 /2]|d 2 B K d /d(logK d ) 2 | = (∆K d ) 2 (1−T /K d )/[2(K d +T ) 3 ]. The ratio of the logK d to the K d second derivative contribution is |1−T /K d |/2, which vanishes at T = K d and is less than one over most of the relevant range of concentrations.
Σ of the {φ i } estimates as

B. REDUCED RETENTION CAUSED BY WASHING
Aptamer-target complexes may dissociate during washing. We show here that the only effect of this non-equilibirum process on the retention curves is to reduce r.
We can evaluate its effect on non-cooperative binding using the dynamic model where B (k + ,k − ) (T,t) is the concentration of bound ligands at time t, L 0 is the total ligand concentration, k + and k − are the association and dissociation rates, t b is the duration of the binding reaction, and Θ(·) is the Heaviside step function. We assume, as was the case for this study, that the experiments are performed in great target excess so T , the unbound target concentration, is essentially constant. Also, because the wash volumes were large, rebinding during the wash period is ignored. Solving Eq. S1 for t > t b determines the solution binding at time t = t b +t w , where t w is the total wash duration, to be is the equilibrium binding fraction given by Eq. 1. As long as binding phase has progressed to equilibrium [i.e., t b (k − +k + T ) −1 ], the only non-equilibrium effect will be that of washing, and this will simply introduce a T -independent factor e −k − tw < 1 that will reduce r. 2 To analyze solid-phase binding we refer to Eq. 5 and assume, as seems likely, that the time-scale of the U B 1 transition is much slower than that of the subsequent B n B n+1 transitions. Hence, we can use an adiabatic approximation wherein the different bound states are assumed to be in equilibrium with each other so that B(t) ≡ n B n (t) ∝ B 1 (t). Therefore, Eq. S2 is still applicable with the replacement of k + and k − by the TSD-dependent cooperative rate constants. Thus, in this case as well, washing only reduces r and does not change the shape of the binding curve.

C. ASYMPTOTIC DDC MODEL
To set the stage we first analyze non-cooperative solution binding using partition functions in a manner that generalizes to the cooperative immobilized case. headingSolution binding We assume that a ligand and a target must be located within the same "binding volume" v to bind, that the total volume V 0 equals the "reference volume" (1 liter for molar concentration units), and that the concentration of unbound targets is T . Since the ligands do not interact, we can consider one individually. It can be located in one of the binding volumes, which is the number of its (discrete) spatial configuration states. Including the correct Boltzmann counting factor, the targets can be in Since the complete configuration space of unbound states is the outer product of the ligand and target spaces and there is an unbound microstate corresponding to each configuration, there are Ω L V ×Ω T V unbound microstates. Without loss of generality we can assign each of them zero free-energy, so the partition function of the unbound states is The partition function of the bound states also factors into a product of ligand-and target-dependent terms, but the latter is now more complicated. As before, the ligand can be in any one of the Ω L V binding volumes. However, only those target configurations having at least m ≥ 1 targets in the binding volume containing the ligand can contribute to the Boltzmann sum over the bound microstates, and the contribution will depend on m. Therefore, we begin by counting the number of ways of having m randomly distributed targets in this binding volume. It is easy to do this by working "backwards" from the corresponding probability: This is the Poisson distribution P(m;µ V ), where µ V = T N a v is the mean number of targets per binding volume. 3 Since this probability is the ratio between the number of target configurations satisfying the criterion, (Ω T V ) m , and the total number of target configurations, Ω T V , Thus, the number of configuration states having m targets colocalized with the ligand is Ω L V ×(Ω T V ) m . For each such configuration state, there are ways of selecting one of the colocalized targets for binding to the ligand. Therefore, there are Ω L V ×(Ω T V ) m ×ω B m bound microstates having m colocalized targets. Summing over m, the total number of bound microstates is All of the bound microstates have the same "intrinsic" binding free-energy for a colocalized ligand and target, ∆G . 4 Therefore, denoting temperature as τ , Boltzmann's constant as k B , and β B ≡ (k B τ ) −1 , the bound state partition function is and the bound/unbound ligand ratio is and ∆S v is the entropy of localization, which corresponds to the v/V 0 prefactor in Eq. S4. This proves Eq. 4 and, in combination with Eq. 3, proves Eq. 1.

Immobilized binding
Denoting the bead 5 concentration (expressed as the fraction of the reaction volume occupied by packed beads) as s and the number of moles of binding volumes on the surface of a liter of packed beads as C s , the ligand has Ω L V possible configuration states in solution and states on the bead surface. The targets are immobilized to the bead surface in one of possible configuration states. Therefore, ignoring (vanishingly small) terms of O(Ω L σ /Ω L V ), the partition function of the unbound states is We calculate the number of ways of having m randomly distributed targets located within a specific binding volume on the bead surface by again working backwards from the probability of this event assuming random placement of the targets on the bead surface. This is P(m;µ σ ), where is the mean number of targets per binding volume on the bead surface. Therefore, the number of ways of having m randomly distributed targets in a specified binding volume is and the number of configuration states having m targets colocalized with the ligand is Ω L σ ×(Ω T σ ) m . For each of these configuration states, there are ways of selecting n of the colocalized targets for binding. Thus, there are Ω L σ ×(Ω T σ ) m ×ω B mn microstates having m ligand-colocalized targets and n ligand-bound targets. Summing over m, n ≤ m, the total number of microstates having n targets bound to a ligand is Each of these bound microstates has residual binding free energy ∆G 1 +(n−1)∆G + , where ∆G 1 and ∆G + are 5 The analysis generalizes to any immobilizing substrate. the residual free-energies for binding the first and subsequent colocalized targets, respectively. Therefore, the partition function for states having n ≥ 1 bound targets is Summing over n, 1 ≤ n, gives the complete bound state partition function where we used Eqs. S3, S6, S8, and S9 and the Poisson distribution's summation property to get the last line. Therefore, thereby proving Eq. 6. Comparing this with Eq. S5, we see that K 0 d would be the same as the solution binding dissociation constant if ∆G 1 = ∆G, i.e., if a ligand were to bind colocalized solvated and immobilized targets with the same affinity.κ, the asymptotic cooperativity constant, is determined by ∆G + . Z B n is proportional to a Poisson distribution with expectation valueκT /s and the probabilities of binding one or two targets are the same whenκT /s = 2. Thus,κ determines the scale of TSDs at which cooperativity becomes significant.

Relationship to intrinsic binding energy
The decomposition of the binding free-energy in terms of ∆G 1 , ∆G + , and the translational entropy component T ∆S v is similar in principle to Jencks' (2) decomposition of binding free-energy into spatial (translational plus rotational) entropy, strain, and the intrinsic components. However, the decomposition used here, which is optimized for the analysis of the sequential binding process of Eq. 5, differs in detail. While ∆G 1 excludes the change in translational entropy, it includes the change in rotational entropy and ∆G + includes any contribution from strain that may be present. However, we expect strain contributions to be small, so that ∆G + is approximately equal to the intrinsic binding free-energy for non-specific binding of the ligand once the first binding has induced any global conformational changes that are needed for non-specific binding.

Number of bound targets
The mean number of targets bound by a ligand conditional on the ligand being bound is the Boltzmann-weighted average of n over the bound states. Noting that −∂Z B n /∂(β B ∆G + ) = (n−1)Z B n and ∂/∂(β B ∆G + ) = −∂/∂ logκ and using Eq. S11 we get

Reparameterization for VTSD experiments
In a VTSD experiment s = s * , and it is convenient to define γ ≡κ K 0 d s * and to rewrite Eq. 6 as ρ VTSD where the dependence of γ on s * is implicit. While cooperativity will increase binding whenever γ > 0, the retention curve will have an inflection point only when γ > 2. This corresponds toκK 0 d /s * > 2, which means that there will be an inflection point in a VTSD experiment when the probability of binding one or two targets is the same when T = K 0 d . The inflection point will be at T infl = (K 0 d /γ)log(γ −1). It is sometimes useful to express Eq. 6 in terms of K 0.5 (K 0 d ,γ), the half-binding dissociation concentration at which andn(κ), the mean number of targets bound to a bound ligand when T = K 0.5 (K 0 d ,γ): n(γ) and K 0.5 (K 0 d ,γ) are often determined more accurately by VTSD data than K 0 d and γ because the first two parameters primarily depend on the (often most accurate) data in the mid-range of binding. The asymptotic DDC equation can reparameterized in terms of these variables using the inverse relations where W (x) is the principal branch of the Lambert W function. 6

D. COMPLETE DDC MODEL
In the complete model we assume that the ligand has N non-overlapping, non-specific binding sites, where N is determined by the size of the target in relation to the ligand and possibly other factors. (While this is a simplificationthe positions of potential non-specific binding may not be fixed-we expect that the error it introduces will be small.) The analysis in this case follows that of the asymptotic case. Ω L σ (Eq. S6) and Ω T σ (Eq. S7) are the same as in the asymptotic case, but (Ω T σ ) m and its dependent quantities are different. We calculate the number of ways that m targets can be colocalized with these binding sites by again working backwards from the probability of the event. Since there are µ σ targets per binding volume and N binding sites per binding volume, the number of targets per binding site is µ σ /N and the probability of having m targets colocalized with the ligand binding sites is Bin(m;µ σ /N,N ), where is the binomial distribution. Therefore, the number of target configurations having m targets colocalized with binding sites of a ligand at a specific location is (Ω T σ ) m = Ω T σ Bin(m;µ σ /N,N ) . 6 We can also parameterize the binding curves using K 0 d and K0.5 by inverting these equations to get γ = −1− is the lower of the two real branches of the Lambert W -function over −exp(−1) < x < 0. ω B mn (Eq. S10) is unchanged from the asymptotic analysis.
Multiplying Ω L σ , (Ω T σ ) m , and ω B mn and summing over m, n ≤ m ≤ N , gives Multiplying this by the Boltzmann factor gives the partition function for states having n bound targets Summing over n, 1 ≤ n ≤ N , and using the binomial distribution's summation property gives Bin n; Using Eqs. S3, S6, and S8 and the definition of K 0 d we get thereby proving Eq. 8. In this case, Z B n is proportional to a binomial distribution with expectation value N κ(T /s)/(N − 1+κT /s). The probabilities of binding one or two targets are the same when κT /s = 2.
As before, we calculate the number of bound targets usinḡ n (K 0 d ,κ,N ) = ∂ logZ B /∂ logκ+1. This gives which, as expected, has the limits lim T /s→0n For analyzing VTSD data alone, it is convenient to define and to rewrite Eq. S13 as This has the same non-cooperative behavior for small κT /s as the asymptotic equation, but it has different behavior for large κT /s, which corresponds to the concurrent binding of N targets. The dissociation concentration is and the mean number of targets bound when T = K 0.5 (K 0 d ,γ,N ) is It is easily verified that all these equations go to their asymptotic counterparts as N → ∞.

E. APPLICATION OF THE DDC BINDING MODEL TO MOBILE TARGETS
Although the DDC binding equations were derived assuming that the targets are immobilized, the same analysis applies when they are not. They only require that the targets be restricted to a two-dimensional (or even one-dimensional) subspace. This follows because target mobility does not affect the counting of the target conformational states, so the Boltzmann sums are not affected: in both cases the sum over bound states contains the term Ω B n that includes the same Poisson (in the asymptotic case) or binomial (in the complete case) distribution, which represents the entropy associated with target localization. However, the physical interpretations are different: in immobilized binding this term comes from the probability that the ligand is located in a region having the specified number of targets; in binding to mobile targets that can move within the lower-dimensional subspace it comes from the entropy change occurring when the specified number of targets are moved into the binding volume of an arbitrarily placed ligand. The key factor in both cases is that the large reduction in entropy that occurs when the ligand is localized from the three-dimensional solution to the lower-dimensional subspace only occurs for the first, but not the subsequent, bindings.

F. COMBINED SPECIFIC AND COOPERATIVE NON-SPECIFIC BINDING
Here we consider the extent to which cooperative non-specific binding will affect the relative affinities of a low-affinity ligand that binds only non-specifically and a high-affinity ligand that possesses a specific binding site. This will depend on whether the conformational restrictions imposed by specific binding preclude or permit additional non-specific binding. If, as seems likely, it were precluded, then the dissociation constant of the high-affinity ligand will be unchanged but that of the low-affinity ligand would be decreased from K 0 d to K * d . Changes as large as those observed in the experiments reported here would decrease enrichment by two ordersof-magnitude. At the other extreme, if specific binding did not diminish additional non-specific cooperative binding, we show below that the DDC equation predicts that the effective affinity of the high-affinity ligand would be increased even more than that of the low-affinity ligand; this would cause up to an N -fold increase in the inverse ratio of the cooperative dissociation constants. Thus, in this case as well, cooperative binding will distort the observed (effective) affinity. This will be particularly problematic in a SELEX experiment if the cooperative effect differs between different specifically binding aptamers: because of this, the ranking of effective affinities, which determines which aptamers will be selected, will not match the ranking of the non-cooperative affinities, which determines the utility of the aptamers in an application.
We analyze the change in effective affinity occurring in the second case by modifying the procedure used in Sec. D. The calculation of Ω B n is unchanged. However, to calculate Z B n , the analog of Z B n , we must treat the configurations that include or do not include the specific binding site separately and account for the difference between the specific and nonspecific binding free-energies. Since there are N −1 n−1 ways of choosing the specific binding site and n−1 non-specific binding sites, N −1 n ways of choosing all n target-bound sites from the N −1 non-specific binding sites, and a total of where ∆G is the residual free-energy for specific binding to a target located in the binding region and Z B n is given by Eq. S12. In the first line, the first term represents states having targets bound to the specific site plus n−1 non-specific sites and the second term represents states having targets bound to n non-specific sites.
Summing over n, 1 ≤ n ≤ N , evaluating the first term using the formula for the mean of the binomial distribution and the second term as before, we get the high-affinity ligand boundstate partition function Dividing by Z U gives where K 0 d and κ are defined as before and is the dissociation constant in the absence of cooperative binding. 7 Table S1. DDC and Langmuir best-fit parameters to the retention data of Fig. 2.