Diffusion of intrinsically disordered proteins within viscoelastic membraneless droplets

Multivalent proteins and nucleic acids form biomolecular condensates, which are so-called membraneless organelles (MLOs), through phase separation Hyman et al. (2014); Banani et al. (2017); Pappu et al. (2023). A primary driving force for phase separation is the sequence-specific multivalent interactions among the components. Notably, intrinsically disordered proteins (IDPs), due to their flexibility and low-complexity sequences, play a crucial role in the formation of these condensates Protter et al. (2018); Wang et al. (2018). The organization of MLOs has been identified as contributing to essential biological processes, such as gene regulation, biochemical reactions, signal transduction, and molecular protection as a stress response. These functions depend on the unique physical properties of the condensates. The presence of flexible IDPs render the interior of the condensate highly dynamic and complex, exhibiting liquid-like properties with viscoelasticity Brangwynne et al. (2009); Elbaum-Garfinkle et al. (2015); Jawerth et al. (2018, 2020); Alshareedah et al. (2021).

Several studies have investigated the dynamic property of proteins within the condensates using experimental and computational methods such as pulsed field gradient NMR Murthy et al. (2019), fluorescence recovery after photobleaching (FRAP) Burke et al. (2015), single-molecule fluorescence microscopy Elbaum-Garfinkle et al. (2015); Fisher and Elbaum-Garfinkle (2020); Kamagata et al. (2021, 2022) and molecular dynamics (MD) simulations Zheng et al. (2020); Benayad et al. (2021); Tejedor et al. (2021, 2022); Ranganathan and Shakhnovich (2022). The diffusivity of proteins within condensates is significantly reduced compared to that in bulk and can be tuned by sequence variation and length. It has also been highlighted that the interface of biomolecular condensates acts as a key element in driving pathological protein aggregation Shen et al. (2023) and molecular exchange between the condensates and bulk Bo et al. (2021). A detailed exploration of the intricate diffusion process of IDPs using physical models could shed light on the properties pivotal for cellular processes. However, a comprehensive understanding of the spatiotemporal conformation and dynamics of IDPs within the condensates remains elusive at the molecular scale.

Here, we perform large-scale MD simulations of spherical droplets formed by mono IDPs to investigate the diffusive dynamics of proteins within the droplets. We demonstrate that proteins exhibit the anomalous diffusion with fluctuating diffusivity attributed to the viscoelastic environment of the droplet. Furthermore, we show that diffusivity varies between the interior and the interface of the droplets, resulting in a non-Gaussian distribution of displacement.

Subdiffusive motion of proteins within the droplet–We performed coarse-grained (CG) MD simulations of droplets formed by 1,000 IDPs, well-studied Fused in Sarcoma low-complexity domains (FUS-LCD), using a residue-revel CG model Dignon et al. (2018); Regy et al. (2021); Tesei et al. (2021) (see Fig. 1A and details in SI). Within the simulation timescales of $20\,\mathrm{\mu s}$, proteins diffuse between the interior and the interface of the droplet several times (Fig. Fig. 1B). To characterize the diffusive dynamics of proteins, the time-averaged mean squared displacement (TAMSD)

of proteins was calculated, where $t$ is the observation time, $\Delta$ is the lag time, and $\vec{r}_{i}$ is the relative coordinates of the center of mass (COM) of the $i$th protein obtained by subtracting the COM position of the droplet. In what follows, we exclusively analyzed proteins that reside within the droplet during the observation time. Figure 1C shows that TAMSDs of each protein exhibit sub-linear increase, $\overline{\delta^{2}_{i}(\Delta;t)}\propto\Delta^{\alpha}$, with $\alpha\sim 0.85$ within the droplets. After $1\,\mathrm{\mu s}$ the TAMSDs show a plateau, indicative of the spatial confinement of proteins within the droplet, where the plateau value of the TAMSDs corresponds to the square of the droplet radius, approximately $20\,\mathrm{nm}$. Note that the ensemble-averaged MSD is consistent with the TAMSDs, indicating that the diffusion process is ergodic (see Fig. S1). The subdiffusive behavior is distinct from the Brownian motion of an isolated protein in a solution, where the TAMSD increases linearly without saturation. The subdiffusive motions are often observed in the molecular diffusion within viscoelastic media, e.g. polymer solutions Mason and Weitz (1995); Mason et al. (1997) and lipid membranes Akimoto et al. (2011); Jeon et al. (2012); Yamamoto et al. (2014); Jeon et al. (2016). We also confirmed similar subdiffusion for another protein (DDX4), under different temperature conditions, and with different model parameters (see Fig. S2).

