
==== Front
J Chem Phys
J Chem Phys
JCPSA6
The Journal of Chemical Physics
0021-9606
1089-7690
AIP Publishing LLC

39225536
5.0223001
10.1063/5.0223001
JCP24-AR-MCM2023-02584
ARTICLES
Biological Molecules and Networks
Direct computations of viscoelastic moduli of biomolecular condensates
Cohen, Banerjee, and Pappu
Note: This paper is part of the JCP Special Topic on Monte Carlo methods, 70 years after Metropolis et al. (1953).

Cohen Samuel R. 1
https://orcid.org/0000-0002-6169-1461
Banerjee Priya R. 2
https://orcid.org/0000-0003-2568-1378
Pappu Rohit V. 1a)

1 Department of Biomedical Engineering and Center for Biomolecular Condensates, James McKelvey School of Engineering, Washington University in St. Louis, St. Louis, Missouri 63130, USA
2 Department of Physics, The State University of New York at Buffalo, Buffalo, New York 14260, USA
a) Author to whom correspondence should be addressed: pappu@wustl.edu
07 9 2024
03 9 2024
03 9 2024
161 9 09510311 6 2024
19 8 2024
© 2024 Author(s).
2024
Author(s)
https://creativecommons.org/licenses/by-nc/4.0/ All article content, except where otherwise noted, is licensed under a Creative Commons Attribution-NonCommercial 4.0 International (CC BY-NC) license (https://creativecommons.org/licenses/by-nc/4.0/).
Biomolecular condensates are viscoelastic materials defined by time-dependent, sequence-specific complex shear moduli. Here, we show that viscoelastic moduli can be computed directly using a generalization of the Rouse model that leverages information regarding intra- and inter-chain contacts, which we extract from equilibrium configurations of lattice-based Metropolis Monte Carlo (MMC) simulations of phase separation. The key ingredient of the generalized Rouse model is a graph Laplacian that we compute from equilibrium MMC simulations. We compute two flavors of graph Laplacians, one based on a single-chain graph that accounts only for intra-chain contacts, and the other referred to as a collective graph that accounts for inter-chain interactions. Calculations based on the single-chain graph systematically overestimate the storage and loss moduli, whereas calculations based on the collective graph reproduce the measured moduli with greater fidelity. However, in the long time, low-frequency domain, a mixture of the two graphs proves to be most accurate. In line with the theory of Rouse and contrary to recent assertions, we find that a continuous distribution of relaxation times exists in condensates. The single crossover frequency between dominantly elastic vs dominantly viscous behaviors does not imply a single relaxation time. Instead, it is influenced by the totality of the relaxation modes. Hence, our analysis affirms that viscoelastic fluid-like condensates are best described as generalized Maxwell fluids. Finally, we show that the complex shear moduli can be used to solve an inverse problem to obtain the relaxation time spectra that underlie the dynamics within condensates. This is of practical importance given advancements in passive and active microrheology measurements of condensate viscoelasticity.

National Institutes of Health https://doi.org/10.13039/100000002 R01NS121114 Air Force Office of Scientific Research https://doi.org/10.13039/100000181 FA9550-20-1-0241 Foundation for the National Institutes of Health https://doi.org/10.13039/100000009 T32 EB028092 crossmark
==== Body
pmcI. INTRODUCTION

Biomolecular condensates are membraneless bodies that form via spontaneous or regulated phase transitions of multivalent proteins and nucleic acids.1–5 Macromolecular phase separation is a major component of the phase transitions that contribute to condensate formation.3,6,7 Condensates are defined by the presence of two or more coexisting phases. Each pair of coexisting phases is delineated by a phase boundary that gives rise to distinct equilibrium and dynamical interphase properties.8 Equilibrium interphase properties refer to differences in concentrations of macromolecules and solutes that are engendered by the equalization of chemical potentials across phase boundaries.8,9 These concentration gradients also create differences in material properties, which we refer to as dynamical interphase properties. They quantify differences in nanoscale dynamics and micrometer-scale rheology across coexisting phases.10–24 Here, we present details of a recent generalization of the Rouse model,25 which enables direct computations of viscoelastic moduli for condensates, including those formed by intrinsically disordered proteins such as prion-like low complexity domains (PLCDs).26

Systematic measurements have yielded detailed descriptions of phase behaviors of the PLCD of the protein hnRNP-A1, which is also referred to as A1-LCD.27–29 Using passive and active microrheology aided by optical traps, Alshareedah et al.26 measured the material properties of a subset of the A1-LCD variants studied by Bremer et al.29 They reported the following findings: A1-LCD and designed variants thereof form dense phases that are viscoelastic materials. The dominance of viscous vs elastic moduli depends on the physical age of condensates, and the aging process is sequence-specific. Each system is defined by a characteristic timescale known as the crossover frequency (ωc). Above the crossover frequency, elastic moduli (G′) are larger than viscous moduli (G″), whereas the converse is true below the crossover frequency. This frequency is lowered by increasing the strengths of cohesive motifs referred to as stickers. Enhancing the strengths of stickers also leads to an increase in the magnitudes of G′ and G″ across the entire frequency range that can be probed. Finally, there is a clear inverse correlation between measured intra-condensate viscosities and csat values. Overall, the measurements were taken to imply that condensates formed by A1-LCD variants are viscoelastic Maxwell fluids. The measured viscoelasticities defy assertions of condensates being purely viscous materials such as Newtonian fluids that form via liquid–liquid phase separation.14,30,31 Instead, phase separation coupled to percolation is the operative process for condensate formation of many intrinsically disordered proteins and nucleic acids.3,4,32–41

PLCDs and other intrinsically disordered proteins are exemplars of linear flexible polymers that have non-ideal configurational statistics in aqueous solvents that are relatively poor solvents for these and related systems.28,42,43 Alshareedah et al. recently generalized the Rouse model to analyze equilibrium configurations generated via lattice-based Metropolis Monte Carlo (MMC) simulations of two-phase systems formed by PLCDs.26 Here, we lay out the details of how the Rouse model was generalized using graph-theoretic descriptions of polymer configurations in dense phases44,45 to compute and analyze sequence-specific viscoelasticities of biomolecular condensates.

II. MMC SIMULATIONS

The simulations we analyze for the computations of viscoelasticity are those of Farag et al.45 They used LaSSI,34 which is a lattice-based simulation engine that relies on MMC sampling and a blend of two versions of the bond fluctuation model.46,47 The coarse-grained model of Farag et al. used one lattice bead per amino acid residue. The solvent is pseudo-implicit because a vacant lattice site is occupied by solvent molecules. This approach accounts for solvent entropy explicitly, but the solvent-mediated energetics are implicit. We analyzed data generated by Farag et al. for two systems that were designated as WT+NLS and allY. The former is the wild-type sequence of A1-LCD that also features a nuclear localization signal (NLS); in the latter, all aromatic residues of WT+NLS were replaced by tyrosine residues. The simulations use O(102) molecules and were performed across a series of different temperatures to map coexistence curves.

Monte Carlo moves on a cubic lattice are accepted or rejected based on the criterion of Metropolis et al.,48 whereby the probability of accepting a move is the min[1, exp(−∆E/kBT)]. Here, ∆E is the change in total system energy of the attempted move, and kBT is the thermal energy. Total system energies were calculated using a nearest neighbor model, whereby any two beads that are within one lattice unit of each other along all three coordinate axes contribute to the total energy of the system. The full set of Monte Carlo moves used to obtain coexistence curves included local moves, multi-local moves, reptation or slithering snake moves, translational moves, and pivot as well as double pivot moves.45 Biases introduced into the move sets, which were designed to enhance sampling, were accounted for and obviated to ensure that detailed balance is preserved.34 Farag et al.45 introduced multi-local and pivot moves that add to the move sets in the original LaSSI engine.34 The multi-local move attempts to move a single bead and its covalently bonded partners by −1, 0, or 1 lattice unit along each coordinate axis. The pivot move chooses a random chain, then chooses a random bead, x, along this chain. The move attempts to pivot every bead in the chain beyond x in the same direction by 90°.

To calculate coexistence curves, Farag et al. performed LaSSI-based MMC simulations at a series of simulation temperatures with several hundred chain molecules in the simulation setup. All simulations used a 120 × 120 × 120 cubic lattice with periodic boundary conditions. The simulations were initialized in a smaller 35 × 35 × 35 cubic lattice, and the box size was changed after initialization. This approach allows for significantly faster equilibration of coexisting phases. However, it precludes the analysis of pre-equilibrium processes, such as condensate coalescence. For each variant, 200 chains, each with 137 beads, were placed in the cubic lattice. Therefore, the total volume fraction of beads was ∼0.016. Computed coexistence curves were compared to the measurements of Bremer et al.,29 and the agreement was found to be very good.45

