A unifying model of T-cell signaling protein condensates in reconstitution experiments

The Smoluchowski aggregation model describes the time evolution of the concentration of particle aggregates under well-mixed conditions. For LAT condensates in reconstitution experiments, aggregation is effected by the crosslinking of LAT molecules via Grb2 and SOS1. To formulate the model, we consider clusters containing $i$ LAT molecules with number density $c_i\left(t\right)$ on the SLB. Two clusters of sizes $i$ and $j$ can bind with the rate $k_{ij}$ to form a new cluster of size $i+j$. By only accounting for two-body binding events, the Smoluchowski aggregation model is written as

The first term on the right-hand side of Equation 1 describes the formation of a cluster of size $i$ through binding of clusters of sizes $i-j$ and $j$, while the second term captures elimination of clusters of size $i$ through the binding to a cluster of size $j$. To solve Equation 1, we impose the initial condition that only monomers are present at $t=0$,

Note that Equation 1 describes irreversible aggregation, i.e. the fragmentation of clusters is not accounted for. We will discuss the role of fragmentation for the formation of LAT condensates at the end of this section.

The coupled system of differential equations in Equation 1 needs to be closed by prescribing the form of the rate constants $k_{ij}$. To obtain expressions for these rate constants, the limits of diffusion-limited cluster aggregation (DLCA) and reaction-limited cluster aggregation (RLCA) are commonly considered (Müller, 1926; Jullien, 1992). DLCA is marked by high binding affinities and slow diffusion, while RLCA results from low binding affinities and fast diffusion. Given the diffusion coefficient of monomeric LAT (D_LAT ≈ 1 µm²/s), diffusion is fast on micrometer length scales compared to the minute-long experimental time scale (Figure 2; Huang et al., 2017a), suggesting RLCA is initially the dominant aggregation mechanism. However, as clusters grow, diffusion slows down to D_LAT < 0.01 µm/s in the LAT-rich phase (Huang et al., 2017a), and we expect a transition from RLCA to DLCA. To capture this transition, the rate constants $k_{ij}$ are written as the product of an encounter rate $k_{ij}^{\mathrm{B}\mathrm{r}}$ and a reaction probability $p_{ij}^{\mathrm{r}\mathrm{x}\mathrm{n}}$ (Moncho-Jordá et al., 2001),

The role of $k_{ij}^{\mathrm{B}\mathrm{r}}$ and $p_{ij}^{\mathrm{r}\mathrm{x}\mathrm{n}}$ is depicted in Figure 4. The Brownian contribution $k_{ij}^{\mathrm{B}\mathrm{r}}$ describes the rate at which clusters of size $i$ and $j$ come within a distance where they could bind, referred to as an encounter (Moncho-Jordá et al., 2001), and $p_{ij}^{\mathrm{r}\mathrm{x}\mathrm{n}}$ describes the probability that such an encounter leads to binding of the two clusters.

Figure 4

Download asset Open asset

In choosing the form of the Brownian and reactive contributions, we follow the rate constants proposed by Moncho-Jordá et al., 2001. Specifically, the reactive contribution takes the form

In this expression, $P$ is the binding probability of two clusters upon collision and $N_{ij}$ is the number of collisions per encounter of clusters of size $i$ and $j$ before they diffuse apart. While the form of $N_{ij}$ is unknown, it is argued by Moncho-Jordá et al., 2004; Moncho-Jordá et al., 2001 that it must be an increasing function of $i$ and $j$ to reflect that larger clusters collide more frequently before diffusing apart if no reaction occurs. Adopting the ansatz from these references, we use the expression

with $N_{11}$ denoting the number of collisions per encounter for two monomers and $\sigma >0$ being a free parameter.

The Brownian contribution to the binding rate $k_{ij}^{\mathrm{B}\mathrm{r}}$ is determined by the diffusion coefficients of the clusters, i.e.

with $D_i$ and $D_j$ denoting the diffusion coefficients of clusters of size $i$ and $j$ (the Brownian kernel in Equation 7 is often referenced in the more general form $k_{ij}\propto (R_i+R_j{)}^{d-2}(D_i+D_j)$, where d denotes the spatial dimension and R_i is a measure of the radius of clusters of size i [Jullien, 1992; Moncho-Jordá et al., 2004; Odriozola et al., 2001]. However, the original derivation makes the provision d > 2 [Krapivsky et al., 2010; Meakin, 1992]. Yet, dimensional analysis and an exit time analysis indicate that the characteristic time indeed takes the form of the rate constant in Equation 7, Krapivsky et al., 2010). This leads to a system-size dependent, non-dimensional prefactor that is omitted here for brevity. To proceed, we note that colloidal aggregation through both DLCA and RLCA leads to fractal clusters of particles (Meakin, 1987b; Robinson and Earnshaw, 1992). In three-dimensions, the diffusion coefficients of such aggregates can be evaluated using Kirkwood–Riseman theory (Kirkwood and Riseman, 1948; Meakin et al., 1985) or the porous sphere model (Debye and Bueche, 1948; Warren, 1994). While extensions to the collective diffusion of particles exist for free-standing lipid membranes (Oppenheimer and Diamant, 2009; Sokolov and Diamant, 2018; Sorkin and Diamant, 2021), they have not been applied to fractal aggregates or adapted to SLBs. Therefore, we make the ansatz that the diffusion coefficients depend on the hydrodynamic radius $R_{\mathrm{H}i}\propto R_0^2{i}^{{\gamma }_{\mathrm{H}}}$, where $R_0$ is the radius of gyration of a monomer and ${\gamma }_{\mathrm{H}}$ is the mobility-mass scaling exponent (Sorensen, 2011; it is well-understood that the fractal dimension depends on whether aggregation is diffusion- or reaction-limited [Meakin, 1992]. Thus, the effective fractal dimension might vary with time or might be ill-defined in our simulations, suggesting that the exponent γ_H could vary over the course of our simulations. However, we show in Appendix 1—figure 1 that our results are not sensitive to the choice of γ_H). We can then use the Evans–Sackmann theory to relate the hydrodynamic radius to the diffusion coefficient of a cluster containing $i$ particles (Evans and Sackmann, 1988),