Anticorrelated motion and non-Gaussian distribution of protein movements–Fractional Brownian motion (FBM) Kolmogorov (1940); Mandelbrot and Van Ness (1968) is a fundamental model of anomalous diffusion characterized by Brownian motion driven by long-term correlated noise. If the noise exhibits anticorrelation, the TAMSD shows a sublinear increase Metzler et al. (2014). Figure 2A shows the normalized form of the displacement autocorrelation function (DAF) of proteins,

where $\langle\cdot\rangle$ indicates ensemble average and $\Delta=0.1\,\mathrm{ns}$. The transition of the DAF from positive to negative values indicates anticorrelated motion. The dip in the DAFs signifies the degree of anticorrelation. Viscoelastic effects are significant on a short timescale and diminish gradually over $\Delta=10^{-2}\,\mathrm{\mu s}$. Conversely, boundary confinement effects become noticeable around $\Delta=10^{-1}\,\mathrm{\mu s}$ and are more pronounced at $\Delta=1\,\mathrm{\mu s}$, a dip of sharp peak of $-0.5$ Weber et al. (2012), where the TAMSDs reach a plateau. Moreover, the normalized DAF exhibits a power-law decay, although it slightly deviates from the free FBM Burov et al. (2011), $C_{\Delta}(t)/C_{\Delta}(0)\propto t^{\alpha-2}$, where $\alpha$ is the same exponent as in $\overline{\delta^{2}_{i}(\Delta;t)}\propto\Delta^{\alpha}$ (Fig. 2A). As expected from the FBM with a DAF exponent of $-1.4$, the MSD should increase sublinearly with $\alpha=0.6$. The slight increase in the MSD exponent might be attributed to factors such as protein entanglements, conformational changes, etc.

To examine the Gaussianity of the displacement, we calculated a normalized propagator $P(\chi,\Delta)$, which represents the probability density function (PDF) of the displacement that the particle diffuses within a lag time $\Delta$. Since no systematic differences were observed in $x,y,z$ directions, we combined all increments from these directions into a single set denoted as $\chi$, i.e. displacements calculated as $x(t+\Delta)-x(t)$, $y(t+\Delta)-y(t)$, and $z(t+\Delta)-z(t)$ are collectively combined and represented as $\chi$. In the case of a normal FBM, the normalized propagator is expected to be Gaussian and consistent for different $\Delta$. Figure 2B shows that the propagators of proteins in the droplet exhibit non-Gaussianity with different $\Delta$ values. On a short time scale, until $\Delta=0.1\,\mathrm{\mu s}$, there are significant deviations from the Gaussian distribution in the tails of the distribution. Such non-Gaussianity has been observed in the diffusion process with fluctuating diffusivity Wang et al. (2012); Chubynsky and Slater (2014); He et al. (2016); Cherstvy and Metzler (2016); Chechkin et al. (2017), i.e., individual proteins within the droplet diffuse non-uniformly. Conversely, on a long time scale, the propagator is close to the Gaussian distribution. The $\Delta=0.1\,\mathrm{\mu s}$ is the same timescale at which the TAMSD converges to a plateau (see Fig. 1C), which again indicates the confinement effect in the droplet.

Fluctuating diffusivity of proteins within the droplet–To quantitatively evaluate the fluctuation of the diffusivity, we calculated the relative standard deviation (RSD) of TAMSD Uneyama et al. (2012); Yamamoto et al. (2017), the square root of the ergodicity breaking parameter  Burov et al. (2011),

In diffusion following a Gaussian process such as Brownian motion, the RSD decays as $t^{-0.5}$, but for non-ergodic processes such as CTRW He et al. (2008); Metzler et al. (2014) or with significant fluctuation in the diffusivity Akimoto and Yamamoto (2016a); Miyaguchi et al. (2016); Akimoto and Yamamoto (2016b), the RSD does not decay to zero even if observed for a long time. In polymeric solutions, the RSD changes from non-Gaussian fluctuations (RSD decays at $t^{-\alpha}(0<\alpha<0.5)$) to Gaussian fluctuations (RSD decays at $t^{-0.5}$) because of polymer entanglement, and the crossover time corresponds to the longest relaxation time Uneyama et al. (2012, 2015). The RSD of FUS-LCD exhibits slow decay at short timescales, then a crossover to Gaussian fluctuations (see Fig. 3).

Conformation and diffusivity of proteins at the interior and the interface of the droplets–From the radial density profile, we define the boundary between the interior and the interface of the droplet by $20\,\mathrm{nm}$ (see Fig. 4A). The residence time distributions for the interior and the interface of the droplet exhibit a power-law distribution with an exponential cutoff $P(t)\propto t^{-\beta}\exp{(-t/\tau)}$ (see Fig. 4B). The power-law exponent $\beta=-1.5$ is consistent with the first-passage time distribution of one-dimensional Brownian motion Karatzas et al. (1991); Godec and Metzler (2016). The cutoff of the residence time of protein within the droplet interior ($\sim 0.2\,\mathrm{\mu s}$) is comparable to the crossover time of the RSD, suggesting that the non-Gaussian fluctuations are due to transitions between residence in the interior and the interface.