We analyzed the dense phases of WT+NLS and allY condensates by extracting snapshots from the LaSSI simulations. The dense phase is taken to be the largest cluster of connected chains. Two molecules are treated as being connected if at least one pair of residues from the molecules in question are nearest neighbors on the cubic lattice, which we define as belonging to any of the 26 sites that surround a given lattice site. This approach is used to quantify the numbers of intra- and inter-chain contacts for molecules in the dense phase.

III. THE ROUSE MODEL FOR VISCOELASTICITY

We first summarize the Rouse model for viscoelasticity. This summary helps us establish how we generalized the model to extend it beyond the formalism deemed to apply only to ideal chains, thus allowing us to compute viscoelastic moduli for non-ideal systems such as condensates. The theory of Rouse describes the near-equilibrium dynamics of random coil polymers in the linear response regime.49 The equations of motion for the Rouse theory can be expressed in matrix form using a Langevin equationζdXdt=−3kBTb2ZX+F.(1)

Here, ζ is the friction coefficient of the medium, kB is the Boltzmann constant, and T is the temperature. F is the random force that satisfies the fluctuation–dissipation relation. X is the column vector of the positions of Kuhn monomers of length b. The pre-factor of 3kBT/b2 comes from treating Kuhn monomers in a bead–spring model [Fig. 1(a)] with mean end-to-end distances that conform to a Gaussian distribution. A key component of the Rouse model is the connectivity matrix, Z, that describes the connectivity between Kuhn monomers.

FIG. 1. Illustration, for completeness, of the well-established and well-known Rouse theory for an ideal chain in solution. (a) Bead-and-spring model for a single chain (N = 7). The beads are Kuhn monomers, each of size b. (b) The central object of our adaptation of the Rouse theory is the graph Laplacian, which is shown here for an ideal chain of seven beads. (c) The spectrum of relaxation times, τp, for a 137-mer chain is calculated from the eigenvalues of the graph Laplacian. A fit of the dimensionless relaxation times to [1/(6π2)](N/p)2 agrees well with the slowest relaxation times. The mode corresponding to τN−1 describes the relaxation time of a single monomer, whereas τ1 describes the relaxation time of the entire chain. (d) Plots of the computed relaxation moduli for an ideal 137-mer. (e) The storage modulus (G′) and the loss modulus (G″) for an ideal chain are plotted against angular frequency. The scaling ∼ ω1/2 in the frequency range, 1/τ1 ≪ ω ≪ 1/τN−1, is observed, with G″ converging to G′. (f) Viscosities are calculated using iωη* = G′ + iG″. The quantities here are normalized to the zero-shear viscosity, η0, and the polymer-mediated viscosity of the solvent, ηs. In (c)–(f), we plot the dimensionless quantities.

We generalize the Rouse model by relaxing the imposition of Gaussian chain statistics, which is achieved by recasting Z. This approach takes advantage of the fact that the connectivity matrix is formally equivalent to a graph Laplacian,44 wherein polymers or beads within a polymer are nodes that are connected to other polymers or beads by edges. We generalized the Rouse model to allow for the nodes to be monomers within a chain or entire chains. The edges are either covalent bonds, a combination of covalent bonds and non-covalent contacts for a real chain, or non-covalent contacts between chains if each node is a single chain. The graph Laplacian is defined as the difference between the degree (D) and adjacency (A) matrices50 such that Z ≡ (D − A). In the degree matrix D, the off-diagonal elements are zero, and the diagonal elements quantify the number of connections for each node. Elements of the adjacency matrix, A, have a value of 1 if two nodes are connected or 0 if they are not connected by either covalent or non-covalent crosslinks. The graph Laplacian is a version of the connectivity matrix introduced by Zimm,51 with the difference being that our version ignores hydrodynamic interactions because we are interested in the contributions made by equilibrium structures of condensates.

To aid in understanding the structure of the graph Laplacian, we first show how this is constructed for the familiar case of an ideal linear chain as introduced by Rouse.25 The graph is undirected and unweighted. For an ideal chain, the diagonal elements of the graph Laplacian have values of 2 or 1, depending on whether the monomers are internal to the chain or are at the ends of the chain. A value of 2 signifies the fact that each monomer is connected to one bead on each side, whereas the two terminal beads have one connection each. For an ideal chain, all off-diagonal elements other than the band diagonal elements are zero. The band diagonal elements have values of −1 for an ideal chain. This sparsity is not the case for real chains.

Equation (1) can be solved by introducing normal coordinates where the positions are independent of one another.25,52 We can rewrite Z as Z ≡ BBT, where B is the unoriented incidence matrix in graph theory. Accordingly, for an undirected graph consisting of N nodes and m edges, the elements of the N × m incidence matrix are set to 1 if there is a node i that is incident upon edge j; otherwise, the elements are set to zero. The incidence matrix is used to rewrite Eq. (1) asζdXdt=−3kBTb2BBTX+F.(2)

Following the approach of Rouse,25 we define a new set of coordinates, r, which we can relate to the original coordinates usingr=BTX.(3)

Multiplying both sides of Eq. (2) by BT givesζdrdt=BTdXdt=−3kBTb2BTBBTX+BTF=−3kBTb2Rr+BTF,(4)

where R ≡BTB. If we only consider the first term in Eq. (4), because we are interested only in fluctuations about equilibrium, the solution to Eq. (1) becomes a solution to an eigenvalue problemRuk=λkuk,(5)

for eigenvectors uk and eigenvalues λk for k = 1, 2, …, m, where m is the total number of springs in the Rouse model or edges in our graph-theoretic representation of the Rouse model. The solution to Eq. (1) is the set of nonzero eigenvalues of either the m × m Rouse matrix, R, or of the N × N graph Laplacian, Z, which are equivalent.44 Note that for N monomers or nodes, there are N − 1 nonzero eigenvalues. This follows from spectral graph theory, where the first eigenvalue of the graph Laplacian is zero for a fully connected graph.53

A. Eigenvalues of the graph Laplacian and their use

We are interested in the eigenvalues of the graph Laplacian,Zup=λpup,(6)

for p = 1, 2, …, N − 1. Following Zimm,51 we diagonalize Z usingQ−1ZQ=Λ,(7)

where Q is the N × N matrix with the eigenvectors up as its columns, and Λ is the diagonal matrix of the eigenvalues λp. The matrix Q can be used to transform the original coordinates X into normal coordinates q according toq=Q−1X.(8)

Taking the time derivative of both sides of Eq. (8) leads todqdt=Q−1dXdt.(9)

By comparison with Eqs. (1) and (8), we haveζdqdt=ζQ−1dXdt=−3kBTb2Q−1ZQQ−1X+Q−1F,(10)

and this simplifies toζdqdt=−3kBTb2Λq,(11)

where we have neglected the contribution of the random force since it is not needed to solve the eigenvalue problem.52 To arrive at Eq. (11), we use QQ−1 = Q−1Q = I, where I is the identity matrix. The normal coordinates are related to their time derivative byζdqpdt=−3kBTb2λpqp,(12)

and this yields the expressionτp=ζb26kBTλp.(13)

The expression in Eq. (13) defines the full distribution of relaxation times, τp.52 Therefore, the distribution of dynamical modes, as captured in the spectrum of relaxation times for an ideal polymer in solution, comes from knowledge of the graph Laplacian in the Rouse formalism. We adapt this feature for real chains and analyze simulations of two-phase systems by constructing graph Laplacians that describe the equilibrium structures within dense phases. Prior to describing this adaptation, we discuss key dynamical properties of an ideal polymer in a solvent that can be extracted by analyzing the spectrum of relaxation times.

B. Relaxation times for an ideal

The use of normal coordinates, as shown in Eq. (11) resolves the motions of the polymer into a series of modes. Each of the modes is defined by a characteristic relaxation time. For an ideal chain of N Kuhn monomers, there are N − 1 non-zero modes. Since A1-LCD is a 137-mer, we first computed the distribution of relaxation times, τp, for an ideal 137-mer chain [Fig. 1(c)]. These relaxation times are inversely proportional to the eigenvalues of the graph Laplacian. The mode corresponding to τN−1 describes the relaxation time of a single monomer, whereas τ1 describes the relaxation time of the entire chain. A fit of the dimensionless relaxation times to [1/(6π2)](N/p)2 agrees well with the slowest relaxation times [Fig. 1(c)]. As noted by Rouse,25 this fit does not work for the fastest modes. Indeed, we find a clear deviation for the higher p modes, which correspond to the dominant contributions made by the elasticity of the ideal chain [Fig. 1(c)].