with the translational drag coefficient given by

Here, ${\eta }_{\mathrm{m}}$ denotes the membrane viscosity, ${\epsilon }_i=R_{\mathrm{H}i}\sqrt{\frac{b_{\mathrm{s}}}{{\eta }_{\mathrm{m}}}}$ is the effective dimensionless radius with $b_{\mathrm{s}}$ being the friction coefficient describing the membrane–substrate interactions, and $K_0$ and $K_1$ are modified Bessel functions of the second kind. The drag coefficient in Equation 9 is dominated by the hydrodynamic drag arising from the membrane viscosity at small scales ($D_i\propto -\ln {\epsilon }_i$) and by the membrane–substrate friction at large scales ($D_i\propto {\epsilon }_i^{-2}$). This form of the diffusion coefficient has also been confirmed experimentally (Kaizuka and Groves, 2004).

Before proceeding, we recall that the form of the aggregation model in Equation 1 does not capture the effects of fragmentation of clusters resulting from the unbinding of crosslinkers (Huang et al., 2016). Accounting for fragmentation would lead to additional terms in Equation 1 that are linear in the concentrations $c_i$ (Krapivsky et al., 2010). While these equations are well-known, the form of the fragmentation kernels is generally not known as they depend on the internal aggregate connectivity. However, at early times, the concentration of all clusters of size $i\ge 2$ is small compared to the monomer concentrations, suggesting that the fragmentation of these clusters does not qualitatively change the dynamics even though it could lead to an effective slowing of cluster growth. At later times, we expect that the high valency of LAT quickly results in the bulk energy of a cluster dominating over its interfacial energy, thus stabilizing larger clusters against fragmentation, consistent with classical nucleation theory (Karthika et al., 2016). Thus, we do not expect fragmentation of clusters to qualitatively alter the results presented in the next section. In Appendix 1, we confirm this assertion for a simple form of the fragmentation kernel.

To analyze the Smoluchowski aggregation model, we numerically integrate the system of ODEs given by Equation 1 with the reaction kernels defined by Equations 4–9 and subject to the initial conditions in Equations 2 and 3. Table 1 shows the parameters used in our results unless stated otherwise. For the friction coefficient $b_{\mathrm{s}}$ and membrane viscosity ${\eta }_{\mathrm{m}}$, we use typical values found in the literature. However, values for the binding probability $P$, collision frequency parameters $N_{11}$ and $\sigma$, and the mobility-mass scaling exponent ${\gamma }_{\mathrm{H}}$ are not known. We focus on varying the parameter $P$ in the following and show in Appendix 1 (Figure 1) that the obtained results do not qualitatively depend on the values of the remaining parameters.

The numerical integration is implemented using MATLAB’s ode45 solver (The MathWorks Inc, 2024). To keep the integration computationally tractable, we solve Equation 1 for $i\le N_{\mathrm{m}\mathrm{a}\mathrm{x}}=4000$ and terminate our simulations once 5% of the initial mass has been lost to clusters of size $i>N_{\mathrm{m}\mathrm{a}\mathrm{x}}$. In the following, we present our results in terms of the mean cluster size, $s=\sum _{i}i{c}_i/\sum _j{c}_j$ and the non-dimensional time $t^{\ast }=\frac{C_{\mathrm{B}\mathrm{r}}{k}_{\mathrm{B}}T{\overline{c}}_{1}t}{4\pi {\eta }_{\mathrm{m}}}$, where $C_{\mathrm{B}\mathrm{r}}$ is the proportionality factor in Equation 7.

Figure 5, left panel, shows the mean cluster size for several values of the binding probability $P$. In all cases, we find an initial lag time without any significant cluster growth. This initial lag time increases with decreasing $P$. The lag time is followed by a rapid increase of the average cluster size that appears independent of $P$. After sufficiently long times, the rate at which the average cluster size grows decreases. The time and mean cluster size at which the growth rate decreases also depend on the value of $P$.

Figure 5

Download asset Open asset

We first seek to understand the origin of the lag time. To this end, consider the concentration of clusters of size $i>1$. The initial condition in Equation 3 suggests that the concentrations $c_i$ increase from zero as clusters of size $i$ form from smaller clusters. Since we consider irreversible aggregation, clusters of size $i$ are continuously converted to larger clusters, suggesting that $c_i$ reaches a maximum and subsequently decreases (see also Krapivsky et al., 2010). In Appendix 1, we find an approximate expression for the time at which the dimer concentration is maximum, i.e.

yielding

In Figure 5, right panel, we rescale time by $t_{\mathrm{c}\mathrm{h}\mathrm{a}\mathrm{r}2}$ and find that the mean cluster size collapses for all values of $P$ at early times and that $t/t_{\mathrm{c}\mathrm{h}\mathrm{a}\mathrm{r}2}\approx 1$ when the cluster size starts to rapidly increase. This suggests that $t_{\mathrm{c}\mathrm{h}\mathrm{a}\mathrm{r}2}$ indeed represents the lag time and further indicates that the onset of aggregation of LAT molecules can only proceed after initial formation of small dimers. This observation arises from the increasing binding probability of an encounter with cluster size (Equations 4–6), accelerating aggregation once sufficiently many dimers have formed.