The regions either in the interior or in the interface of the droplet can influence the structure and dynamics of proteins. The distribution of the radius of gyration $R_{g}$ has a singular peak, and proteins in the interior of the droplet exhibit a slightly larger conformation than those at the interface (see Fig. 4C). This is attributed to the hydrophobicity inside the droplet. To compare diffusive dynamics, we calculated the temporary diffusion coefficient (TDC) Yamamoto et al. (2017); Akimoto and Yamamoto (2017); Yamamoto et al. (2021),

where $d$ is the spatial dimension, $\Delta$ is $0.1\,\mathrm{ns}$, and $t$ is $10\,\mathrm{ns}$. In the TDC analysis, we consider $D(t^{*})$ for a protein interior or interface regions at $t=t^{*}$. The TDC distribution implies that proteins in the interface region diffuse faster than those in the interior of the droplet. This observation is due to the reduced density and the smaller size of proteins at the interface compared to those within the droplet (Fig. 4D). The movement between the interior of the droplet and its surface leads to significant fluctuations in diffusivity (Figs. 2B and 3).

Universality of fluctuating diffusivity in droplet–The universal nature of the diffusive protein behavior within spherical droplets is underscored by MD simulations performed at different temperatures ($320\,\mathrm{K}$), with different protein N-terminal disordered regions of DDX4 (DDX4-DR), and using different force fields (see Fig. S2–S7). Qualitatively, no clear differences are observed between FUS-LCD and DDX4-DR. When droplets form from 1,000 proteins, the density inside the DDX4-DR droplet is lower than that of the FUS-LCD droplet (Fig. S6). Additionally, because the residue length of DDX4-DR is greater than that of FUS-LCD, the droplet size of DDX4-DR becomes larger than that of FUS-LCD. This results in slightly faster diffusivity and a longer relaxation time for the fluctuating diffusivity of DDX4-DR compared to those of FUS-LCD (see Figs. 3, S5, S7). If the observation time scale is sufficiently longer than the scale at which proteins significantly diffuse between the droplet’s surface and its interior, diffusive behavior might be more determined by the size of the protein and the density of the droplet than by the sequence-specific interactions of IDPs.

Discussion– In summary, we performed MD simulations to investigate individual protein diffusion within protein droplets. We revealed subdiffusion within the droplets due to their viscoelastic nature. The diffusivity of proteins differs significantly between the interior and the interface of the droplet. The relaxation time of these fluctuations is equivalent to the residence time of the proteins within the interior of the droplets. Information on protein diffusion within the droplet and residence time in interior and interface regions will lead to an understanding of cell functions such as chemical reactions, transport of nanoparticles, and biomolecular storage using the droplet as a medium.

Performing large-scale MD simulations at the atomic level that are comparable to experimental spatiotemporal scales remains challenging Zheng et al. (2020). Residue-level CG-MD simulations with an implicit solvent are a unique and powerful tool to explore molecular phenomena over tens to hundreds of nanometers and tens of microseconds Dignon et al. (2018); Joseph et al. (2021). However, the diffusion coefficient of proteins within the droplets, as simulated using the CG models, was three orders of magnitude larger than the experimental values Burke et al. (2015); Murthy et al. (2019); Kamagata et al. (2021). This implies that simulations of $20\,\mathrm{\mu s}$ correspond to real-time durations of $20\,\mathrm{ms}$ (see Fig. S8). Fast dynamics of the residue-revel CG models have been also reported in other studies Tejedor et al. (2021, 2022), as these models were parameterized to reproduce the structural and phase behavior rather than the dynamic properties Regy et al. (2021); Joseph et al. (2021). The higher diffusivity of residue-level CG models with implicit solvent is due to the absence of factors such as side-chain entanglements and hydrodynamic interactions. One method to match the simulated diffusion coefficients with experimental values is by adjusting the friction coefficient in the Langevin equation. Increasing the friction coefficient by factors of 10 to 1,000 resulted in a decrease of the diffusion coefficient from 200 to $4.0\,\mathrm{\mu m^{2}/s}$, respectively, eventually reaching the same order of magnitude as the experimental values. Note that increasing the friction coefficient does not alter the structural property of $R_{g}$ and the density of the droplets, and the physical mechanism underlying the diffusion process remains qualitatively consistent.