Following Rouse,25 we use the relaxation times, τp, to compute the complex shear modulus in terms of the storage and loss moduli usingG′ω=ϕkTNb3∑p=1N−1ω2τp21+ω2τp2,(14)

andG″ω−ωηs=ϕkTNb3∑p=1N−1ωτp1+ω2τp2.(15)

Here, ϕ is the volume fraction of polymers in solution, ω is the angular frequency, and ηs is the viscosity of the solvent. We do not need to know the solvent viscosity since we compute the loss modulus dereferenced against (ωηs) in Eq. (15), and all quantities on the right-hand side are known. For illustration, we set ϕ, b, and ζ to be unity, and T is in simulation units with k set to unity. These choices are made because the quantities are pre-factors that do not materially alter any of the conclusions. We use a simulation temperature of 53, which corresponds to 296.8 K for A1-LCD,28,45 thus closely matching experimental conditions. Using the same values for the constants, we also compute the relaxation modulus for an ideal chain, which is the time-domain analog of the complex shear modulus,Gt=ϕkTNb3∑p=1N−1exp−t/τp.(16)

The relaxation modulus scales as ∼tτN−1−1/2 in the regime between τN−1 and τ1. Analysis of the eigenvalues for the ideal chain reveals well-established behaviors for the frequency dependence of the storage and loss moduli. In the low-frequency regime, G′ < G″, and the storage modulus G′ scales as ∼ω2, whereas the loss modulus G″ scales linearly with ω. In the high-frequency regime, G′ > G″, and while G′ plateaus, G″ decreases as ∼ ω−1. Between the regimes where G′ < G″, and G′ > G″, there is a regime where G′ = G″. In this regime, 1/τ1 ≪ ω ≪ 1/τN−1, both G′ and G″ scale as ∼ω1/2.

The complex shear moduli illustrate the well-established result that even an ideal polymer in a dilute solution has finite elastic and viscous moduli. Given the N − 1 finite relaxation times, Rouse proposed that the viscoelastic properties of ideal polymer solutions conform to a generalized Maxwell model.25,54 In the Rouse formalism, Maxwell elements representing each of the relaxation modes contribute on the order of kBT to the overall storage modulus at high frequencies. One can then calculate complex viscosities using iωη* = G′ + iG″. Plots of these normalized quantities reproduce the results of Rouse for an ideal chain.25

Next, we asked how the inclusion of the effects of non-bonded interactions in the form of contacts between non-nearest neighbors along a linear chain influences the storage and loss moduli of an ideal chain. These additions, which are relatively sparse, are made to ensure that the contacts do not significantly perturb the ideal chain statistics. This computational experiment is intended to be illustrative and is designed to explore how controlled deviations from ideal chain behavior influence equilibrium dynamics. As shown in Fig. 2, we consider a 137-mer. For an ideal chain, the off-diagonal density of the graph Laplacian is ρ = 0.0146. This non-zero density is entirely due to bonded interactions. We randomly added non-bonded connections to this system by increasing the density of off-diagonal elements. The new result shows that as random contacts are added, the intermediate region, where G′ = G″ and the moduli scale as ∼ω1/2, shrinks such that the crossover region corresponding to the equality of the moduli becomes a crossover point. Adding random non-bonded contacts to an ideal chain yields frequency-dependent profiles for the complex shear moduli that are concordant with a conventional Maxwell model.55 The changes to the viscoelastic responses observed by the addition of random, non-bonded interactions motivated computations of graph Laplacians matrices by analyzing contact patterns from the LaSSI-based MMC simulations of phase equilibria for WT+NLS and allY versions of A1-LCD.

IV. GRAPH LAPLACIANS EXTRACTED FROM MMC SIMULATIONS FOR PLCDS

Analysis of the LaSSI simulations showed that chains are non-ideal and more expanded within dense phases when compared to the non-ideal dilute phase.45 To adapt the Rouse model to model non-ideal systems such as condensates, we extracted graph Laplacians matrices directly from lattice-based MMC simulations of the A1-LCD system [Fig. 3(a)].

FIG. 2. Assessments of how non-bonded contacts impact computed moduli of an ideal chain. (a) Starting with the 137-mer ideal chain, we randomly add contacts by increasing the off-diagonal density of the graph Laplacian starting with the ideal chain at ρ = 0.0146. We compute the dynamical moduli (solid lines) and compare the results to the ideal chain (dashed lines) for (b) ρ = 0.015, (c) ρ = 0.017, and (d) ρ = 0.019. The intermediate region in the dynamical moduli disappears as random contacts are added.

A. The single-chain graph

We first tested an approach that we termed the single-chain graph [Fig. 3(b)]. Here, the nodes represent residues on a single chain within the condensate, and the edges are defined by bonded and non-bonded intra-chain contacts. This represents a generalization of the Rouse model. The graph Laplacian matrix was calculated as Z ≡ BBT, where the unoriented incidence matrix B has elements Bij. The elements Bij are set to 1 if there is a node i incident upon edge j; otherwise, Bij equals zero. For N nodes and m edges, the incidence matrix has size N × m. This allows us to consider the effect of edge-weighting the graph, although we used unweighted graphs because previous work26 showed that weighting does not have a material impact on our findings.

For distinct snapshots in the simulations, chains located within the interiors of dense phases were analyzed for the presence of non-bonded intra-chain contacts. For each snapshot of the dense phase, we construct a graph Laplacian for every single chain and compute the eigenvalues of the graph Laplacian. Each chain has 136 relaxation modes that yield finite relaxation times. These eigenvalues are then used to compute the average relaxation time by averaging over the 136 modes. These averaged relaxation times are then used to compute the snapshot-specific moduli. We then compute the ensemble averages of these moduli using the snapshot-specific moduli.

For the single-chain graph used to capture the dynamics within a condensate, the Rouse time, τ1, which is the slowest relaxation mode, is longer than for an ideal chain in a dilute solution. This is due to internal friction imparted by the non-bonded interactions. However, when compared to what we observe by adding random, non-bonded contacts at higher off-diagonal densities [Figs. 2(c) and 2(d)], the intermediate regime where the moduli scale as ∼ω1/2 is largely preserved, although the storage and loss moduli are unequal in this regime. This highlights the sparsity of intra-chain contacts for A1-LCD molecules within the dense phases. The sparsity of intra-chain contacts comes from the fact that chains form expanded conformations within dense phases.45

FIG. 3. The single-chain graph and comparisons of computed and measured moduli. (a) Snapshot of a single chain extracted from MMC simulations of the dense phase of A1-LCD: WT+NLS. Each residue is a Kuhn monomer. (b) The graph Laplacian accounts for both bonded and non-bonded intra-chain contacts. The inset indicates a region of the graph Laplacian showing the types of intra-chain contacts that form in dense phases. (c) The corresponding dynamical moduli show an intermediate response scaling as ∼ω1/2. The Rouse time is longer than the Rouse time, τ1, for the ideal chain due to internal friction imparted by the non-bonded interactions. (d) Comparisons to measurements following a rescaling to match the measured crossover frequency show that the single-chain graph overestimates the experimental storage and loss moduli. The computed moduli were analyzed over 30 snapshots across three replicates by averaging the pth relaxation times over all chains in the dense phase.

Next, we compared the moduli computed using the single-chain graph to the moduli measured using passive microrheology with optical trapping.26 To put the comparisons on a quantitative footing, we extracted the crossover frequency ωc from the experimental data and scaled the computed crossover frequency to match the measured value. Rescaling the moduli involves using a single-point constraint such that the value of the crossover of the computed moduli is set equal to the value from the experiment. We multiply all frequency values of the computed moduli by the crossover frequency of the experiment and divide it by the crossover frequency of the computed moduli. We also computed moduli by the value of the experimentally measured moduli at the crossover and divided this by the value of the computed moduli at the crossover. This rescaling procedure brings the computations into the same currency as the experiments. This is essential since the lattice simulations do not have any representation of solvent-mediated interactions. Importantly, the rescaling based on information regarding the crossover frequency from the measurement and nothing more allows us to ask if this is sufficient for recapitulating the frequency-dependent values of the moduli across the entirety of the experimentally accessible range.

Comparisons of the computed and measured storage and loss moduli show that the single-chain graph overestimates moduli across the frequency range that was probed experimentally. While there is reasonable agreement with the loss modulus in the intermediate frequency range, the storage modulus computed using the single-chain graph is an overestimate when compared to the measurements. These deviations were surprising given the excellent agreement between measured and computed phase diagrams for A1-LCD and numerous designed variants thereof.45 This prompted us to appreciate that the measurements, which construct complex shear moduli from autocorrelation functions of the motions of beads within condensates, investigate the material properties of the dense phase as a whole and not just the contributions of individual chains. Furthermore, the concentrations of polymers within dense phases have been measured to be above the overlap concentration with the measured concentrations corresponding to the semidilute regime.45 In this regime, the differences between intra- and inter-chain interactions vanish, and chain segments, referred to as thermal blobs,56 move under the influence of inter-segment interactions, which may be within the same chain but are more likely to be between different chains. Additionally, dense phases may be viewed as percolated networks.4,34,35,57 Based on these considerations, we reasoned that the network of chains making up the dense phase can be treated as a single chain leading to what we refer to as the collective graph.