To obtain the scaling of the characteristic time with the binding probability $P$, we note that for $P\ll 1$, μ is approximately independent of $P$. From this, we immediately find the result

This shows that the scaling of $t_{\mathrm{c}\mathrm{h}\mathrm{a}\mathrm{r}2}$ is determined by the monomer concentration and the rate at which monomers can be converted to dimers when $P\ll 1$. The validity of this assumption is discussed in detail below.

For $t>t_{\mathrm{c}\mathrm{h}\mathrm{a}\mathrm{r}2}$, the average cluster size increases rapidly. In this regime, cluster growth is determined by both the reactive and diffusive contributions to $k_{ij}$ and an expression for the growth rate is currently not known. However, we observe a second crossover to a lower growth rate with the time of this crossover increasing with decreasing $P$. The lower growth rate is associated with the transition to diffusion-limited cluster aggregation. In this regime, the growth rate of $s$ follows a power law relation,

This result can be obtained from the scaling ansatz described by Krapivsky et al., 2010 and the scaling of the diffusion coefficient for large clusters in Equation 8 as $D_i\propto i^{-2{\gamma }_{\mathrm{H}}}$. Furthermore, the cluster size at which this transition occurs can be determined by finding $i_{\mathrm{c}\mathrm{r}\mathrm{o}\mathrm{s}\mathrm{s}}$ such that $p_{i_{\mathrm{c}\mathrm{r}\mathrm{o}\mathrm{s}\mathrm{s}}{i}_{\mathrm{c}\mathrm{r}\mathrm{o}\mathrm{s}\mathrm{s}}}^{\mathrm{r}\mathrm{x}\mathrm{n}}\approx 1$, yielding the crossover cluster size for $P\ll 1$ as (Moncho-Jordá et al., 2001)

In Appendix 1 (Figure 2), we rescale the mean cluster size by $i_{\mathrm{c}\mathrm{r}\mathrm{o}\mathrm{s}\mathrm{s}}$, showing that the slowing of the growth rate is indeed determined by the transition to DLCA.

In summary, the Smoluchowski aggregation model showed that the lag time can be explained by a low binding probability of LAT monomers. This low binding probability originates from LAT binding with the help of a Grb2 and SOS1 as crosslinkers: Aggregation does not only require two LAT molecules to come within binding distance on the membrane, it also requires the crosslinkers to bind to LAT in the correct stoichiometric. ratio, leading to an effectively low binding probability in experiments. Furthermore, the binding probability is affected by the average number of phosphorylated sites, allosteric effects (Huang et al., 2016), and the binding of monomeric Grb2 or Grb2:SOS1 dimers to LAT instead of the required Grb2:SOS1:Grb2 trimer (Huang et al., 2016), effectively blocking tyrosine residues. These parameters may be tuned evolutionarily and dynamically, in response to various stimuli. Lastly, we discussed above that the model in Equation 1 assumes binding to be irreversible. In Appendix 1, we consider a simple form of fragmentation that scales inversely with aggregate size to reflect the increased stability of larger aggregates. Here, we find that the onset of rapid aggregate growth remains dominated by aggregation and is not qualitatively affected by fragmentation.

In a recent study, White et al., 2025 proposed a model that captures the time delay between TCR–pMHC binding and LAT condensate nucleation, and suggested that this mechanism constitutes a kinetic proofreading step. As opposed to reconstitution experiments, their work also accounts for the balance of kinases and phosphatases and TCR–pMHC binding kinetics. This allows them to show that the time delay observed in live cells is modulated by the slow phosphorylation of tyrosine residue Y132 (Zeng et al., 2021; McAffee et al., 2022), which preferentially binds PLCγ1. This mechanism likely acts in concert with the time delay described here. Specifically, our findings suggest that before Y132 phosphorylation, the effective binding probability $P$ is small and phase separation proceeds slowly. In contrast, upon Y132 phosphorylation, $P$ increases and condensation can proceed faster. This change in aggregation dynamics is necessary for the proofreading mechanism proposed by White et al., 2025. Before phosphorylation of Y132, condensation must be slow to avoid premature condensation, but once phosphorylation reaches a sufficient level, condensation should take off rapidly. Through this mechanism, downstream responses are suppressed initially but amplified upon successful proofreading.

To construct the field-theoretic model, we employ the insights obtained from the Smoluchowski aggregation model. Specifically, we found that we require sufficiently many dimers to form before nucleating the growth of large clusters. As discussed in the Smoluchowski aggregation model section, a threshold number of small clusters instead of dimers is needed in the case of reversible binding. This is equivalent to requiring a sufficiently large average number of bonds per LAT molecule. Once sufficiently many bonds have formed, cluster growth proceeds rapidly and macroscopic phase separation is observed. To incorporate the above idea into a field-theoretic model, we introduce the LAT density, $L(\mathit{𝒙},t)$, and the bond density per LAT molecule, $\varphi \left(\mathit{x},t\right)\in \left[0,\frac{\nu }{2}\right]$, where $\nu$ is the average valency per LAT molecule. To model the evolution of the two fields, we propose two coupled dynamical equations.

Since LAT is conserved, the evolution of $L\left(\mathit{x},t\right)$ can be related to a flux ${\mathit{j}}_L\left[\varphi \left(\mathit{x},t\right),L\left(\mathit{x},t\right)\right]$, as

