Protein aging dynamics: A perspective from non-equilibrium coarse-grained models

Cellular compartments without membranes, such as stress granules,^[1] Cajal bodies,^[2] and nucleoli,^[3] play essential roles in regulating intracellular biochemical activities.^[4–6] For example, stress granules form through LLPS in response to cellular stress,^[7] while the phase separation of ribosomes and nucleoli is involved in the regulation of transcription and translation.^[8,9] In recent years, it has been widely recognized that these membrane-less organelles form through liquid–liquid phase separation (LLPS).^[10–12] During LLPS, biomolecules aggregate via weak multivalent interactions to form phase-separated condensates.^[13,14] These interactions include, but are not limited to, hydrophobic interactions,^[15] electrostatic interactions,^[16] hydrogen bonds,^[17] and π–π stacking.^[18] Through these interactions, molecules can condense into highly concentrated phase-separated regions under specific environmental conditions, facilitating the efficient execution of various biochemical processes within the cell.

However, the lifespan of cells is limited, and over time, proteins may undergo misfolding and aberrant aggregation. Such misfolding and aggregation not only disrupt normal biological functions but also lead to the destabilization of membrane-less organelle structures,^[19] further accelerating cellular aging. This process is closely associated with various neurodegenerative diseases,^[20–22] such as Alzheimer’s disease and Parkinson’s disease. Specifically, within membrane-less organelles, certain proteins, such as FUS^[23] and TDP-43,^[24] when abnormally aggregated, can significantly impair cellular function and become critical factors in disease pathogenesis. During the aging process of biomolecular condensates, their physical properties may transition from a fluid-like state with high mobility to a more viscous gel-like or even solid state.^[25,26] These studies suggest that this liquid-to-solid transition plays a pivotal role in the hardening process of protein aggregation.

In existing studies, the aging of condensates is typically categorized into gelation and vitrification processes. Misfolded proteins may form condensate nuclei, which can further develop into insoluble amyloid aggregates. Mutations in the FUS protein and other RNA-binding proteins (RBPs), associated with amyotrophic lateral sclerosis (ALS) and Alzheimer’s disease, have been found to accelerate this hardening process.^[25] While the association between amyloid fibrils and diseases has long been recognized, advances in structural biology techniques, such as cryo-electron microscopy (cryo-EM) and solid-state nuclear magnetic resonance (NMR) spectroscopy, have provided deeper insights into the relationship between fibril structures and disease mechanisms. For example, Ray et al. identified domain-level interactions that drive the liquid-to-solid transition of α-synuclein.^[27] Pantoja-Uceda et al. investigated the process by which the low-complexity domain (LCD) of TDP-43 transitions from a liquid state to an amyloid structure under low pH (salt-free or RNA-free) conditions, revealing distinct roles of the LCD region in the liquid-to-solid transition.^[28] Using fluorescence confocal microscopy, they observed that the formation of amyloid fibrils begins at the interface of liquid droplets.^[29] These experimental studies provide critical insights into the dynamic mechanisms underlying cellular aging processes.

In recent years, computational simulations have become a significant complement to experimental techniques and have made remarkable progress in the high-resolution study of protein dynamics and self-assembly.^[30] Recent studies have shown that all-atom simulations investigating the formation of Aβ dimers in solution and neuron models provide a theoretical basis for understanding oligomer-mediated neurotoxicity mechanisms.^[31] Coarse-grained molecular dynamics simulations, due to their ability to simplify complex structures while retaining core structural features, have been widely applied in studies of cellular aging processes.^[32,33] Various types of sticker-and-spacer models have been developed and have demonstrated varying degrees of success in the study of folding proteins and linear multivalent proteins.^[34,35] For instance, coarse-grained simulations of FUS-PLD condensates have shown that homotypic and heterotypic interactions can control the morphology of condensates.^[36] By regulating the interaction strength between domains, the morphology of the condensates can be transitioned from a disordered mixed assembly to an ordered, core-shell structured fibrous assembly.^[33] Moreover, coarse-grained models have also been used to study solvent exclusion effects, the rigidity and flexibility of biopolymer scaffolds, and the impact of sticker contact lifetimes in proteins on the viscoelastic properties and aging dynamics of condensates.^[37]

