Generic self-stabilization mechanism for biomolecular adhesions under load

Mechanical loading generally weakens adhesive structures and eventually leads to their rupture. However, biological systems can adapt to loads by strengthening adhesions, which is essential for maintaining the integrity of tissue and whole organisms. Inspired by cellular focal adhesions, we suggest here a generic, molecular mechanism that allows adhesion systems to harness applied loads for self-stabilization through adhesion growth. The mechanism is based on conformation changes of adhesion molecules that are dynamically exchanged with a reservoir. Tangential loading drives the occupation of some states out of equilibrium, which, for thermodynamic reasons, leads to association of further molecules with the cluster. Self-stabilization robustly increases adhesion lifetimes in broad parameter ranges. Unlike for catch-bonds, bond rupture rates can increase monotonically with force. The self-stabilization principle can be realized in many ways in complex adhesion-state networks; we show how it naturally occurs in cellular adhesions involving the adaptor proteins talin and vinculin.


I. INTRODUCTION
Multicellular organisms are held together by complex biomolecular adhesion structures. For decades, cellular adhesions have been extensively studied, both, as a paradigm for fundamental biophysical mechanisms, as well as to understand their essential biological function. Yet, a very fundamental property of cell-matrix adhesions remains mysterious -their ability to adapt their size to the mechanical load. How is it possible -from a physics perspective -that the strength of adhesions increases under load?
Cell-matrix adhesions, also called focal adhesions, are crucial for cell physiology [1,2], cell motility [3], cancer metastasis [4,5], and development [6,7]. The structures consist of transmembrane integrins and adaptor proteins that connect the force-generating actomyosin cytoskeleton with the extracellular matrix. Accordingly, focal adhesions have been likened to a "molecular clutch" [8][9][10][11][12][13]. In a pioneering experiment by Riveline et al. [14], it was shown that local application of centripetal forces to adherent cells induces focal adhesion growth [14]. Moreover, the size of focal adhesions is proportional to the load [15]. The regulatory network associated with focal adhesions is complex and several biological processes have been proposed to play a role for adhesion stability. These include mechanosensitive activation of integrins [11], catch-bond behavior of integrins and vinculinactin binding [16], non-linear mechanical response of unfolded proteins, and downstream signaling, e.g., mediated by the adaptor protein p130Cas [17,18]. During the past years, the pivotal role of the adaptor protein talin for adhesion maturation has also been established [19,20]. Talin directly transmits forces by binding with its glob-ular head domain to integrin, while its rod domain links to F-actin [21]. Under stretch, conformation changes in talin occur, leading to an unfolding of protein domains and to the exposure of cryptic binding sites for vinculin [22][23][24]. Vinculin, in turn, further recruits Factin and thereby strengthens the linkage [25][26][27]. Moreover, other adhesion types such as adherens junctions are also capable of a load adaptation based on unfolding and recruitment of further molecules [28]. In spite of the considerable amount of theoretical approaches [29][30][31][32][33][34][35][36][37] and pioneering work combining modeling and experiment [12,38,39], the physical principle underlying loadadaptation of focal adhesions remains largely not understood.
We propose here a generic mechanism through which molecular adhesions can harness mechanical load for adapting their size and stability without active feedback. A minimal model is employed, which relies on a combination of unfolding of adhesion molecules under force [39,40], the dynamical exchange of molecules with a reservoir [29,30], and possibly the recruitment of additional molecules that stabilize the unfolded conformations, all under the constraint of thermodynamical consistency. The key property of our molecular statenetworks is that some states are separated from the reservoir and driven out of equilibrium. Under tangential stress, the state occupations shift, leading to a growth of the adhesion under increasing load, until the adhesion ultimately fails at very high loads. This non-equilibrium mechanism of adaptation is simple and robust. We show how the mechanism is naturally realized by the talinvinculin system at focal adhesions and perhaps in other bioadhesions. Moreover, the proposed mechanism is general enough that it could possibly be implemented in en- and β + (h). Molecules can enter or leave the adhesion cluster with rate constants given by γ ± . b) Mean number of bonds in steady state as a function of loading force. Symbols show simulation results and lines correspond to approximate analytical results. In the basic model, an increased load on the adhesion reduces the number of bonds (red). In contrast, an increased load produces a growth of the number of bonds in the generalized model (blue). c) Illustration of the generalized adhesion model incorporating molecule unfolding and refolding with rates δ ± a,b (h) as well as a molecule-exchange with the reservoir. Mechanical load drives the system out of equilibrium, shifting the state occupations. See also SI movies 1,2. gineered biomimetic systems.