We assume that the flux ${\mathit{j}}_L\left[\varphi \left(\mathit{x},t\right),L\left(\mathit{x},t\right)\right]$ can be obtained as the gradient of a free energy functional $\mathcal{F}\left[\varphi \left(\mathit{x},t\right),L\left(\mathit{x},t\right)\right]$ (Hohenberg and Halperin, 1977), but with a mobility that also depends on both densities, i.e.

To construct the free energy functional, we employ an analytic expansion that is linear in $\varphi$ and quartic in $L$. This is the lowest-order model exhibiting phase separation regulated by the bond density $\varphi$. With this choice, $L$ takes the generic form

where $f[\varphi (\mathit{𝒙},t),L(\mathit{𝒙},t)]$ is the bulk area free energy density, $\mathrm{\Lambda }$ denotes the line tension, $\delta L=L-\overline{L}$ is the LAT density with respect to some constant reference value $\overline{L}$, and ${\alpha }_i$ and ${\beta }_i$ are parameters of the free energy density. To reduce the number of parameters, we postulate that the free energy should take the characteristic form shown in Figure 6a, where the red dots indicate the minima of the free energy density $f$ at fixed $\varphi$. Specifically, we require that there is only one minimum for $\varphi \le {\varphi }_{\mathrm{t}\mathrm{h}}$, for some threshold value of the bond density, ${\varphi }_{\mathrm{t}\mathrm{h}}$, and two minima for $\varphi >{\varphi }_{\mathrm{t}\mathrm{h}}$. As shown in Appendix 2, enforcing these constraints reduces the permissible free energy density to

Figure 6

Download asset Open asset

with the new parameter ${\tilde{\alpha }}_4={\alpha }_4+{\beta }_4{\varphi }_{\mathrm{t}\mathrm{h}}$. In addition, the parameters in Equation 20 are subject to the constraints

Substituting Equations 18 and 20 into Equation 17 yields the flux of LAT as

and using this result in Equation 16 leads to the evolution equation for $L\left(\mathit{x},t\right)$,

where we have now omitted the explicit dependencies on space and time for clarity.

In the context of the Smoluchowski aggregation model, we discussed that the diffusion coefficient decreases with the cluster size and used the Evans–Sackmann theory to quantify this. This effect has also been experimentally measured in the LAT-rich phase (Huang et al., 2017a), motivating the dependence of the diffusion coefficient on the bond and LAT densities in Equation 17. However, in the field-theoretic model, there does not exist a measure of the size of an aggregate such that the Evans–Sackmann theory cannot be applied. Therefore, we propose a phenomenological model of the diffusion coefficient,

with $r$ describing the change in diffusivity, and $\omega$ being the width of the transition region where the diffusivity changes from $D_0$ to $D_{0}r$. An example of Equation 28 is shown in Figure 6b. We note that the dependence of the expression in Equation 28 on $\overline{L}$ and ${\varphi }_{\mathrm{t}\mathrm{h}}$ ensures that the diffusion coefficient is not decreased before phase separation occurs, i.e. $D\left[\varphi <{\varphi }_{\mathrm{th}},L<\overline{L}\right]\approx D_0$.

Equations 16–28 establish the governing equation for the LAT density given the bond density per LAT $\varphi$. We propose that the dynamics of $\varphi$ follow a reaction–diffusion equation,

This equation is expressed in terms of the bond density $\varphi L$ instead of the bond density per LAT to allow for an easier interpretation. The first term on the right-hand side of Equation 29 describes that a diffusing LAT molecule carries $\varphi$ bonds along with it, leading to the flux of bonds $\varphi {\mathit{j}}_L$. The remaining two terms describe the binding of two LAT molecules via a crosslinker from the bulk fluid and the corresponding unbinding reaction. We note that the rate constants $k_1$ and $k_{-1}$ also describe the effects of the bulk concentrations of Grb2 and SOS1 on the binding propensity of LAT. In here, we assume that the bulk fluids act as reservoirs for Grb2 and SOS1 such that we can assume $k_1$ and $k_{-1}$ to be constant in time.

To solve the field-theoretic model given by Equations 26 and 29, we derive their non-dimensional forms in Appendix 2. Using asterisks and hats to describe non-dimensional fields and parameters, respectively, this yields the non-dimensional flux

and the dynamical equations

The non-dimensional parameters in Equations 30–32 are given by

For ease of notation, we will omit the asterisk symbol in the following. After non-dimensionalization, we are left with the parameters $r$, $\omega$, ${\varphi }_{\mathrm{t}\mathrm{h}}$, $\overline{L}$, ${\hat{\beta }}_3$, ${\hat{\beta }}_4$, ${\hat{k}}_1$, and ${\hat{k}}_{-1}$, characterizing the diffusion coefficient, free energy landscape, and reaction rates, as well as the initial LAT density $L_0$. In the following, we hold the values of $r$ and $\omega$ constant and investigate the role of the remaining parameters.

We solve Equations 30–32 on a square domain subject to periodic boundary conditions using the Dedalus Project (Burns et al., 2020), a Python package implementing pseudospectral methods for partial differential equations. Using a Fourier basis and Dedalus’ second-order Crank–Nicolson–Adams–Bashforth (CNAB2) time stepper, we evolve the dynamical equations subject to the initial conditions,

where $L_0$ is the average initial LAT density and $\epsilon \left(\mathit{x}\right)\ll 1$ is a spatially varying, uniformly distributed noise term with zero mean. We run the spatially resolved simulations for a fixed time $T_{\mathrm{s}\mathrm{i}\mathrm{m}}$ (see Table 2) and record the observed patterns at regular intervals. For the phase diagrams in Figures 10 and 12, we only consider the final recorded snapshots of the simulations.