In protein aging studies, Garaizar et al. developed a time-dependent coarse-grained algorithm capable of dynamically simulating the structural diversity during the aging process of proteins.^[38] This method revealed how completely disordered PLD structures spontaneously transition into ordered structures, and through the study of a coarse-grained model based on FUS proteins, it explored the impact of aging on the formation of the multicomponent structures of condensates and the temporal changes in their material properties. This dynamic aging algorithm has now been applied to various studies investigating protein aging processes.^[39–41] However, the algorithm still has limitations. For example, it fails to predict the amyloid fibril structures.^[42,43] Therefore, we have improved the existing model by replacing the single coarse-grained aging particles with patchy particles. Two binding sites were introduced at two sides of the aging particles, with the binding directions partially constrained, thereby facilitating the occurrence of a chaining effect in the aging particles. The layout of this paper is as follows. Section 2 introduces the model and simulation methods. Section 3 provides an in-depth analysis of the structural and dynamic changes in condensates pre- and post-aging. Section 4 offers a comprehensive summary of the study and outlines prospects for future research.

In this study, we coarse-grained proteins into polymer chains consisting of 50 beads. Based on the sticker-and-spacer model framework, we assign two types of sticker beads to each chain, with three beads of each type, denoted by green and yellow, respectively (see Fig. 1). The remaining beads serve as spacer beads connecting the sticker beads and are represented in pink. The distribution of the sticker beads (denoted as A and B hereafter) along the polymer chain follows a regular pattern, where a sticker bead appears every seven beads in the chain. In each polymer chain, we position a sticker bead every 7 particles, with possible aging patchy particles located at the fourth particle of each chain. The model setup is consistent with the characteristic distribution and ratio of tyrosine and arginine residues in the actual FUS protein, which exhibits a diblock-like distribution. The model has been successfully employed to predict FUS protein phase separation in our previous work.^[44] Over time, localized high-concentration regions emerge within protein aggregates, creating conditions conducive to aging. Low-complexity aromatic-rich kinked segments (LARKS) have been identified as aging-related motif, primarily located in the N-terminal and central regions of FUS.^[45] To model the aging process of LARKS, we set the condition that aging occurs when four aging particles aggregate and their distance is less than 2.0 nm. This parameter choice is based on previous simulation studies,^[32] which indicated that four short peptides represent the minimal unit for modeling the aging process. Upon aging, the particles change color from pink to sky blue; we present a schematic representation of the coarse-grained model after aging, as shown in Fig. 1.

In the coarse-grained protein solution, the primary driving force for phase separation is the weak multivalent interactions between the two types of stickers.^[44] In the coarse-grained simulations, the protein solution is represented as linear chains embedded in an implicit solvent environment. Each pair of beads is connected by an elastic harmonic potential, with the potential energy function expressed as follows:

The bond length r_b = 4.5 nm and the bond stiffness k=4.0899ε0/rb2, where ε₀ is the energy unit. To simulate the specific one-to-one interaction between the two types of stickers, an attractive interaction potential is introduced between sticker-A and sticker-B, which is expressed as follows:^[44,46]

Here, U₀ = −20.7ε₀, the cutoff radius r₀ = 2.0 nm, and the energy unit is ε₀. To avoid the situation where multiple beads are simultaneously connected to a single bead, a softened cutoff Lennard–Jones potential is introduced between sticker-A and sticker-B, expressed as

where ε = 0.621ε₀, λ = 0.68, σ = 3.5 nm, and the cutoff radius r_c = 5 nm. Additionally, during the early stage of phase separation, the aging process is not involved; thus, the interaction between the aging beads, the orientational beads, and the spacer beads is identical. A uniform repulsive potential is applied to spacer beads and their interactions with other beads.^[14,46] The only difference is that the energy unit of the repulsive potential is modified to ε = 0.00001ε₀, with all other parameters remaining unchanged.

To facilitate the formation of aging patchy particles in the simulation,^[47] we introduced an angular potential between the particles at the two binding sites and the aging particles. This angular potential ensures that during the aging process, the two binding sites remain positioned at the two ends of the aging particle, partially restricting the spatial binding locations of the aging particle. The functional form of the angular potential is as follows:^[38]