B. The collective graph for computing graph Laplacians

In the collective graph, the nodes are individual chains within the dense phase [Fig. 4(a)]. For a pair of molecules, an edge is established if at least one pair of residues, one from each chain, are nearest neighbors on the cubic lattice. The nearest neighbor belongs to any of the 26 sites that surround a given lattice site. The graph constructed in this way is undirected, meaning that the edges are symmetric. In the dense phase, the number of inter-chain contacts per chain is at least an order of magnitude larger than the number of intra-chain contacts per chain [Fig. 4(b)]. Using this approach, it was previously shown that the internal organization of condensates formed by A1-LCD and other molecules corresponds to a small-world network.45 Given this observation, the coarse-graining procedure used to generate graph Laplacians for the collective graph, whereby individual chains are replaced by a node [Fig. 4(c)], has precedent in the literature. Specifically, Song et al. used block renormalization to demonstrate that complex networks, including those with small-world topology, have self-similar, scale-free properties.58 This justifies our coarse-graining, whereby each node is now a single chain. The hub-and-spoke nature of the internal organization within dense phases gives rise to graph Laplacians that are differently dense or sparse [Fig. 4(d)].45

FIG. 4. In the collective graph, the entire polymer network is treated as a single chain. (a) We analyze dense phases from LaSSI simulations of WT+NLS. (b) Chains in the dense phase have a high likelihood of forming inter-chain interactions.45 (c) A coarse-grained model of the dense-phase network was developed by treating each chain as a node. An undirected edge between nodes indicates that at least one pair of residues between the chains forms a contact. (d) The graph Laplacian corresponding to the dense phase in panel (a). The insets indicate regions with different sparsities.

For each chain in the network, which is represented as a node in the graph, the number of edges associated with the node is determined by the number of intermolecular contacts involving the node in question. The number of edges is distinct for each node and for each snapshot. We have also constructed the graphs using the same criterion from Farag et al., with the graphs for the collective graph defined by sticker–sticker interactions rather than any residue–residue interaction. There was no significant difference in the moduli, which suggests that edge-weighting the graph Laplacians by the number of intermolecular contacts would have no effect. Instead, sequence specificity comes from the interaction model used for the MMC simulation. The calculations generate sequence-specific moduli because the ensembles are sequence-specific and so are the graph Laplacians. To make this point, Fig. 4(c) shows a sub-graph depicting the fact that in the collective graph, inspired by the block renormalization procedure of Song et al., each chain is modeled as a node on the graph, and edges between chains are determined by the inter-chain, intermolecular interactions. Therefore, neither all nodes nor all monomers are equivalent. Instead, there is sequence specificity and heterogeneity.

We compared the moduli, computed using the collective graph to the measured moduli, using the same rescaling of frequencies that was used in the comparison of the single-chain graph [Fig. 5(a)]. Inclusion of the inter-chain contacts and treating the network of chains as a single chain leads to almost perfect agreement between the measured and computed loss moduli and good agreement across the intermediate- and high-frequency ranges for the storage moduli. This suggests that the actual network is more compliant than is anticipated by the single-chain graph, meaning that the network can more readily deform in response to internal stresses.19 This point is made by comparing the relaxation modulus computed using the single-chain vs collective graphs [Fig. 5(b)].

FIG. 5. Performance of the collective graph, the computed relaxation modulus, and the mixture of two graphs. (a) The storage and loss moduli computed using the collective graph are plotted alongside the experimental values for WT+NLS. (b) The relaxation modulus shown here is computed using both the single-chain and collective graphs. The collective graph leads to a shortening of the relaxation times as the network is less stiff and more compliant. In (a) and (b), 30 snapshots were analyzed over three replicates. The error bands indicate the standard deviation. (c) A linear combination of the moduli computed using the single-chain and collective graphs reproduces the experimental values with a coefficient of ≈0.88 for the collective graph.

C. Mixture of graphs

Returning to the comparison between computed and measured moduli, we note that the storage moduli in the measurements deviate from the computed ones at low frequencies. Alshareedah et al. speculated that deviations between computations and measurements in the low-frequency regime might have been due to larger errors in the storage moduli in this frequency range. This seemed like an intuitive explanation since the moduli were derived from passive microrheology and there are fewer data points from which to extract the long-time, low-frequency portions of the moduli. We revisited the totality of the data for all the A1-LCD variants studied by Alshareedah et al. and noted that the error bars in the low-frequency range are ∼5%–10% larger than the error bars in the high frequency range. Given the consistency of the size of the error bars across all variants, it appears that deviations between the collective graph and measurements are due to a systematic underestimation of the storage moduli by the collective graph. Moreover, since the deviating behavior appears only in the storage modulus, it does not seem likely to be an artifact of the Fourier transform of the positional autocorrelation function.

While the collective graph does a good job of recapitulating the elasticity of the network at intermediate and high frequencies, it underestimates the elasticity at low frequencies. Individual polymers are not coarse-grained beads freely moving around one another. While this coarse-graining is valid for describing the network, it appears that on longer timescales the motions combine attributes of the collective and single-chain graph. This hypothesis is based on the observation that the measured values for storage moduli lie between those of the single-chain and collective graphs at low frequencies. We asked if a mixture of the two graphs can be parameterized to capture the totality of the measured frequency dependence of storage and loss moduli. In the mixture of graphs, the measured storage modulus is postulated to be a linear combination of the storage moduli computed from the collective (cm) and single-chain (sm) models,Gexp′ω=aGcm′ω+bGsm′ω.(17)

This is subject to the constraint that the coefficients a and b are positive scalars and that a + b = 1. To determine the values of the coefficients, we performed an optimization using gradient descent. The learning rate, γ, was 0.06, and we iterated up to 2000 times or until convergence, with a numerical tolerance of 10−6. All data were first normalized to ensure the reliability of the procedure. This was performed by scaling and centering the data such that they have a mean of zero and a standard deviation of one. For each iteration, the error was computed as the difference between the linear combination as defined earlier and the normalized target (experimental) data. The coefficient a was updated according toan+1=an−γd.(18)

Here, d is the arithmetic mean of the error multiplied by the normalized test (cm) data. The coefficient a was then reset to the original scale by multiplying the value by the ratio of the standard deviation of the unnormalized experimental data to the standard deviation of the data from the collective graph. Using this procedure, we obtain a ≈ 0.88 and b = 1 − a.

The optimized values of the coefficients a and b were used in Eq. (17) to compute the storage and loss moduli, and the frequency dependence was compared to the measured moduli. Note that the estimation of the coefficients did not use the measured loss moduli as a target. The mixture model generates profiles that are in very good agreement with the measurements for both the storage and loss moduli [Fig. 5(c)]. The physical interpretation is that viscoelastic moduli of dense phases are influenced by a combination of the inter-chain contacts that give rise to a percolated, small-world network and the rigidity that arises from the segments being part of a polymer. Accordingly, threadlike molecules, to use the phrasing of Rouse,25 form networks through a collection of inter-chain contacts that make and break across a range of timescales. The long-time motions are hindered by the internal friction of a single chain,59,60 and this diminishes the compliance of the network. As a result, both inter-chain interactions, which are dominant in semidilute solutions, and intra-chain contacts, which are captured in the single-chain graph, jointly, albeit disproportionately, contribute to the totality of the frequency-dependent viscoelasticity of A1-LCD condensates. The precise parsing of the contributions from the single-chain vs collective graphs is likely to be system-specific.

Before concluding this section, we note the following: Below the crossover frequency, the storage modulus is smaller than the loss modulus. However, this does not imply that the dense phases are purely viscous materials that, via physical aging such as a glass transition61 or conversion to solids,62 become dominantly elastic.63 Such binary characterizations of early time condensates as being viscous fluids that age and convert to elastic solids gloss over the fact that even on timescales when the loss modulus is larger than the storage modulus, the two moduli are within an order of magnitude of one another [Fig. 5(c)] and the materials are viscoelastic. This negates the characterization of liquid-like condensates as purely viscous Newtonian fluids.16,30,61,64,65 Instead, the finite storage modulus gives dense phases a network-like structure,37,45,63 and this will have a direct influence on the mechanical forces transduced, sensed, and exerted by condensates.24,66–71