In this section, we present the results obtained by solving the model proposed given by Equations 30–32. We begin by considering the time before phase separation can occur, i.e. when $\varphi <{\varphi }_{\mathrm{t}\mathrm{h}}$. During this time, the diffusive flux ${\mathit{j}}_L$ is small and the LAT density remains approximately homogeneous. We can then solve Equation 32 to determine the time when $\varphi$ reaches ${\varphi }_{\mathrm{t}\mathrm{h}}$. In the flux-free case, Equation 32 is the Riccati equation, which can be cast as a linear second-order differential equation. By solving this differential equation subject to Equation 38, we obtain

where we introduced the dissociation constant

and assumed $K_{\mathrm{D}}>0$. Figure 7a plots the bond density $\varphi$ against the non-dimensional time $t^{\ast }={\hat{k}}_1{L}_{0}t$ for different values of the dissociation constant $K_{\mathrm{D}}$ according to Equation 39. At early times, the majority of binding sites are vacant, leading to a rapid increase in bond density. As the bond density increases, the number of available sites decreases, while unbinding additionally counters binding. This effects a slower rate of bond formation until binding and unbinding are balanced and a steady-state is reached. The steady-state bond density decreases with increasing dissociation constants.

Figure 7

Download asset Open asset

The dependence of the steady-state bond density on the dissociation constant suggests there exists a maximum value of $K_{\mathrm{D}}$, denoted by ${\overline{K}}_{\mathrm{D}}$, such that the threshold bond density ${\varphi }_{\mathrm{t}\mathrm{h}}$ cannot be reached for $K_{\mathrm{D}}>{\overline{K}}_{\mathrm{D}}$. Indeed, we find that the dissociation constant $K_{\mathrm{D}}$ must satisfy

for $\varphi$ to exceed ${\varphi }_{\mathrm{t}\mathrm{h}}$ as $t\to \mathrm{\infty }$. Otherwise, $\varphi$ equilibrates at a value $\varphi <{\varphi }_{\mathrm{t}\mathrm{h}}$ and phase separation cannot occur. If Equation 42 is satisfied, we obtain the time required to reach ${\varphi }_{\mathrm{t}\mathrm{h}}$ from Equation 39 as

In Figure 7b, the non-dimensional time ${\hat{k}}_1{L}_0{\tau }_{\mathrm{t}\mathrm{h}}$ is plotted against the dissociation constant for different values of the threshold bond density ${\varphi }_{\mathrm{t}\mathrm{h}}$, commonly assumed to satisfy ${\varphi }_{\mathrm{t}\mathrm{h}}\approx 1$ (George Barisas, 2003; Nag et al., 2009; Nag et al., 2012). For decreasing values of $K_{\mathrm{D}}$, the non-dimensional time ${\hat{k}}_1{L}_0{\tau }_{\mathrm{t}\mathrm{h}}$ approaches a constant value. Thus, ${\tau }_{\mathrm{t}\mathrm{h}}\propto {\hat{k}}_1$, implying that ${\tau }_{\mathrm{t}\mathrm{h}}$ is limited by the forward rate of the reaction for small $K_{\mathrm{D}}$. As $K_{\mathrm{D}}$ approaches ${\overline{K}}_{\mathrm{D}}$ from below, ${\tau }_{\mathrm{t}\mathrm{h}}$ increases rapidly due to the increasing rate of unbinding. This shows that the onset of phase separation is modulated by the effective binding constants ${\hat{k}}_1$ and ${\hat{k}}_{-1}$ as well as the initial LAT density (Equation 41), which is consistent with the Smoluchowski aggregation model (Equation 13) while providing additional insight into the lag times observed in experiments.

Lastly, we plot the change of ${\tau }_{\mathrm{t}\mathrm{h}}$ as a function of the average number of binding sites per LAT molecule for ${\varphi }_{\mathrm{t}\mathrm{h}}=1$ in Figure 7c. For every value of $K_{\mathrm{D}}$, ${\tau }_{\mathrm{t}\mathrm{h}}$ is normalized by ${\tau }_{\mathrm{t}\mathrm{h}}$ when the number of binding sites is maximum, i.e. $\nu =4$. Here, we find that reducing the number of available binding sites can increase ${\tau }_{\mathrm{t}\mathrm{h}}$ several-fold. Therefore, an additional time delay can occur when LAT is not efficiently phosphorylated, effectively reducing the number of binding sites. This could, for example, result from phosphatase pressure leading to dephosphorylation of tyrosine residues. Similarly, it has been shown that phosphorylation of tyrosine residue Y132 and subsequent binding of PLCγ1 significantly accelerates LAT condensation (Zeng et al., 2021; McAffee et al., 2022).

Thus far, we examined the behavior of LAT prior to the onset of phase separation. We now turn to analyzing the dynamics that govern phase separation. However, while our aim is to construct a minimal model, it does involve a number of (non-dimensional) parameters. Namely, the initial LAT concentration $L_0$, the parameters describing the diffusivity $\omega$ and $r$, the average number of binding sites per molecule $\nu$, the threshold bond density ${\varphi }_{\mathrm{t}\mathrm{h}}$, the free energy parameters ${\hat{\beta }}_3$ and ${\hat{\beta }}_4$ as well as the rate constants ${\hat{k}}_1$ and ${\hat{k}}_{-1}$. This relatively large number of parameters makes it infeasible to fully explore the parameter space in this article. Therefore, we focus on investigating the effects of the free energy parameters ${\hat{\beta }}_3$ and ${\hat{\beta }}_4$ and the initial LAT density $L_0$. The values of the remaining parameters are listed in Table 2. We note here that the value of $r$ is larger than measured in experiments (Huang et al., 2017a) but of the same order of magnitude. However, smaller values of $r$ lead to numerical instabilities in our model, motivating the choice in Table 2.