II. RESULTS
Molecular adhesion model. We consider a generic adhesion system consisting of N molecules that form harmonic bonds between two planar, parallel surfaces, see Fig. 1. The molecule extension, i.e., the difference between actual length and rest length, is denoted by h and the spring constant by κ. The bottom surface is fixed in space and a constant tangential loading force F is exerted on the top, leading to a time-dependent tangential shift s. The model is two-dimensional and forces normal to the surface are not considered. The adhesion-molecule number N can vary by exchange with a reservoir, as determined by the rate constants γ ± . Individual molecules undergo stochastic transitions between the adhesion cluster and the reservoir, and all transition rates are chosen to fulfill detailed balance when F = 0, i.e., in thermal equilibrium, see Methods IV, which avoids unphysical energy injection that could produce an apparent motor-like behavior.
In a first, basic model, see Fig. 1a, molecules from the reservoir associate reversibly with the adhesion system via the state a, in which they have not yet formed a bond between the upper and lower surfaces. The state where a bond has formed is denoted by b. The bond formation rate is β + (h), which is maximal when the extension of the adhesion molecule equals an optimal binding distance |h| = b . For bond formation, it is assumed that the extension h fluctuates thermally with magnitude σ = k B T /κ, where k B T is the thermal energy scale. For bond dissociation, we focus on slip-bond dynamics with rupture rates β − (h) that increase exponentially with bond extension, see Secs. IV A-B.
In a second, generalized model, see Fig. 1c, the molecules can undergo abrupt conformational changes to unfolded states, denoted by a subscript u, such that a u are the unfolded, unbound states and b u the unfolded, bound states. The overall number of bound molecules, irrespective of their conformational state, is N B = N b + N bu , and the corresponding overall number of unbound molecules is N A = N a + N au . Details of the reaction rates are given in Methods IV A-C. We assume that mechanical relaxations occur instantaneously and viscous damping is neglected, so that the sum of the forces borne by the bonds equals the applied load F at all times. Stochastic bond dynamics is simulated with an exact algorithm, see Supplementary Information II.
Self-stabilization of macromolecular adhesions.
We first consider the basic adhesion model, in which molecules do not change conformation, see Fig. 1a  a stationary state of persistent tangential sliding of the surfaces. For these adhesions, the mean number N B of bonds always decreases monotonically with F , see Fig. 1b. Therefore, increasing load on adhesions consisting of simple molecules promotes adhesion failure characterized by rupture of all bonds. Next, we consider the generalized adhesion model consisting of molecules that can undergo an unfolding transition under force, see Fig. 1c and SI movie 2. For simplicity, we assume that unfolding only entails an increase in the rest length, while the elastic properties remain unchanged. Remarkably, the mean number N B of bonds now initially grows with increasing load F , see Fig. 1b, which prevents early adhesion failure. This striking effect, which we call "self-stabilization", is the central finding of this work. We want to emphasize that this effect crucially depends on the exchange of "folded" molecules (state a) with the reservoir, so that the transitions a → b → b u lead to a depletion of a molecules, which can be replenished from the reservoir, thereby leading to an increase of the adhesion size. This mechanism is explained in detail below. The simulation results can be corroborated with a mean-field approximation, which can be solved analytically [33]. The stationary distributions of molecules with extension h in the bound and unbound states are denoted by n b (h) and n a (h), respectively. For the basic adhesion model without molecule unfolding, a drift-reaction equation is assumed where the average sliding velocity of the adhesion v = ṡ stretches the molecules, such that Only stationary solutions with ∂ t n b (h) = 0 are considered. The total number of molecules in the adhesion is obtained as where the extension of the unbound molecules obeys a Gaussian distribution, n a (h) ∝ N (0, σ 2 a ). The non-linear equations are solved by expanding the distributions for small absolute values ofṽ = v/(k β σ b ), where the binding constant k β is employed as time unit and the extension variance of bound molecules σ 2 b as length unit, e.g., . The asterisk ( * ) here and in the following denotes equilibrium quantities calculated with F = 0. Using the additional assumption that the optimal molecule extension for binding, b , is much smaller than the typical length fluctua- and because, due to symmetry under reversal of the force direction, to leading orderṽ ∝ F . For the general case, where b σ b does not hold, the first non-vanishing correction to the equilibrium solution for the bonds N B can be shown to still be of the orderṽ 2 and strictly negative, see Supplementary Information III. Thus, tangential force reduces the number of bonds and thereby always destabilizes simple adhesions consisting of molecules that do not undergo conformation changes, as expected intuitively.
To support the effect of self-stabilization in the generalized model shown in Fig. 1c with analytical theory, we supplement Eq. (1) by two additional equations for binding and unfolding transitions, see Supplementary Information IV. For |ṽ| 1, the overall number of bound molecules, N B = N b + N bu , can be expanded as The leading contribution is linear inṽ with a positive coefficient N B,1 (obtained from numerical analysis), as required for self-stabilization! The second-order coefficient N B,2 can be positive or negative, see Supplementary Information Fig. S3. Hence, analytical models confirm the existence of a self-stabilization regime where the number of adhesion bonds initially increases with load.