where K = 900 represents the angular stiffness, and θ₀ = 180° is the equilibrium angle. Once the aging particles are formed, the particles at the binding sites experience a strong intermolecular attraction. The interaction between the binding sites is set to be the same as the interaction potential between sticker-A and sticker-B, but with U₀ = −41.4ε₀, ensuring an irreversible binding. According to the above model, various types of binding may form between aging particles, as illustrated in Fig. 2, which could include complete linear binding of aging beads, pairwise binding, and combinations of single particles with chains of three beads. It is worth noting that the selection of parameters in our coarse-grained model primarily aligns with realistic physical scenarios. Specifically, before aging, the stickers form reversible bonds, whereas after aging, irreversible bonds are established between the aging beads. Given that the coarse-grained beads are larger than actual residues or motifs, it is challenging to quantitatively compare our simulation parameters with experimentally measured interaction parameters between residues or motifs. A quantitative comparison would require more precise all-atom simulations.

In the simulation, we established a model consisting of 500 coarse-grained protein chains, placed within a cubic simulation box with dimensions of 300 nm×300 nm×300 nm, and periodic boundary conditions applied in all directions. All simulations were conducted under the NVT ensemble, with the temperature controlled at 300 K using a Langevin thermostat. The mass of each bead was set to m = 3534.3 ag. Additionally, the simulations used a damping coefficient of τ = m/γ = 125 ns and a time step of dt = 2.5 ns. Simulation data were recorded every 2.5 ns, with a total simulation time of 25 ms. In the case of studying the internal diffusion with varying rates, we used a damping coefficient τ′ = m/γ = 1.5 ns for the internal molecules, while keeping the other parameters unchanged. These simulations were performed using the LAMMPS molecular dynamics software (version from June 23, 2023),^[48] and all analyses were based on data from five independent simulations.

Proteins driven by weak multivalent interactions undergo phase separation in the simulation, ultimately forming a coexistence of dense and dilute phases.^[49] The simulation results reveal that the generated condensates adopt a nearly spherical structure, as shown in Fig. 3(a). The black arrows indicate the diameter of the condensates. With the introduction of the aging mechanism, condensate size noticeably decreases, as shown in Fig. 3(b). We present snapshots illustrating the liquid–liquid phase separation and aging of the system in supporting information (SI). A comparison indicates that the diameter of the condensates after aging is significantly smaller than before, reflecting that the aging process causes the internal structure of the condensates to become more compact.

Figure 3(c) shows the snapshot of aging beads within the condensate post-aging. To further investigate the behavior of aging beads, we zoomed in on the interior of the condensate to focus on the binding of aging beads. The results indicate that aging beads primarily undergo linking, leading to the formation of more compact structures within the condensate. To quantify the impact of aging on the internal structure of the condensate, we compared the density distributions inside the condensates pre- and post-aging. Figure 3(d) displays the internal density distributions of the condensates, centered around their center of mass. We observed that, within a small distance range, the density distributions of both the pre- and post-aging condensates are parabolic, and have wide interface widths. Upon comparing the density distributions of the stabilized condensates, we found that the internal density of the condensate significantly increases after aging. This suggests that the aging process makes the condensate more compact, which may contribute to enhancing its stability. Additionally, we examined the density distribution of condensates containing only aging particles, as shown in the inset of Fig. 3(d). It was found that aging particles tend to concentrate in the interior of the condensate, with the density distribution gradually decreasing with increasing distance from the center.

To further understand the microscopic mechanisms behind the increased internal density of the condensate during aging, we analyzed the number of bonds formed between sticker-A and sticker-B on each protein chain. Since our phase separation system is primarily driven by weak multivalent interactions between sticker-A and sticker-B, as well as strong interactions between aging particles, we calculated the number of bonds formed between the two types of stickers on each protein chain pre- and post-aging, along with the number of bonds formed by aging particles, as shown in Fig. 3(e). The results indicate that over time, the number of bonds between aging particles gradually increases, which in turn promotes the increase in the number of bonds between sticker-A and sticker-B. This suggests that aging not only involves the linking of aging particles but also enhances the number of bonds formed due to weak multivalent interactions, further driving the increase in the internal density of the condensate.

Furthermore, in biological systems, condensates often exchange material with their external environment.^[13,50] In our simulation, we compared the number of chains exchanging with the external environment in two scenarios (see Fig. 3(f)). The results demonstrate that the aging process significantly reduces the exchange between the condensate and the external environment, a phenomenon consistent with observations from other research.^[20,51,52]