V. SPECTRUM OF RELAXATION TIMES

Measured moduli of condensates suggest that these materials belong to the same class as Maxwell fluids.10,12,26,61 A Maxwell element features a spring and a dashpot in series.55,72 In coarse-grained descriptions based on appropriate constitutive equations, a Maxwell fluid with a single element will generate a frequency-dependent profile akin to what has been measured for A1-LCD-based systems.55,61,73 This has been taken to mean that condensates can be reduced to a single Maxwell element and that the presence of a single crossover frequency implies that condensates are defined by a “single relaxation time.”61 However, the number of relaxation modes is governed by the number of non-zero eigenvalues of the graph Laplacian. Therefore, each mode can be thought of as a Maxwell element.25 Accordingly, if dilute or semidilute polymer solutions are dominantly viscous at long times, then they are best described as generalized Maxwell fluids,63 whereby Maxwell elements, one per mode, are assembled in parallel.74 While the relaxation spectrum, hτ, cannot be measured directly, it can be extracted by solving an inverse problem because it is related to the storage and loss moduli viaG′ω∼∫−∞∞ω2τ21+ω2τ2hτd⁡ln⁡τ,(19)

andG″ω∼∫−∞∞ωτ1+ω2τ2hτd⁡ln⁡τ.(20)

Note that the continuous relaxation spectrum is related to the discrete spectrum byG′ω∼∫−∞∞ω2τ21+ω2τ2hτd⁡ln⁡τ=∑p=1N−1gpω2τp21+ω2τp2,(21)

andG″ω∼∫−∞∞ωτ1+ω2τ2hτd⁡ln⁡τ=∑p=1N−1gpωτp1+ω2τp2.(22)

For N monomers or chains, there are N−1 Maxwell modes giving finite relaxation times, and gp,τp describe the weights and relaxation times, respectively, for the pth mode. We do not have a priori knowledge of the weights, but we have the measured moduli, and we also have the computed moduli that match the measurements via the mixture of graphs. Accordingly, we can use either the computed or measured moduli as joint inputs and solve for hτ. We solve this inverse problem using a nonlinear Tikhonov regularization method.75 We use the substitution hτ=expHτ so that hτ>0. We then minimize the cost function given byVλ=∑iGm/c′ωi−G′ωi;HτGm/c′ωi2+∑iGm/c″ωi−G″ωi;HτGm/c″ωi2+λ∫−∞∞d2Hτdτ22d⁡ln⁡τ.(23)

Here, m/c implies measured (m) or computed (c) moduli, and we use one or the other, but not both. In the current setup, we use the measured moduli. Equation (23) can be rewritten compactly as Vλ=σ2+λη2. The first term, σ2, denotes the mean-square error between the measured or computed moduli, Gm/c′ωi and Gm/c″ωi, and the moduli inferred from the optimization procedure, G′ωi;Hτ and G″ωi;Hτ, respectively, for the ith frequency. Inferred moduli appear in the optimization because we are estimating the relaxation spectra while minimizing the deviations from the measured moduli. The second term consists of the norm of the curvature of Hτ and the regularization parameter, λ, which controls the smoothness of Hτ. To obtain the relaxation spectrum, we used the approach developed by Takeh and Shanbhag76 that first uses a least-squares method to find the Hτ while minimizing the cost function and then determines the optimal λ from the L-curve in the log–log plot of η vs σ. We then compute hτ=expHτ. As a consistency check, we used the computed hτ to solve the forward problem and extract the measured moduli from the estimated relaxation spectra.

Figure 6 shows the relaxation spectra that are derived from the measured moduli. The relaxation time spectrum, hτ, shows a broad, continuous distribution for the WT+NLS and allY systems. Accordingly, condensates formed by these systems are generalized Maxwell fluids featuring system-specific distributions of relaxation modes. Each relaxation spectrum has a distinct peak that corresponds to the system-specific crossover time in the dynamical moduli. Note that this crossover time shifts to longer values for the allY system compared to the WT+NLS, and this is consistent with the downshift of the crossover frequency. Our analysis shows that measured moduli can be used to extract the spectrum of relaxation modes for condensates. This highlights the superiority of direct measurements of dynamical moduli over measurements of single-chain properties within condensates. This does not need the use of optical traps, although the methodology affords several advantages. Instead, light scattering or particle tracking measurements can be adapted to quantify complex shear moduli,26,77,78 and these measured moduli can be used to infer the relaxation spectra.

FIG. 6. A continuous spectrum of relaxation times indicates that a generalized Maxwell model best describes dense phases of PLCDs. (a) We calculated the spectrum for the WT+NLS from the experimentally measured moduli in Fig. 5(a) following a nonlinear regularization strategy. As a consistency check, we solve the forward problem exactly from the knowledge of hτ, thereby recovering the storage and loss moduli. (b) We also compute the spectrum from the experimentally measured moduli for the allY variant following the same approach. The arrow in each panel indicates the peak in the relaxation time spectrum corresponding to the crossover in the dynamical moduli.

VI. DISCUSSION

In this work, we have provided details of a recently introduced generalization of the Rouse model that incorporates different types of graph Laplacians whose eigenvalues are used to enable direct computations of viscoelastic moduli of dense phases formed by PLCDs. The graph Laplacians are computed using a graph-theoretic representation of lattice-based MMC simulations of condensate-forming systems. Two types of graph Laplacians were constructed: the single-chain graph that accounts only for intra-chain contacts, whereas the collective graph accounts for inter-chain contacts. Overall, the collective graph generates better agreement with the measured moduli. However, this model is not perfect, and a mixture of the two graphs with dominant contributions from the collective graph explains the totality of measured, frequency-dependent moduli. The collective graph, which is in accord with the dense phase being a semidilute solution,45 generates a compliant network. This model conceptualizes condensates as a collection of blob-sized segments that interact with one another to generate a percolated network. However, the motions of the segments are not unhindered, and the long-time dynamics require a proper accounting of the polymeric nature, as evidenced by the need for a mixture of graphs to account for the totality of the frequency dependence of storage and loss moduli [Fig. 5(c)]. At long timescales, the collective graph overestimates the compliance of the network, and the slithering of single chains and the sparsity of intra-chain contacts needs to be accounted for to generate the correct, frequency-dependent moduli as measured using microrheology.

Building on the foundational work of Rouse25 and aided by the analysis of MMC simulations of two-phase systems, we demonstrate explicitly that the presence of a single crossover frequency is not sufficient to assert that there is a single relaxation time within a condensate.61 Instead, the computations show that there is a spectrum of relaxation times. Given the congruence between measured moduli and expectations for a Maxwell fluid, we propose that the dense phases may be viewed as generalized Maxwell fluids comprising a collection of Maxwell elements being assembled in parallel.25,63 This contrasts with there being a single Maxwell element that describes the entire system. Accordingly, condensates that are dominantly viscous are best described as generalized Maxwell fluids that are viscoelastic. The latter point is noteworthy because the storage moduli cannot be ignored even at long times, and both the storage and loss moduli contribute to the dynamics of condensates across the entire frequency range.

The methods introduced here show how dynamical moduli of dominantly viscous, albeit viscoelastic materials can be computed from simulations of two-phase systems formed by disordered proteins, which are becoming increasingly tractable both at the coarse-grained and fine-grained levels.28,38,45,59,65,79–86 However, for MMC simulations, we do not have information about timescales, nor do we have access to hydrodynamic interactions that are mediated by solvent contributions. While molecular dynamics simulations do not have such issues, the timescales that are accessible cannot span the range of relaxation modes that define dynamics within condensates.59,65,86 Accordingly, analysis of simulations, without any input from experiments, can be useful for extracting robust insights that enable relative comparisons, but the computed values of absolute moduli will be unreliable. Knowledge of the crossover frequency is more than sufficient to put the computed moduli on realistic time and energy scales. Importantly, the computations in concert with microrheology measurements, deployed across a range of systems, have the promise of enhancing our understanding of the connections between driving forces for condensate formation and the material properties of condensates.

ACKNOWLEDGMENTS

This work was supported by the St. Jude Research collaborative on the Biology and Biophysics of RNP granules (P.R.B. and R.V.P.), the US Air Force Office of Scientific Research (Grant No. FA9550-20-1-0241 to R.V.P.), and the US National Institutes of Health (No. R01NS121114 to R.V.P. and the training Grant No. T32 EB028092 that provided partial support to S.R.C.). We thank Mina Farag for the simulations of the A1-LCD systems and Ibraheem Alshareedah and Tanja Mittag for helpful discussions. The code for computing viscoelastic properties from graph Laplacians matrices is available on GitHub at https://github.com/Pappulab/material_properties.

AUTHOR DECLARATIONS

Conflict of Interest

The authors have no conflicts to disclose.