Mechanism of self-stabilization.
To obtain deeper insight into the necessary ingredients for selfstabilization, and to demonstrate the essential contribution of the reservoir, we compare several variants of the basic models, see Fig. 2 -such as a simple adhesion model with fixed molecule number (model I), a model comprising a molecule reservoir and therefore adhesions of variable size (model II), a model with a fixed number of molecules that can unfold under force (model III), and a model combining unfolding molecules with a molecule reservoir (model IV). Some results from models II and IV are also shown in Fig. 1. For models I-III, the mean number of bonds N B decreases with force, see Fig. 2b. For model II with variable system size, even the number of molecules in the adhesion decreases with force, leading to an earlier adhesion failure on average. Thus, neither a variable adhesion size nor molecule unfolding alone result in self-stabilization. In model IV, which combines a variable adhesion size with molecule unfolding, both the total number of molecules N and the number of bonds N B initially grow with increasing load on the adhesion, Figs. 2a,b. The increased number of bonds improves load sharing among the molecules. One consequence is a significant reduction of the sliding motion of the adhesion, Fig. 2c. The inset in Fig. 2b shows that the bound fractions of molecules, N B /N , as a function of F collapse onto a single master curve for all models (I-IV). Hence, self-stabilization results from force-induced growth of the adhesion and not from changes of the rupture properties of individual molecules. This is the key difference to established catch-bond models, where individual molecules exhibit an increase of bond lifetime within limited force regimes. Figure 2d illustrates the underlying mechanism of selfstabilization. The tangential load F causes a continuous molecular-state turnover with recurring stretch, unfolding, and rupture of molecules along the transitions a → b → b u → a u . Overall, the load increases the occupation of state a u . Meanwhile, the state a, representing unbound, folded molecules, is in contact with the reservoir and molecules are replenished here, which allows a concurrent increase of the overall molecule number. The recruitment of molecules from the reservoir crucially depends on the intermediate state a u not to be equilibrated, see Suppl. Figs. S6,S8, and to rebind and contribute to load sharing again via a u → b u . It is important to note that without the reservoir, self-stabilization cannot occur because it requires growth of the adhesion cluster. Generically, the principle behind self-stabilization is that molecular-state occupation statistics are driven out of equilibrium in a way that results in further influx of molecules from a reservoir. Self-stabilization is a robust mechanism. To be an efficient and universally applicable mechanism, selfstabilization must be able to compensate load changes in non-stationary conditions and should not depend on a fine-tuning of parameters. To investigate these aspects, we consider step-like load changes for different model parameters. After a load jump from F 1 = 0 to F 2 = F , adhesion clusters either dissociate quickly or reach a new non-equilibrium steady state. Self-stabilization after a load jump depends on the strength of the molecule exchange with the reservoir, which is controlled by the values of γ + and γ − , Fig. 3a. Exemplary trajectories for the number of bonds in the folded and unfolded states, b and b u , are shown in Fig. 3b. The most-likely rupture forces are higher in self-stabilizing than in non-self-stabilizing adhesions, and grow with increasing reservoir-exchange rates (γ ± > k β ), see Supplementary Fig. S4, because adhesion molecules can quickly be recruited from or released into the reservoir. The self-stabilization mechanism can thus also work under dynamic load conditions.
To measure adhesion lifetimes, we simulate systems consisting of few molecules, because lifetimes increase rapidly with the number of molecules. Load-jump simulations are carried out for a reservoir-exchange rate ratio γ = 1, which leads to adhesions with N * ≈ 10 molecules in equilibrium. Lifetime is measured as the time from the force jump, F 1 → F 2 , to the rupture of the last adhesion bond. While lifetimes of adhesions with no reservoir coupling decrease monotonically with force, i.e., show a pure slip-bond behavior, lifetimes of self-stabilizing adhesions exhibit a maximum at non-vanishing forces, see Fig. 3c. This lifetime maximum becomes more pronounced for increasing γ + and also depends on the rate-constant ratio k δ /k β , see Supplementary  rameter choices, we vary the unfolding length ∆, which determines the width of the energy barrier between the native and the unfolded molecule state. For adhesions without reservoir coupling, γ ± = 0, the unfolding length ∆ does not have a large impact on the adhesion lifetime, see Fig. 3d(i). However, a significant increase in adhesion lifetimes is observed for a broad range of unfolding lengths ∆ > 0 if the reservoir coupling is strong, see Fig. 3d(ii). Self-stabilization is less effective for very small or very large values of ∆, where native and unfolded state become indistinguishable or the unfolded state becomes inaccessible, respectively.
The principle of selfstabilization can be realized in a large variety of molecular-state networks -as long as they allow for nonequilibrated states that are increasingly populated under mechanical stress. As stated in the introduction, talin unfolding and its interaction with other focal-adhesion proteins have been identified as major contributors to adhesion formation. We study here a model of selfstabilization that incorporates the binding of an adaptor protein such as vinculin to talin, see Fig. 4a. The talin rod domain contains 11 cryptic vinculin binding sites. Under load, subdomains of the rod successively unfold, thereby enabling vinculin recruitment which blocks talin refolding and promotes focal-adhesion growth [41,42]. Talin unfolding typically starts at forces around 5 pN [22,25,41] with the R3 domain, by which two vinculin binding sites are exposed. In our model, we focus on this first unfolding transition. Six additional vinculin-bound states b u, and a u, are introduced as both sites in the R3 domain can be occupied independently, see Fig. 4b. Rate constants for binding and unbinding of vinculin to unfolded talin are denoted by λ ± and are assumed to be the same for all transitions. Other model parameter values are estimated according to experimental results [21,23,40], see Supplementary Information I.
In the vinculin-bound states, talin refolding is blocked.   Thus, vinculin binding generates talin states that are not in direct contact with the reservoir and can therefore be driven out of equilibrium. Simulation results for different values of the parameter λ = λ + /λ − are shown in Fig. 4c. For λ = 0, all unbound states are equilibrated with the reservoir and no self-stabilization occurs. For λ > 1, the number of bonds increases with force, which not only stabilizes the adhesion but also translates into a reduction of adhesion sliding, see Fig. 4d. Given the high affinity with a dissociation constant in the range 10 −7 − 10 −8 M of unfolded talin for vinculin, we expect the vinculin binding constant to be larger than unity [40,[43][44][45]. Moreover, we conjecture that the remaining cryptic binding sites in talin that open at higher forces extend the demonstrated self-stabilization effect to larger loads. In summary, the minimal adhesion model can be applied to specific biological systems like the integrin adhesome, where vinculin binding to unfolded talin domains results in a reaction network containing molecular states that are populated by force application -and thus leads to self-stabilization. Analytical results of a simplified, corresponding mean-field model is given in Supplementary Information VI.

III. DISCUSSION
Our theoretical model reveals a strikingly simple mechanism for a counter-intuitive load-response of adhesions, in which tangential mechanical load can result in adhesion enhancement -instead of the ubiquitous adhesion weakening and rupture. This self-stabilization relies on molecular-conformation state-networks that are driven out of equilibrium by a mechanical load. By shifting the state occupations, the load causes a net influx of adhesion molecules from a surrounding reservoir. Notably, this self-stabilization does not require extra chemical energy, but the non-equilibrium conditions produced by the load suffice.
The primary motivation for our theoretical work are experimental observations on integrin-based focal adhesions, which display size adaptation to the applied load in planar cell cultures [14,15]. Different focal adhesion stabilization mechanisms presumably act in parallel, including actin polymerization, transcription regulation, integrin activation, and conformation changes of the adaptor protein talin. Contrasting to this complexity, we predict that adhesion self-stabilization emerges naturally in systems that merely incorporate the unfolding transition of the adhesion molecules, like talin, and a mechanism which prevents rapid bulk-exchange of unfolded states, e.g., due to vinculin binding. Since the proposed models respect basic physical constraints, such as the detailed balance in thermal equilibrium, no energy is artificially injected into the system. Additional chemical driving, e.g., through the Rap1-GTP-interacting adaptor molecule [46] or phosphorylation of vinculin or paxillin [47,48], can provide an additional layer of biological control over the suggested adhesion-stabilization mechanism.
Mechanosensitive conformation changes of adhesionlinked proteins and subsequent recruitment of additional molecules are recurring motifs in many fundamental adhesion structures besides focal adhesions, for instance adherens junctions [41,49,50] and hemidesmosomes [51]. Motor proteins also undergo mechanosensitive conformation changes and can form dynamical ensembles. Therefore, we expect that the suggested non-equilibrium mechanism for self-stabilization can help to decipher many physiological and pathophysiological processes controlled by mechano-chemical factors, and may even allow novel designs of bio-inspired, artificial adhesion systems.

A. Binding
Different states are allowed to have their own spring constants and corresponding quantities are denoted by subscripts, e.g., κ a and σ a = k B T /κ a . The probability per unit time for an unbound molecule to bind at an extension h is given by the function where k β is the intrinsic binding rate constant and b is a constant. Similar binding rate expressions have been used previously [52]. The total binding rate is obtained by integration over h so that where Erf(x) denotes the error function and the dimen- Binding of unfolded molecules via β + u (h) is defined analogously with the extension h u of unfolded molecules.

B. Unbinding
The rupture rate β − (h) is defined by the detailedbalance condition in thermal equilibrium [53], where b =˜ b + ln(σ a /σ b ) is the effective binding affinity. The rupture rate results as The rupture rate of unfolded bonds, β − u (h), is defined analogously with the extension h u of unfolded molecules.

C. Unfolding and refolding
The unfolding and refolding reaction is modeled as the transition between two local energy minima separated by a single barrier. The distance to the barrier is denoted by ∆ 1 for unfolding and by ∆ 2 for refolding. Their sum is equal to the total unfolding length ∆. The unfolding rates are thus defined as where˜ f = f /(k B T ) is a constant energy contribution for the conformation change. The reverse rates are given by The ratio of unfolding to refolding rate of bonds is given by the energy change of a bond going from extension h to h − ∆ as For unbound molecules, the total unfolding probability per time and bond is given by δ + a = k δ exp (−˜ f ) . The total refolding rate is δ − a = k δ .

COMPETING INTERESTS STATEMENT
The authors declare no competing interests.

Supplementary Discussion I. PARAMETERS
The minimal model for self-stabilization presented in the main text is motivated by focal cell-matrix adhesions. Model parameters were chosen accordingly. In focal adhesions, mechanical force is transmitted along single or multiple proteins linked in series. In the latter case, the overall spring constant κ b will be dominated by the element with the smallest spring constant. The stiffness of cellular adhesion proteins lies in the order of pN nm −1 [54,55]. The employed value of 0.25 pN nm −1 is similar to values used for previous adhesion models [56][57][58][59]. The values of κ a,b used for the talin simulations were chosen to match the experimentally measured unfolding and refolding behavior, as explained below. The thermal energy scale k B T and the spring constant κ determine the mean thermal fluctuation length σ, which is used as a length unit. Experimentally, turnover of talin is significantly faster than turnover of integrin [60]. The binding and unbinding rates of talin are determined by the intrinsic rate k β , the optimal binding distance b , and the constant energy contribution b . The intrinsic binding rate constant k β sets the time unit t 0 . The value of the binding distance b agrees with the value employed in Ref. [58] and is smaller than the fluctuation length. The strength of individual bonds is determined by F 0 = k B T / b ≈ 4 pN, as the rupture rate depends on force F as ∝ exp(F/F 0 ). The constant b can be understood as an effective affinity parameter. Its value is chosen to be rather low to fix the adhesion cluster size at equilibrium, which is proportional to exp( b /k B T ), compare Refs. [52,61].
The unfolding and refolding rates depend on k δ , f , ∆ and ∆ 1,2 with ∆ = ∆ 1 + ∆ 2 . The energy contribution f determines the ratio between folded and unfolded molecules in equilibrium. It shifts the energy barrier between both states that are separated by the distances ∆ 1 and ∆ 2 . For the simulations of talin molecules, the unfolding and refolding rate constants are chosen in accordance with the experimental results in Ref. [23]. For a time unit t 0 = 1 s, unfolding occurs at a rate of 0.015 s −1 at zero force. The folding and unfolding rates intersect at a force of 5 pN, the value at which talin R3 domain unfolding and refolding is observed [23,40,41]. The cross-linking by additional adaptor proteins is described by making use of two rate constants λ ± . In the general simulation model, the value of these constants lies within the same order of magnitude as the intrinsic rate constants k δ and k β . For integrin-based adhesions, this reaction is realized by the association of vinculin to unfolded talin domains. The molecular interactions of talin and vinculin are highly complex. At low forces, vinculin binding at unfolded talin domains strengthens the adhesion [41,62,63]. The recruitment of vinculin is proposed to act as a negative feedback loop that stabilizes the force acting on the complex [40]. Additionally, the vinculin-talin interaction also depends on the direction of forces [64]. In our model, we assume linking-rate constants with a value comparable to the intrinsic binding rate k β . In recent experiments addressing the association of vinculin with the talin R3 domain, vinculin binding has been observed with a rate on the order of 10 −1 s −1 when sufficient tension was applied to talin to induce unfolding [40]. The exchange of molecules with a reservoir is governed by the two rate constants γ + and γ − . Their ratio γ is proportional to the number of molecules in the adhesion-cluster at equilibrium, i.e., when F = 0. For comparability of the results of different models, the rate for adding molecules was increased from γ + toγ + = γ + (1 + exp(− f /k B T )) for the model without unfolding, and similarly when linking is included. This choice ensures the same number of molecules and bonds at F = 0. As described in the main text, the molecule-exchange with the reservoir plays a fundamental role for the self-stabilization mechanism. A strong connection to the reservoir allows the adhesion to grow with increasing force. The construction of molecular state models for talin in integrin-based adhesions poses a challenge due to the large number of states and due to the fact that addition and removal of an adhesion molecule requires several steps in vivo [65,66]. Furthermore, several pathways have been suggested for the recruitment of talin to integrin-based adhesions [67]. Values for the halftime of fluorescence recovery after photobleaching (FRAP) experiments for talin in focal adhesions are on the order of seconds or tens of seconds, depending on the cell type and substrate stiffness [48,68]. The exchange of single molecules can therefore be expected to occur on timescales of seconds or even faster. Hence, the values for γ ± are chosen such, that they ensure a frequent exchange of molecules with the reservoir, but still allow to observe the processes within the cluster.

Supplementary Discussion II. STOCHASTIC SIMULATIONS
Molecule-state trajectories are simulated with the Gillespie algorithm [69,70]. Its basis is provided by the stochastic formulation of chemical kinetics where the probability that a reaction i will happen in the infinitesimal time interval δτ is determined by the product of the available reactants or reactant pairs N i and a parameter µ i . Assuming j ∈ [1, . . . , m] available reactions, the function P (i, τ ) describes the probability that the reaction i is the first one to occur after a waiting time τ in the next infinitesimal time interval. For 0 ≤ τ < ∞ we have Starting from an initial configuration, the time for the next reaction and the type of reaction are drawn according to P (i, τ ) repeatedly. Force balance is restored instantaneously after each event.
In response to the tangential external force, adhesion clusters either reach a quasi-steady-state or dissociate. The timeand ensemble averages are calculated for trajectories that do not lead to complete dissociation during the simulation time. If not stated otherwise, 50 trajectories are tracked for > 10 6 single-bond transitions for the measurements of steady-state quantities, during which the tangential force is held constant. To measure adhesion lifetime, more than 200 cluster trajectories are simulated until their dissociation.

Supplementary Discussion III. ADHESION MODEL WITHOUT UNFOLDING
The following mean-field equations are used to approximate the dynamics of the stretch-dependent state-occupation numbers in the basic adhesion model without molecule unfolding where we omitted explicit t-dependence in our notation. Since the unbound molecules in state a are assumed to relax quickly mechanically, their extensions obey a Gaussian distribution. Equation (12) describes the mean total number of unbound molecules. The mean number of bonds n b (h) evolves according to Eq. (13). The drift term in (13) accounts for the average relative velocity v = ṡ between the adhesion planes that is due to the tangential force F . An expansion of n b (h) as withṽ = v/(k β σ b ) solves the steady state conditions ∂ ∂t N a = 0 and where p a (h) is a Gaussian function with zero mean and variance σ 2 a and the prime denotes the derivative with respect to the extension h. The first correction n b1 (h) is an odd function, so that its integral vanishes. The first non-vanishing correction is given by the term for j = 2. Up to second order, the integrated steady state solution for the mean number of bonds N B = N b is given by have been used and Erfc(x) denotes the complementary error function. This first correction ∝ v 2 is negative since the complementary error function for an argument x > 0 is bounded by The force-balance equation F = ∞ −∞ n b (h)κ b h dh connects the external force with the resulting mean velocity v. The lowest-order contribution (j = 0) to the force balance vanishes because n b0 (h) is symmetric. For the first non-vanishing correction, we findF The bracketed term on the right hand side of Eq.  The numbers of unbound molecules obey the following approximate rate equations The distributions of molecule numbers with extension h experience a drift due to the velocity v of the upper boundary and obey The steady-state solution for v = 0, i.e., the equilibrium solution, is given by Integration yields the total number of bonds in equilibrium as To investigate the case F > 0, the state distributions are expanded in powers of the average velocityṽ as above, see Supplementary Discussion III. The connection of the state a with the reservoir enforces N aj = 0 for j > 0. The solutions of the above equations cannot be found with an exact iterative formula. Instead, a numerical approach is used to find the corrections N au,j , which can then be inserted back into the rate equations to solve for n bj (h) and n bu,j (h). Up to first order inṽ, the mean number of bonds remains unaffected or increases with the velocity, as demonstrated for different unfolding steps ∆ and optimal binding lengths b in Fig. S3a. The first-order correction vanishes for ∆ → 0 and ∆ → ∞. The second-order correction can be both positive or negative, Fig. S3b. Strong self-stabilization is found for˜ b 1.

IV.i. Rupture behavior and adhesion cluster lifetimes
Application of a sudden force jump from F 1 = 0 to F 2 > 0 either leads to a quick dissociation of the adhesion or the adhesion system relaxes to a non-equilibrium steady state. For large adhesion clusters (N > 20), the lifetimes of those adhesions that reach the non-equilibrium steady-state usually exceeds the finite simulation time. The fraction of those adhesion clusters that rupture almost immediately, Φ rupt , is shown in Fig. S4a as a function of F 2 for different reservoir rates γ + . Without reservoir connection, when no self-stabilization occurs, the fraction of ruptured adhesion clusters sharply increases at a threshold force. Self-stabilization broadens the curves and increases the typical forces at which rupture occurs. The numerical results are well-fitted by a shifted and scaled error function. The fraction of ruptured adhesion clusters after sufficiently long simulation times represents the cumulative probability distribution for cluster dissociation after a sudden force application. Therefore, the force value F 2 at which Φ rupt (F 2 ) = 0.5 corresponds to the highest rupture probability, see Fig. S4b. The most probable rupture force increases with increasing reservoir rate γ + .  application of a force jump. Small adhesion clusters (N * ≈ 10) with different values for the reservoir rates γ ± are studied. The small molecule-numbers allow a direct measurement of the lifetimes by simulating the clusters until the last bond dissociates. For an adhesion system without molecule unfolding, realized by setting k δ = 0, the lifetime decreases monotonically with increasing force for any value of γ ± , see Fig. S5a. The equilibrium lifetime is largest for adhesion clusters without reservoir connection, realized by γ ± = 0. Figs S5b,c display lifetimes of adhesion models with unfolding molecules. The lifetime of adhesions without molecule exchange with a reservoir decrease monotonically with increasing force (black bullets). However, when the system is coupled to a molecule reservoir, so that the self-stabilization mechanism takes effect, the adhesion lifetime curves have a maximum at non-zero, finite forces. A substantial lifetime increase through self-stabilization is realized when both the reservoir rates and the intrinsic binding rate k β are large compared to the unfolding rate k δ .
IV.ii. Special case: no unfolding in the unbound state a To investigate how cyclic fluxes along the state network a − b − b u − a u − a affect self-stabilization, we set δ ± a = 0. Thereby, the cycle in the single-molecule transition-diagram is broken. Physically, this modification means that unfolding is only allowed when the molecule is bound between both planes. It should be emphasized that this model variant still allows the emergence of cyclic fluxes in the high-dimensional continuous state space spanned by the extensions of the molecules. In Fig. S6, results from the new model variant are compared with results from model IV defined in the main text. For δ ± a = 0, the increase of the adhesion molecule number with force is significantly stronger than for the model IV with cyclic flux. Thus, self-stabilization is enhanced in the absence of cyclic flux, see Fig. S6a,b. The relative velocity of the two planes bounding the adhesion is reduced accordingly, see Fig. S6c. The enhanced self-stabilization can be attributed to a stronger accumulation of molecules in the state a u , from which molecules can only escape via state b u in this model variant. Note that the force value at which first rupture events are observed does not increase greatly, see black, vertical lines in Fig. S6a. A comparison with the inset in Fig. S6b shows, that the force per bond, F/N B , is strongly reduced because of the increased number of bonds N B .  Fig. 2) and adhesion clusters with δ ± a = 0. When unbound unfolding and refolding is inhibited, the self-stabilization effect is more pronounced.