In many systems where phase separation occurs in proteins, the aging process of the condensate and the structures formed during aging have been the focus of research. Experimental studies have shown that proteins preferentially form fibrous structures on the surface of condensates.^[29,41] However, the mechanisms underlying the preferential aging phenomena at the surface of the condensates are not yet fully understood. In our simulations, density analysis revealed that the condensate’s interior exhibits a high and uniform density, while the exterior has a broader interface width with a corresponding decrease in density (see Fig. 3(d)). A recent study revealed that the uneven network structure formed by multivalent interactions within condensates significantly influences the distribution of acidic regions and the diffusion behavior of peptides.^[53] Furthermore, other research has shown that the highly crowded intracellular environment also plays a crucial role in affecting diffusion rates.^[54] Based on these observations, we propose the hypothesis that the high density within the condensate results in slower diffusion inside, while the broader interface on the exterior leads to faster diffusion at the surface. We hypothesize that the difference in diffusion rates between the interior and exterior of the condensate may drive the aging process to occur preferentially on the surface, promoting the formation of fibrous structures.

To test this hypothesis, we present schematic diagrams of two different aging models in Figs. 4(a) and 4(e). In Fig. 4(a), the diffusion rate within the condensate is homogeneous and represented by a uniform orange background. In this case, as the aging process progresses, Figs. 4(b)–4(d) show the distribution of aging particles within the condensate. The light blue particles represent those that have undergone aging. Over time, we observe that aging particles initially form within the interior of the condensate and gradually age outward toward the periphery. In this case, the number of aging particles achieves a uniform and random distribution throughout the condensate. For the model with heterogeneous diffusion rates, the dark orange regions represent areas with higher density and slower diffusion rates, while the lighter areas indicate regions with lower density and faster diffusion rates. We set the radius of the low diffusion region to be at a distance of 35 nm from the center of mass of the condensate. This value was determined based on the density distribution of the condensate after aging, where the density decreases to half of its maximum value (see Fig. 3(d)). Figures 4(f)–4(h) depict the spatial distribution of aging particles within the condensate when there is a heterogeneous diffusion rate. We observe that, as the simulation progresses, aging particles first concentrate at the outer layers of the condensate. This observation further supports our hypothesis that aging preferentially occurs at the surface of the condensate.

To quantitatively analyze the distribution of aging particles, we calculated the density distribution of aging particles under both homogeneous and heterogeneous diffusion rates. In Fig. 5(a), we present the density distribution of aging particles under homogeneous diffusion conditions. The results indicate that, under these conditions, aging particles are more likely to preferentially form in the interior of the condensate due to the high concentration. As the number of aging particles increases, the density within the condensate gradually rises. In Fig. 5(b), we show the situation when there are different diffusion rates within the condensate, where aging particles tend to concentrate at the outer layers. As time progresses, the number of aging particles increases and becomes more concentrated at the interface of the condensate. Due to the faster diffusion rate at the interface of the condensate, the probability of collisions between aging particles in the interfacial region increases, leading to more frequent collisions and aggregation of aging particles at the interface. This suggests that when the diffusion rates are unevenly distributed, aging particles are more likely to appear first on the surface of the condensate.

Furthermore, we calculated the probability distribution of the chains formed by aging particles after forming irreversible bonds for three different times. Figure 5(c) shows that, as time increases, aging particles are able to form more irreversible bonds, leading to the formation of longer chains with more aging particles. For example, it can be observed in Fig. 5(c) that 20 aging particles can link together to form irreversible bonds. In Fig. 5(d), due to the heterogeneous diffusion rate inside the condensate, shorter chains of aging particle aggregates are formed. This result indicates that the different diffusion rates within the condensate not only affect the distribution of aging particles but may also disrupt the chaining of aging particles.

In the above discussions, we compared the microscopic structural differences of the condensates pre- and post-aging. In this section, we further examine the macroscopic dynamic properties of the condensate, specifically its diffusion behavior and mechanical response.