Figure 8 shows an example trajectory obtained from our simulations with snapshots taken at evenly distributed times. In the top row, the color coding indicates the LAT density $L$ while in the bottom row, the color coding indicates the bond density per LAT $\varphi$. At early times, both the LAT and bond density remain homogeneous while the bond density increases due to binding of Grb2 and SOS1 from the bulk fluid. Once the bond density reaches ${\varphi }_{\mathrm{t}\mathrm{h}}$ between the second and third frame, we observe the onset of phase separation and LAT-rich islands form in an LAT-depleted majority phase. As a result of the reduced diffusion coefficient in the LAT-rich phase, the observed morphology remains near-stationary, as also observed in experiments. However, the density of LAT increases between the final two frames in the LAT-rich phase. Furthermore, we observe that the bond density per LAT molecule is larger in the LAT-rich domains. This is a result of the quadratic dependence on the LAT density in the binding term of Equation 29.

Figure 8

Download asset Open asset

While Figure 8 shows an example trajectory for a choice of parameters that leads to LAT-rich, non-circular domains in an LAT-depleted majority phase, our proposed model can yield all other experimentally observed morphologies. Figure 9 shows the remaining morphologies, namely near-circular LAT-rich and LAT-depleted droplets, non-circular LAT-depleted inclusions as well as a bicontinuous phase.

Figure 9

Download asset Open asset

To understand how different parameters affect the resulting morphologies, we systematically varied the non-dimensional free energy parameters ${\hat{\beta }}_3$ and ${\hat{\beta }}_4$, the parameters associated with the cubic and quartic terms in the free energy in Equation 20, respectively. The results are shown in Figure 10 for different initial LAT densities $L_0$. Here, the blue background indicates the permissible parameter range based on the constraints posed by Equations 22 and 24. We characterize each simulation result based on the final time step and indicate the morphology using the symbols in Figure 9. In addition, we use a cross (×) symbol to indicate that no phase separation is observed by the end of the simulation.

Figure 10

Download asset Open asset

We begin by considering Figure 10b where the initial LAT density is $L_0=3.0$. In this case, phase separation is observed for all admissible parameter combinations. Furthermore, for all values of ${\hat{\beta }}_4$, small values of ${\hat{\beta }}_3$ lead to LAT-rich inclusions whereas large values of ${\hat{\beta }}_3$ give rise to LAT-depleted inclusions. This is a result of the change of location of the smallest free energy minimum at some $\delta L>0$ for ${\hat{\beta }}_3<0$ to some $\delta L<0$ for ${\hat{\beta }}_3>0$. In between these two regimes, we observe a narrow band with a bicontinuous phase. With increasing values of ${\hat{\beta }}_3$, we also observe the transition from non-circular LAT-depleted inclusions to circular LAT-depleted inclusions. This is a result of the free energy minimum associated with the LAT-depleted phase shifting toward lower densities, consequently forming smaller inclusions due to mass conservation. For these smaller inclusions, line tension becomes dominant, driving them toward circular shapes.

An analogous transition to circular LAT-rich inclusions is not observed for ${\hat{\beta }}_3<0$, despite the form of the free energy density in Equation 20 suggesting symmetry about ${\hat{\beta }}_3=0$ at equilibrium. Thus, we first emphasize that our results are not at equilibrium. Furthermore, the binding terms in the reaction–diffusion equation, Equation 29, break the symmetry of $\delta L$ with respect to $\delta L=0$. This also suggests that the symmetry with respect to ${\hat{\beta }}_3$ is perturbed. Specifically, the formation of LAT-rich inclusions leads to increased binding in those regions, thus stabilizing them. In contrast, when ${\hat{\beta }}_3>0$ the formation of LAT-depleted inclusions leads to reduced binding in these regions while the density in the LAT-rich phase is also lower than when ${\hat{\beta }}_3<0$, further reducing binding. This lower binding then favors the formation of circular LAT-depleted inclusions.

In Figure 10a and c, initial densities of $L_0=2.8$ and $L_0=3.2$ are considered, respectively. We find that in both cases, there exist parameter combinations at which phase separation is not observed. Phase separation does not occur at small and large values of ${\hat{\beta }}_3$ for $L_0=2.8$ and $L_0=3.2$, respectively. These results suggest that the constraint of LAT mass conservation does not make the free energy density minima energetically favorable solutions, thus inhibiting phase separation. For $L_0=2.8$, we now also observe the emergence of circular LAT-rich droplets not seen in the other two phase diagrams, and when $L_0=3.2$, we only observe LAT-depleted inclusions. The latter are non-circular for small values of ${\hat{\beta }}_3$ and become circular as ${\hat{\beta }}_3$ increases. For both Figure 10a and c the arguments for the transition from LAT-rich to LAT-depleted and from circular to non-circular inclusions follow the same reasoning as discussed for Figure 10b. If we now consider a fixed set of parameter values, ${\hat{\beta }}_3=-0.26$ and ${\hat{\beta }}_4=-1.98$, highlighted in red in Figure 10, we observe that increasing the density changes the morphology from circular LAT-rich inclusions to non-circular LAT-depleted inclusions to circular LAT-depleted inclusions. The same trend can be observed in reconstitution experiments, as shown in Figure 11 (Huang et al., 2016, Figure 2c).

Figure 11