Supplementary Discussion V. ADHESION MODEL WITH UNFOLDING AND CATCH BONDS
In the models considered so far, a slip-bond behavior is assumed, i.e., the increase of the single-bond rupture rate with the applied extension h is monotonic, see sec. IV B. Thus, so far, single-bond lifetimes decreased when tension increased. Bonds that become longer-lived when tension increases are called catch bonds [45]. A number of molecular bonds in cellular adhesions have been described as catch-slip bonds [71][72][73][74]. These interactions behave like a catch bond up to some force threshold while the slip-bond behavior takes effect at higher forces. To further test the generality of the self-stabilization mechanism, simulations of different adhesion models are performed where the bonds behave like pure catch bonds. The binding and unbinding rate of folded bonds are changed to The rates for unfolded bonds are defined analogously with the extension h u . Note that the difference to the slip bond rates given in sec. IV A-B lies only in one sign in each exponential function. Thus, detailed balance still holds in equilibrium at F = 0. The steady-state simulation results are shown in Fig. S7. Again, only for model IV with a combination of unfolding and association of new molecules from the reservoir, a pronounced increase in the number of bonds is observed for small forces, see Fig. S7b. Figure S7d shows that the self-stabilization is caused by an accumulation of unbound molecules. Beyond this regime, at high forces, the pure catch-bond dynamics leads to a separation of the bond distribution into two subpopulations. Firstly, one has few, rather static molecules that carry most of the tension. Therefore, these molecules have a large extension and long lifetime. Secondly, a large number molecules form transient bonds with low, symmetrically distributed extensions. As a result, the cluster stops moving, as can be seen in the velocity-force curve in Fig. S7c for high forces. unbound molecules (N

Supplementary Discussion VI. ADHESION MODEL WITH UNFOLDING AND CROSS-LINKING
We consider an extended adhesion model that also contains the effect of adaptor protein binding. Figure S8a shows the corresponding single-molecule state-transition diagram with linking rates λ ± . Linking the load-bearing molecules with adaptor proteins to the adhesion prevents a dissociation of the molecules from the adhesion. Thereby, the linking generates states that are not directly connected to the molecule reservoir. For an approximate analytical description of the system, the set of equations from Supplementary Discussion IV is extended to account for the linked states. The unbound states obey The equations determining the evolution of the extension-dependent state distributions read As before, an expansion of the distributions for small speeds |ṽ| yields the first corrections to the equilibrium distribution. Figure S8 shows a comparison between steady-state results obtained in simulations and in the analytical approximation. Figure S9 shows simulation results for a range of parameters.