Author Contributions

Samuel R. Cohen: Conceptualization (equal); Investigation (equal); Methodology (equal); Software (equal); Validation (equal); Visualization (equal); Writing – original draft (equal); Writing – review & editing (equal). Priya R. Banerjee: Data curation (supporting); Funding acquisition (supporting); Investigation (supporting); Methodology (supporting); Project administration (supporting); Writing – review & editing (supporting). Rohit V. Pappu: Conceptualization (equal); Formal analysis (equal); Funding acquisition (equal); Investigation (equal); Methodology (equal); Project administration (equal); Resources (equal); Supervision (equal); Validation (equal); Visualization (equal); Writing – original draft (equal); Writing – review & editing (equal).

DATA AVAILABILITY

The data that support the findings of this study are available from the corresponding author upon reasonable request.
==== Refs
REFERENCES

1. S. F. Banani, H. O. Lee, A. A. Hyman, and M. K. Rosen, “Biomolecular condensates: Organizers of cellular biochemistry,” Nat. Rev. Mol. Cell Biol. 18 , 285–298 (2017).10.1038/nrm.2017.7 28225081
2. Y. Shin and C. P. Brangwynne, “Liquid phase condensation in cell physiology and disease,” Science 357 , eaaf4382 (2017).10.1126/science.aaf4382 28935776
3. R. V. Pappu, S. R. Cohen, F. Dar, M. Farag, and M. Kar, “Phase transitions of associative biomacromolecules,” Chem. Rev. 123 , 8945–8987 (2023).10.1021/acs.chemrev.2c00814 36881934
4. J.-M. Choi, A. S. Holehouse, and R. V. Pappu, “Physical principles underlying the complex biology of intracellular phase transitions,” Annu. Rev. Biophys. 49 , 107–133 (2020).10.1146/annurev-biophys-121219-081629 32004090
5. J. Berry, C. P. Brangwynne, and M. Haataja, “Physical principles of intracellular organization via active and passive phase transitions,” Rep. Prog. Phys. 81 , 046601 (2018).10.1088/1361-6633/aaa61e 29313527
6. C. P. Brangwynne, P. Tompa, and R. V. Pappu, “Polymer physics of intracellular phase transitions,” Nat. Phys. 11 , 899–904 (2015).10.1038/nphys3532
7. H. Walter and D. E. Brooks, “Phase separation in cytoplasm, due to macromolecular crowding, is the basis for microcompartmentation,” FEBS Lett. 361 , 135–139 (1995).10.1016/0014-5793(95)00159-7 7698310
8. M. R. King, K. M. Ruff, A. Z. Lin, A. Pant, M. Farag, J. M. Lalmansingh, T. Wu, M. J. Fossat, W. Ouyang, M. D. Lew et al. , “Macromolecular condensation organizes nucleolar sub-phases to set up a pH gradient,” Cell 187 , 1889–1906.e24 (2024).10.1016/j.cell.2024.02.029 38503281
9. M. R. King, K. M. Ruff, and R. V. Pappu, “Emergent microenvironments of nucleoli,” Nucleus 15 , 2319957 (2024).10.1080/19491034.2024.2319957 38443761
10. I. Alshareedah, M. M. Moosa, M. Pham, D. A. Potoyan, and P. R. Banerjee, “Programmable viscoelasticity in protein-RNA condensates with disordered sticker-spacer polypeptides,” Nat. Commun. 12 , 6620 (2021).10.1038/s41467-021-26733-7 34785657
11. I. Alshareedah, G. M. Thurston, and P. R. Banerjee, “Quantifying viscosity and surface tension of multicomponent protein-nucleic acid condensates,” Biophys. J. 120 , 1161–1169 (2021).10.1016/j.bpj.2021.01.005 33453268
12. I. Alshareedah, A. Singh, S. Yang, V. Ramachandran, A. Quinn, D. A. Potoyan, and P. R. Banerjee, “Determinants of viscoelasticity and flow activation energy in biomolecular condensates,” Sci. Adv. 10 , eadi6539 (2024).10.1126/sciadv.adi6539 38363841
13. M. T. Wei, S. Elbaum-Garfinkle, A. S. Holehouse, C. C. Chen, M. Feric, C. B. Arnold, R. D. Priestley, R. V. Pappu, and C. P. Brangwynne, “Phase behaviour of disordered proteins underlying low density and high permeability of liquid organelles,” Nat. Chem. 9 , 1118–1125 (2017).10.1038/nchem.2803 29064502
14. N. O. Taylor, M. T. Wei, H. A. Stone, and C. P. Brangwynne, “Quantifying dynamics in phase-separated condensates using fluorescence recovery after photobleaching,” Biophys. J. 117 , 1285–1300 (2019).10.1016/j.bpj.2019.08.030 31540706
15. S. Elbaum-Garfinkle, Y. Kim, K. Szczepaniak, C. C.-H. Chen, C. R. Eckmann, S. Myong, and C. P. Brangwynne, “The disordered P granule protein LAF-1 drives phase separation into droplets with tunable viscosity and dynamics,” Proc. Natl. Acad. Sci. U. S. A. 112 , 7189–7194 (2015).10.1073/pnas.1504822112 26015579
16. H. Zhang, S. Elbaum-Garfinkle, E. M. Langdon, N. Taylor, P. Occhipinti, A. A. Bridges, C. P. Brangwynne, and A. S. Gladfelter, “RNA controls PolyQ protein phase transitions,” Mol. Cell 60 , 220–230 (2015).10.1016/j.molcel.2015.09.017 26474065
17. A. Ghosh, D. Kota, and H.-X. Zhou, “Shear relaxation governs fusion dynamics of biomolecular condensates,” Nat. Commun. 12 , 5995 (2021).10.1038/s41467-021-26274-z 34645832
18. D. Kota and H.-X. Zhou, “Macromolecular regulation of the material properties of biomolecular condensates,” J. Phys. Chem. Lett. 13 , 5285–5290 (2022).10.1021/acs.jpclett.2c00824
19. Y. Shen, F. S. Ruggeri, D. Vigolo, A. Kamada, S. Qamar, A. Levin, C. Iserman, S. Alberti, P. S. George-Hyslop, and T. P. J. Knowles, “Biomolecular condensates undergo a generic shear-mediated liquid-to-solid transition,” Nat. Nanotechnol. 15 , 841–847 (2020).10.1038/s41565-020-0731-4 32661370
20. S. Boeynaems, A. S. Holehouse, V. Weinhardt, D. Kovacs, J. Van Lindt, C. Larabell, L. Van Den Bosch, R. Das, P. S. Tompa, R. V. Pappu, and A. D. Gitler, “Spontaneous driving forces give rise to protein–RNA condensates with coexisting phases and complex material properties,” Proc. Natl. Acad. Sci. U. S. A. 116 , 7889–7898 (2019).10.1073/pnas.1821038116 30926670
21. C. W. Pak, M. Kosno, A. S. Holehouse, S. B. Padrick, A. Mittal, R. Ali, A. A. Yunus, D. R. Liu, R. V. Pappu, and M. K. Rosen, “Sequence determinants of intracellular phase separation by complex coacervation of a disordered protein,” Mol. Cell 63 , 72–85 (2016).10.1016/j.molcel.2016.05.042 27392146
22. R. S. Fisher and S. Elbaum-Garfinkle, “Tunable multiphase dynamics of arginine and lysine liquid condensates,” Nat. Commun. 11 , 4628 (2020).10.1038/s41467-020-18224-y 32934220
23. S. Roberts, T. S. Harmon, J. L. Schaal, V. Miao, K. Li, A. Hunt, Y. Wen, T. G. Oas, J. H. Collier, R. V. Pappu, and A. Chilkoti, “Injectable tissue integrating networks from recombinant polypeptides with tunable order,” Nat. Mater. 17 , 1154–1163 (2018).10.1038/s41563-018-0182-6 30323334
24. L. P. Bergeron-Sandoval, S. Kumar, H. K. Heris, C. L. A. Chang, C. E. Cornell, S. L. Keller, P. Francois, A. G. Hendricks, A. J. Ehrlicher, R. V. Pappu, and S. W. Michnick, “Endocytic proteins with prion-like domains form viscoelastic condensates that enable membrane remodeling,” Proc. Natl. Acad. Sci. U. S. A. 118 , e2113789118 (2021).10.1073/pnas.2113789118 34887356
25. P. E. Rouse, “A theory of the linear viscoelastic properties of dilute solutions of coiling polymers,” J. Chem. Phys. 21 , 1272–1280 (1953).10.1063/1.1699180
26. I. Alshareedah, W. M. Borcherds, S. R. Cohen, A. Singh, A. E. Posey, M. Farag, A. Bremer, G. W. Strout, D. T. Tomares, R. V. Pappu et al. , “Sequence-specific interactions determine viscoelasticity and aging dynamics of protein condensates,” Nat. Phys. (in press) (2024).10.1038/s41567-024-02558-1
27. A. Molliex, J. Temirov, J. Lee, M. Coughlin, A. P. Kanagaraj, H. J. Kim, T. Mittag, and J. P. Taylor, “Phase separation by low complexity domains promotes stress granule assembly and drives pathological fibrillization,” Cell 163 , 123–133 (2015).10.1016/j.cell.2015.09.015 26406374
28. E. W. Martin, A. S. Holehouse, I. Peran, M. Farag, J. J. Incicco, A. Bremer, C. R. Grace, A. Soranno, R. V. Pappu, and T. Mittag, “Valence and patterning of aromatic residues determine the phase behavior of prion-like domains,” Science 367 , 694–699 (2020).10.1126/science.aaw8653 32029630
29. A. Bremer, M. Farag, W. M. Borcherds, I. Peran, E. W. Martin, R. V. Pappu, and T. Mittag, “Deciphering how naturally occurring sequence features impact the phase behaviours of disordered prion-like domains,” Nat. Chem. 14 , 196–207 (2022).10.1038/s41557-021-00840-w 34931046
30. N. Taylor, S. Elbaum-Garfinkle, N. Vaidya, H. Zhang, H. A. Stone, and C. P. Brangwynne, “Biophysical characterization of organelle-based RNA/protein liquid phases using microfluidics,” Soft Matter 12 , 9142–9150 (2016).10.1039/C6SM01087C 27791212
31. J. P. Brady, P. J. Farber, A. Sekhar, Y.-H. Lin, R. Huang, A. Bah, T. J. Nott, H. S. Chan, A. J. Baldwin, J. D. Forman-Kay, and L. E. Kay, “Structural and hydrodynamic properties of an intrinsically disordered region of a germ cell-specific protein on phase separation,” Proc. Natl. Acad. Sci. U. S. A. 114 , E8194–E8203 (2017).10.1073/pnas.1706197114 28894006
32. T. Mittag and R. V. Pappu, “A conceptual framework for understanding phase separation and addressing open questions and challenges,” Mol. Cell 82 , 2201–2214 (2022).10.1016/j.molcel.2022.05.018 35675815
33. J. Wang, J.-M. Choi, A. S. Holehouse, H. O. Lee, X. Zhang, M. Jahnel, S. Maharana, R. Lemaitre, A. Pozniakovsky, D. Drechsel et al. , “A molecular grammar governing the driving forces for phase separation of prion-like RNA binding proteins,” Cell 174 , 688–699.e16 (2018).10.1016/j.cell.2018.06.006 29961577
34. J.-M. Choi, F. Dar, and R. V. Pappu, “Lassi: A lattice model for simulating phase transitions of multivalent proteins,” PLoS Comput. Biol. 15 , e1007028 (2019).10.1371/journal.pcbi.1007028 31634364
35. J.-M. Choi, A. A. Hyman, and R. V. Pappu, “Generalized models for bond percolation transitions of associative polymers,” Phys. Rev. E 102 , 042403 (2020).10.1103/PhysRevE.102.042403 33212590
36. R. B. Davis, T. Kaur, M. M. Moosa, and P. R. Banerjee, “FUS oncofusion protein condensates recruit mSWI/SNF chromatin remodeler via heterotypic interactions between prion-like domains,” Protein Sci. 30 , 1454–1466 (2021).10.1002/pro.4127 34018649
37. G. M. Wadsworth, W. J. Zahurancik, X. Zeng, P. Pullara, L. B. Lai, V. Sidharthan, R. V. Pappu, V. Gopalan, and P. R. Banerjee, “RNAs undergo phase transitions with lower critical solution temperatures,” Nat. Chem. 15 , 1693–1704 (2023).10.1038/s41557-023-01353-4 37932412
38. Y. Dai, M. Farag, D. Lee, X. Zeng, K. Kim, H.-i. Son, X. Guo, J. Su, N. Peterson, J. Mohammed et al. , “Programmable synthetic biomolecular condensates for cellular control,” Nat. Chem. Biol. 19 , 518–528 (2023).10.1038/s41589-022-01252-8 36747054
39. M. Kar, F. Dar, T. J. Welsh, L. T. Vogel, R. Kühnemuth, A. Majumdar, G. Krainer, T. M. Franzmann, S. Alberti, C. A. M. Seidel et al. , “Phase-separating RNA-binding proteins form heterogeneous distributions of clusters in subsaturated solutions,” Proc. Natl. Acad. Sci. U. S. A. 119 , e2202222119 (2022).10.1073/pnas.2202222119 35787038
40. C. Lan, J. Kim, S. Ulferts, F. Aprile-Garcia, S. Weyrauch, A. Anandamurugan, R. Grosse, R. Sawarkar, A. Reinhardt, and T. Hugel, “Quantitative real-time in-cell imaging reveals heterogeneous clusters of proteins prior to condensation,” Nat. Commun. 14 , 4831 (2023).10.1038/s41467-023-40540-2 37582808
41. M. Kar, L. T. Vogel, G. Chauhan, S. Felekyan, H. Ausserwöger, T. J. Welsh, F. Dar, A. R. Kamath, T. P. J. Knowles, A. A. Hyman et al. , “Solutes unmask differences in clustering versus phase separation of FET proteins,” Nat. Commun. 15 , 4408 (2024).10.1038/s41467-024-48775-3 38782886
42. E. W. Martin and T. Mittag, “Relationship of sequence and phase separation in protein low-complexity regions,” Biochemistry 57 , 2478–2487 (2018).10.1021/acs.biochem.8b00008 29517898
43. M. Kar, A. E. Posey, F. Dar, A. A. Hyman, and R. V. Pappu, “Glycine-rich peptides from FUS have an intrinsic ability to self-assemble into fibers and networked fibrils,” Biochemistry 60 , 3213–3222 (2021).10.1021/acs.biochem.1c00501 34648275
44. B. E. Eichinger, “Configuration statistics of Gaussian molecules,” Macromolecules 13 , 1–11 (1980).10.1021/ma60073a001
45. M. Farag, S. R. Cohen, W. M. Borcherds, A. Bremer, T. Mittag, and R. V. Pappu, “Condensates formed by prion-like low-complexity domains have small-world network structures and interfaces defined by expanded conformations,” Nat. Commun. 13 , 7722 (2022).10.1038/s41467-022-35370-7 36513655
46. I. Carmesin and K. Kremer, “The bond fluctuation method: A new effective algorithm for the dynamics of polymers in all spatial dimensions,” Macromolecules 21 , 2819–2823 (1988).10.1021/ma00187a030
47. J. S. Shaffer, “Effects of chain topology on polymer dynamics: Bulk melts,” J. Chem. Phys. 101 , 4205–4213 (1994).10.1063/1.467470
48. N. Metropolis, A. W. Rosenbluth, M. N. Rosenbluth, A. H. Teller, and E. Teller, “Equation of state calculations by fast computing machines,” J. Chem. Phys. 21 , 1087–1092 (1953).10.1063/1.1699114
49. G. D. J. Phillies, “Linear response theory,” in Elementary Lectures in Statistical Mechanics (Springer, New York, 2000), pp. 365–371.10.1007/978-1-4612-1264-5_33
50. M. E. J. Newman, Networks: An Introduction (Oxford University Press, 2010).
51. B. H. Zimm, “Dynamics of polymer molecules in dilute solution: Viscoelasticity, flow birefringence and dielectric loss,” J. Chem. Phys. 24 , 269–278 (1956).10.1063/1.1742462
52. W. L. Peticolas, “Introduction to the molecular viscoelastic theory of polymers and its applications,” Rubber Chem. Technol. 36 , 1422–1458 (1963).10.5254/1.3539650
53. F. R. K. Chung, Spectral Graph Theory (American Mathematical Society, 1994).
54. M. Doi and S. F. Edwards, The Theory of Polymer Dynamics (Oxford University Press, 1988).
55. J. K. Wróbel, R. Cortez, and L. Fauci, “Modeling viscoelastic networks in Stokes flow,” Phys. Fluids 26 , 113102 (2014).10.1063/1.4900941
56. M. Rubinstein and R. H. Colby, Polymer Physics (Oxford University Press, 2003).
57. T. S. Harmon, A. S. Holehouse, M. K. Rosen, and R. V. Pappu, “Intrinsically disordered linkers determine the interplay between phase separation and gelation in multivalent proteins,” eLife 6 , 30294 (2017).10.7554/eLife.30294
58. C. Song, S. Havlin, and H. A. Makse, “Self-similarity of complex networks,” Nature 433 , 392–395 (2005).10.1038/nature03248 15674285
59. N. Galvanetto, M. T. Ivanović, A. Chowdhury, A. Sottini, M. F. Nüesch, D. Nettels, R. B. Best, and B. Schuler, “Extreme dynamics in a biomolecular condensate,” Nature 619 , 876–883 (2023).10.1038/s41586-023-06329-5 37468629
60. A. Soranno, B. Buchli, D. Nettels, R. R. Cheng, S. Müller-Späth, S. H. Pfeil, A. Hoffmann, E. A. Lipman, D. E. Makarov, and B. Schuler, “Quantifying internal friction in unfolded and intrinsically disordered proteins with single-molecule spectroscopy,” Proc. Natl. Acad. Sci. U. S. A. 109 , 17800–17806 (2012).10.1073/pnas.1117368109 22492978
61. L. Jawerth, E. Fischer-Friedrich, S. Saha, J. Wang, T. Franzmann, X. Zhang, J. Sachweh, M. Ruer, M. Ijavi, S. Saha et al. , “Protein condensates as aging Maxwell fluids,” Science 370 , 1317–1323 (2020).10.1126/science.aaw4951 33303613
62. A. Patel, H. O. Lee, L. Jawerth, S. Maharana, M. Jahnel, M. Y. Hein, S. Stoynov, J. Mahamid, S. Saha, T. M. Franzmann et al. , “A liquid-to-solid phase transition of the ALS protein FUS accelerated by disease mutation,” Cell 162 , 1066–1077 (2015).10.1016/j.cell.2015.07.047 26317470
63. S. Biswas and D. A. Potoyan, “Molecular drivers of aging in biomolecular condensates: Desolvation, rigidification, and sticker lifetimes,” PRX Life 2 , 023011 (2024).10.1103/PRXLife.2.023011
64. A. A. Hyman, C. A. Weber, and F. Jülicher, “Liquid-liquid phase separation in biology,” Annu. Rev. Cell Dev. Biol. 30 , 39–58 (2014).10.1146/annurev-cellbio-100913-013325 25288112
65. D. Sundaravadivelu Devarajan, J. Wang, B. Szała-Mendyk, S. Rekhi, A. Nikoubashman, Y. C. Kim, and J. Mittal, “Sequence-dependent material properties of biomolecular condensates and their relation to dilute phase conformations,” Nat. Commun. 15 , 1912 (2024).10.1038/s41467-024-46223-w 38429263
66. F. Yuan, H. Alimohamadi, B. Bakka, A. N. Trementozzi, K. J. Day, N. L. Fawzi, P. Rangamani, and J. C. Stachowiak, “Membrane bending by protein phase separation,” Proc. Natl. Acad. Sci. U. S. A. 118 , e2017435118 (2021).10.1073/pnas.2017435118 33688043
67. F. Yuan, C. T. Lee, A. Sangani, J. R. Houser, L. Wang, E. M. Lafer, P. Rangamani, and J. C. Stachowiak, “The ins and outs of membrane bending by intrinsically disordered proteins,” Sci. Adv. 9 , eadg3485 (2023).10.1126/sciadv.adg3485 37418523
68. M. L. Gardel, K. E. Kasza, C. P. Brangwynne, J. Liu, and D. A. Weitz, “Chapter 19 mechanical response of cytoskeletal networks,” in Methods in Cell Biology (Academic Press, 2008), pp. 487–519.
69. Y. Shin, Y.-C. Chang, D. S. W. Lee, J. Berry, D. W. Sanders, P. Ronceray, N. S. Wingreen, M. Haataja, and C. P. Brangwynne, “Liquid nuclear condensates mechanically sense and restructure the genome,” Cell 175 , 1481–1491.e13 (2018).10.1016/j.cell.2018.10.057 30500535
70. A. R. Strom, R. J. Biggs, E. J. Banigan, X. Wang, K. Chiu, C. Herman, J. Collado, F. Yue, J. C. Ritland Politz, L. J. Tait et al. , “HP1α is a chromatin crosslinker that controls nuclear and mitotic chromosome mechanics,” eLife 10 , e63972 (2021).10.7554/eLife.63972 34106828
71. B. Gouveia, Y. Kim, J. W. Shaevitz, S. Petry, H. A. Stone, and C. P. Brangwynne, “Capillary forces generated by biomolecular condensates,” Nature 609 , 255–264 (2022).10.1038/s41586-022-05138-6 36071192
72. W. M. Lai, D. Rubin, and E. Krempl, “Chapter 8—Non-Newtonian fluids,” in Introduction to Continuum Mechanics, 4th ed., edited by W. M. Lai, D. Rubin, and E. Krempl (Butterworth-Heinemann, 2010), pp. 443–509.
73. R. Xiao, H. Sun, and W. Chen, “An equivalence between generalized Maxwell model and fractional Zener model,” Mech. Mater. 100 , 148–153 (2016).10.1016/j.mechmat.2016.06.016
74. E. Wiechert, “Gesetze der elastischen Nachwirkung für constante temperatur,” Ann. Phys. 286 , 335–348 (1893).10.1002/andp.18932861011
75. J. Honerkamp and J. Weese, “A nonlinear regularization method for the calculation of relaxation spectra,” Rheol. Acta 32 , 65–73 (1993).10.1007/BF00396678
76. A. Takeh and S. Shanbhag, “A computer program to extract the continuous and discrete relaxation spectra from dynamic viscoelastic measurements,” Appl. Rheol. 23 , 24628 (2013).10.3933/applrheol-23-24628
77. T. G. Mason and D. A. Weitz, “Optical measurements of frequency-dependent linear viscoelastic moduli of complex fluids,” Phys. Rev. Lett. 74 , 1250–1253 (1995).10.1103/PhysRevLett.74.1250 10058972
78. T. G. Mason, K. Ganesan, J. H. van Zanten, D. Wirtz, and S. C. Kuo, “Particle tracking microrheology of complex fluids,” Phys. Rev. Lett. 79 , 3282–3285 (1997).10.1103/PhysRevLett.79.3282
79. G. Tesei, T. K. Schulze, R. Crehuet, and K. Lindorff-Larsen, “Accurate model of liquid–liquid phase behavior of intrinsically disordered proteins from optimization of single-chain properties,” Proc. Natl. Acad. Sci. U. S. A. 118 , e2111696118 (2021).10.1073/pnas.2111696118 34716273
80. J. A. Joseph, A. Reinhardt, A. Aguirre, P. Y. Chew, K. O. Russell, J. R. Espinosa, A. Garaizar, and R. Collepardo-Guevara, “Physics-driven coarse-grained model for biomolecular phase separation with near-quantitative accuracy,” Nat. Comput. Sci. 1 , 732–743 (2021).10.1038/s43588-021-00155-3 35795820
81. G. L. Dignon, W. Zheng, R. B. Best, Y. C. Kim, and J. Mittal, “Relation between single-molecule properties and phase behavior of intrinsically disordered proteins,” Proc. Natl. Acad. Sci. U. S. A. 115 , 9929–9934 (2018).10.1073/pnas.1804177115 30217894
82. J. Wang, D. S. Devarajan, A. Nikoubashman, and J. Mittal, “Conformational properties of polymers at droplet interfaces as model systems for disordered proteins,” ACS Macro Lett. 12 , 1472–1478 (2023).10.1021/acsmacrolett.3c00456 37856873
83. M. Farag, W. M. Borcherds, A. Bremer, T. Mittag, and R. V. Pappu, “Phase separation of protein mixtures is driven by the interplay of homotypic and heterotypic interactions,” Nat. Commun. 14 , 5527 (2023).10.1038/s41467-023-41274-x 37684240
84. I. Alshareedah, M. M. Moosa, M. Raju, D. A. Potoyan, and P. R. Banerjee, “Phase transition of RNA–protein complexes into ordered hollow condensates,” Proc. Natl. Acad. Sci. U. S. A. 117 , 15650–15658 (2020).10.1073/pnas.1922365117 32571937
85. A. S. Holehouse, G. M. Ginell, D. Griffith, and E. Böke, “Clustering of aromatic residues in prion-like domains can tune the formation, state, and organization of biomolecular condensates,” Biochemistry 60 , 3566–3581 (2021).10.1021/acs.biochem.1c00465 34784177
86. R. M. Welles, K. A. Sojitra, M. V. Garabedian, B. Xia, W. Wang, M. Guan, R. M. Regy, E. R. Gallagher, D. A. Hammer, J. Mittal, and M. C. Good, “Determinants that enable disordered protein assembly into discrete condensed phases,” Nat. Chem. 16 , 1062 (2024).10.1038/s41557-023-01423-7 38316988