Firstly, to assess the diffusion characteristics of the condensate, we calculated the mean squared displacement (MSD) under different conditions. Figure 6(a) shows the trend of MSD as a function of time for the condensate under three conditions. The results indicate that the aged condensates exhibit a lower diffusion rate. This is because the aging process makes condensates more compact. The originally loose network structure transforms into a dense network structure, significantly reducing the mobility of the condensates. To further investigate the impact of aging on the mechanical properties of the condensate, we employed the shear stress relaxation modulus G(t) to describe the mechanical response of the condensate.^[55,56] In systems with isotropic interactions, the shear stress relaxation modulus can be calculated using the six independent components of the pressure tensor. Its expression is shown in Refs. [57,58],

where σ_αβ(t) represents the off-diagonal components of the stress tensor, V denotes the volume, k_B is the Boltzmann constant, and T is the temperature. N_αβ = σ_αα − σ_ββ represents the difference in normal stress. In Fig. 6(b), we present the evolution of the shear stress relaxation modulus of the condensate over time under three conditions. The red curve represents the condensate pre-aging, which exhibits a relatively fast relaxation rate, indicating that the condensate prior to aging exhibits liquid-like behavior and can quickly reach mechanical equilibrium. In contrast, the other two curves correspond to the aged condensates, showing a slower decay rate of the shear stress relaxation modulus. This suggests that the aged condensate requires a longer timescale to relax into mechanical equilibrium. Our simulation results are highly consistent with previous simulation studies, further validating their reliability and accuracy.^[40] These results indicate that aging enhances the solid-like behavior of the condensate, and its mechanical response becomes slower compared to the pre-aged state.

To quantify the viscosity characteristics of the condensate under different conditions, we further calculated the shear viscosity η. Based on the shear stress relaxation modulus, we performed time integration and used the approach from the literature to estimate the viscosity.^[59] On short timescales, viscosity can be obtained through numerical integration; whereas on longer timescales, we fit the shear stress relaxation modulus using the Maxwell model and derive the viscosity through analytical integration. Therefore, the viscosity is obtained by adding two different effective phases^[40]

where η(t₀) refers to the viscosity calculated at short timescales, while GM(t)=∑i=1MGiexp(−t/τi) represents the long-timescale component fitted using the Maxwell model.^[55,60] Here, t₀ is defined as the time point at which all molecular oscillations within G(t) have decayed and the function transitions into a strictly monotonic decay. Using this method, we can obtain viscosity values under different simulation conditions. As shown in Fig. 6(c), the calculation results indicate that the viscosity of the condensate significantly increases after aging, which aligns with our expected conclusion that aging leads to an increase in viscosity.

In this study, we implemented a non-equilibrium aging algorithm by modeling aging particles as patchy particles with specific binding properties, capable of forming one-to-two bindings and creating linear structures. This approach enhances the description of the aging dynamics of protein condensates.

During aging, the morphology of condensates evolves from sparse spheres to more compact structures. At the macroscopic level, aging leads to a significant reduction in condensate volume and an increase in the number of reversible bonds between particles. Microscopically, the exchangeability of condensates with the external environment decreases, indicating enhanced condensate closure. We also investigated the distribution of aging particles, finding that heterogeneous diffusion rates promote the preferential accumulation of particles at the interface. Aging particles tend to appear first in regions with slower diffusion, suggesting a surface-preferential aging mechanism. Furthermore, the dynamic and mechanical properties of condensates before and after aging were analyzed. Aged condensates showed slower diffusion and higher viscosity, resulting from the liquid-to-solid transition. These findings indicate that aging alters both the microscopic structure and macroscopic rheological properties of condensates. In conclusion, this study provides a new perspective on condensate aging, suggesting that it preferentially occurs at the surface due to diffusion differences, offering valuable insights into the aging process at the molecular level.

In future studies, we will further investigate the orientation of aging particles and the rigidity of particles after binding, aiming to explore the formation of fibrous structures within the condensates and study the underlying dynamics of their formation. Recent experimental and simulation studies have revealed that mutations in proteins associated with ALS can enhance aggregation and lead to the formation of misfolded proteins, triggering disease onset.^[61,62] Drawing on these findings, our simulations can be extended to explore the connections between condensate aging, structural changes, and functional consequences. This approach will provide deeper insights into the mechanisms by which aging-induced alterations in condensates contribute to pathological protein aggregation and dysfunction, offering potential pathways for understanding and addressing neurodegenerative diseases.