Download asset Open asset

Next, we investigate how the rate constants in the reaction diffusion equation, Equation 29, affect the phase diagram in Figure 10. To this end, we consider an unbinding rate constant of ${\hat{k}}_{-1}=\frac{2}{3}\times {10}^{-5}$, twice the value used for Figure 10 and plot the corresponding phase diagrams in Figure 12. These phase diagrams reveal an overall similar structure as those in Figure 10. Most notably, however, the region where phase separation does not occur expanded in Figure 12a and c compared to Figure 10a and c. This is a result of the lower equilibrium bond density due to an increased dissociation constant, making phase separation less favorable. These observations show that changes in binding rates can affect whether phase separation occurs in experiments even when $\varphi >{\varphi }_{\mathrm{t}\mathrm{h}}$. Furthermore, comparison of Figures 10 and 12 shows that the morphologies are also impacted by the rate constants: Throughout all phase diagrams, we find small shifts in the location of phase boundaries between different morphologies. This emphasizes that the observed morphologies can also be affected by the binding parameters.

Figure 12

Download asset Open asset

The model employed in this section accounts for noise in the initial LAT concentration, but does not include Langevin noise. Accounting for Langevin noise might shift the apparent phase boundaries in Figures 10 and 12 and the onset of phase separation in metastable regions (Binder, 1987; Umantsev, 2016), as well as introduce small-scale density fluctuations (Rogers et al., 1988). However, such effects would not qualitatively affect the phase separation behavior described in this section.

Next, we discuss the physical interpretation of the free energy parameters ${\hat{\beta }}_3$ and ${\hat{\beta }}_4$, defined in Equations 33 and 34. ${\hat{\beta }}_3$ and ${\hat{\beta }}_4$ are associated with the cubic and quartic terms of the free energy, respectively, and both multiply the difference $\varphi -{\varphi }_{\mathrm{t}\mathrm{h}}$ (Equation 20). Thus, they indicate how the free energy landscape changes with changing bond densities. The parameter ${\hat{\beta }}_3$ skews the free energy landscape toward either LAT-rich or LAT-depleted minima, indicating that one becomes favorable over the other. This can be understood as cooperative phase separation due to crosslinking—if ${\hat{\beta }}_3<0$ is negative, cooperativity promotes clustering whereas the opposite is the case if ${\hat{\beta }}_3>0$. Therefore, we expect ${\hat{\beta }}_3$ to be affected by, for example, allosteric effects and the LAT valency. In contrast, the parameter ${\hat{\beta }}_4$ affects how the depth and location of the free energy density minima changes (symmetrically) with changing bond densities. For example, ${\hat{\beta }}_4>0$ may capture effects such as the inability of a condensate to restructure its bond configuration to yield higher densities, thus getting kinetically arrested. Similarly, if ${\hat{\beta }}_4<0$, the densities of the energy minima become more dissimilar, indicating the formation of very dense LAT-rich clusters at the expense of a more dilute LAT-depleted phase. Therefore, we expect that ${\hat{\beta }}_4$ is affected by crosslinker properties such as length, flexibility, bond energy and lifetime. Furthermore, ${\hat{\beta }}_3$ and ${\hat{\beta }}_4$ might also be influenced by the lipid composition, which can affect LAT–LAT interactions (Chung et al., 2021; Wang et al., 2023). Finally, we note that the above discussion is purely qualitative and that the inhomogeneous LAT and bond densities, as well as the constraint of mass conservation, render a mapping of the phenomenological free energy parameters to measurable physical parameters of the system challenging. Thus, a deeper understanding of such mappings requires further experimental and theoretical studies to enable predictive modeling of condensate behavior from molecular details.

In this section, we describe results that complement Sec. 2 of the main text. In Sec. A, we derive an approximate expression for the characteristic time of dimer formation, and in Sec. A, we discuss additional simulation results.

Supplementary results

The Smoluchowski aggregation model with the reaction kernels proposed in Odriozola et al., 2001 depends on a number of system parameters (cf. Sec. 2 of the main text). While we do not conduct a complete sweep of the parameter space in this article, we show below that the results described in the main text are sufficiently general and insensitive to the choice of parameter values. To this end, we vary the initial monomer density ${\overline{c}}_1$ (Appendix 1—figure 1a); the mobility-mass scaling ${\gamma }_{\mathrm{H}}$ (Appendix 1—figure 1b); the scaling coefficient of the number of collisions per encounter $\sigma$ (Appendix 1—figure 1c); and the number of collisions per encounter of two monomers $N_{11}$ (Appendix 1—figure 1d). The graphs in Appendix 1—figure 1 show the mean cluster size $s$ plotted against time normalized by the characteristic time of dimer formation $t_{\mathrm{c}\mathrm{h}\mathrm{a}\mathrm{r}2}$ (see Equation 49 and Sec. 2.2 of the main text). In all cases, scaling by $t_{\mathrm{c}\mathrm{h}\mathrm{a}\mathrm{r}2}$ collapses the results at early times and the onset of rapid growth occurs at $t/t_{\mathrm{c}\mathrm{h}\mathrm{a}\mathrm{r}2}\approx 1$. This confirms that $t_{\mathrm{c}\mathrm{h}\mathrm{a}\mathrm{r}2}$ indeed characterizes the onset time of rapid cluster growth. Furthermore, the qualitative behavior at later times of the aggregation process is also unaffected by the parameter choice.

Appendix 1—figure 1

Download asset Open asset

Appendix 1—figure 2

Download asset Open asset

Next, we recall the definition of the crossover cluster size $i_{\mathrm{c}\mathrm{r}\mathrm{o}\mathrm{s}\mathrm{s}}={\left(\frac{1}{P{N}_{11}}\right)}^{\frac{1}{2\sigma }}$ (cf. Sec. 2.2 of the main text). At this cluster size, aggregation becomes diffusion-limited, i.e. $p_{i_{\mathrm{c}\mathrm{r}\mathrm{o}\mathrm{s}\mathrm{s}}{i}_{\mathrm{c}\mathrm{r}\mathrm{o}\mathrm{s}\mathrm{s}}}^{\mathrm{r}\mathrm{x}\mathrm{n}}\approx 1$. In Appendix 1—figure 2, we vary the scaling coefficient of the number of collisions per encounter $\sigma$ (Appendix 1—figure 2a); the number of collisions per encounter of two monomers (Appendix 1—figure 2b); and the binding probability $P$ (Appendix 1—figure 2c). Normalization of the average cluster size by $i_{\mathrm{c}\mathrm{r}\mathrm{o}\mathrm{s}\mathrm{s}}$ shows that $s/i_{\mathrm{c}\mathrm{r}\mathrm{o}\mathrm{s}\mathrm{s}}\approx 10$ when aggregation becomes diffusion-limited. In Appendix 1—figure 2a, the case of $\sigma =0.5$ exhibits different behavior with the crossover being more gradual and occurring at a smaller normalized cluster size. In this case, $p_{ii}^{\mathrm{r}\mathrm{x}\mathrm{n}}$ changes slowly with $i$, leading to $i_{\mathrm{c}\mathrm{r}\mathrm{o}\mathrm{s}\mathrm{s}}=100$. This suggests that the scaling of the diffusion coefficient changes from a logarithmic to a power-law dependence before $i_{\mathrm{c}\mathrm{r}\mathrm{o}\mathrm{s}\mathrm{s}}$ is reached (see Sec. 2.1 of the main text and Evans and Sackmann, 1988). Then, the power-law dependence dominates the reaction rate, rendering it diffusion-limited before $p_{ii}^{\mathrm{r}\mathrm{x}\mathrm{n}}\approx 1$ such that the crossover cannot be observed at $i_{\mathrm{c}\mathrm{r}\mathrm{o}\mathrm{s}\mathrm{s}}$.

Fragmentation

In the main text and the previous sections, we assumed aggregation to be irreversible and, consequently, that aggregates grow indefinitely. If aggregation is instead reversible, aggregates can continuously break-up into smaller aggregates. In this case, the Smoluchowski aggregation model takes the form Krapivsky et al., 2010

(50) ${\displaystyle \begin{array}{c}\frac{\mathrm{d}{c}_i}{\mathrm{d}t}=\frac{1}{2}\sum _{j=1}^{i-1}{k}_{(i-j)j}{c}_{i-j}{c}_j-c_i\sum _{j=1}^{\mathrm{\infty }}{k}_{ij}{c}_j+\sum _{j=1}^{\mathrm{\infty }}{F}_{ij}{c}_{i+j}-\frac{c_i}{2}\sum _{j+k=i}{F}_{jk},i=1,2,...,\end{array}}$

where $F_{ij}$ is the fragmentation kernel that determines the rate at which an aggregate of size $i+j$ breaks up into aggregates of size $i$ and $j$. The fragmentation kernel requires detailed knowledge of the connectivity of the aggregates formed, but such details are not available within the Smoluchowski aggregation model framework. To nonetheless obtain a qualitative understanding of the relevance of fragmentation for our results, we propose the use of a power-law fragmentation kernel of the form

(51) ${\displaystyle \begin{array}{c}{F}_{ij}=\frac{F}{{(i+j)}^{\omega }},\end{array}}$

for some constants $F$ and $\omega \ge 0$. This kernel, originally introduced in Family et al., 1986, reflects that aggregates become more stable with increasing size.

We simulate this model for the parameters listed in Table 1 of the main text. For simplicity, we further choose $\omega =1$ in the following and vary the fragmentation parameter $F$. The corresponding average aggregate sizes, $s=\sum _{i}i{c}_i/\sum _j{c}_j$, are plotted in Appendix 1—figure 3, showing that $s$ reaches a steady-state where aggregation and fragmentation balance. We further observe that the early-time behavior is independent of fragmentation, suggesting that our findings regarding the onset of rapid aggregate growth described in the main text are also valid for reversible aggregation.

Appendix 1—figure 3

Download asset Open asset

To understand this result further, Appendix 1—figure 4 shows the concentration of monomers, dimers, trimers, and tetramers. As expected, the early time behavior is unaffected by fragmentation and the concentrations of all shown N-mers are increased at long times compared to the irreversible case. Importantly, accounting for fragmentation does not significantly shift the peak of the largest dimer concentration used to calculate the characteristic time $t_{\mathrm{c}\mathrm{h}\mathrm{a}\mathrm{r}2}$ in Sec. A. These findings can be understood as follows: The aggregation term is quadratic in the concentrations while the fragmentation is linear, and, at early times, $c_i\gg c_{i+1}$. Thus, aggregation initially dominates over fragmentation. We note that this result hinges on the assumption that the fragmentation rate is sufficiently small. However, Appendix 1—figure 3 indicates that for larger values of $F$ than considered here, the average aggregate size is exceedingly small such that they are not relevant for the present case.

Appendix 1—figure 4

Download asset Open asset

This section provides additional derivations for the field-theoretic model discussed in Sec. 3 of the main text. We first derive the permissible form of the free energy in Sec. B, and subsequently non-dimensionalize the governing equations of the field-theoretic model.