SpiDec: Computing binodals and interfacial tension of biomolecular condensates from simulations of spinodal decomposition

Abstract

Phase separation of intrinsically disordered proteins (IDPs) is a phenomenon associated with many essential cellular processes, but a robust method to compute the binodal from molecular dynamics simulations of IDPs modeled at the all-atom level in explicit solvent is still elusive, due to the difficulty in preparing a suitable initial dense configuration and in achieving phase equilibration. Here we present SpiDec as such a method, based on spontaneous phase separation via spinodal decomposition that produces a dense slab when the system is initiated at a homogeneous, low density. After illustrating the method on four model systems, we apply SpiDec to a tetrapeptide modeled at the all-atom level and solvated in TIP3P water. The concentrations in the dense and dilute phases agree qualitatively with experimental results and point to binodals as a sensitive property for force-field parameterization. SpiDec may prove useful for the accurate determination of the phase equilibrium of IDPs.

Keywords: binodal; biomolecular condensates; interfacial tension; phase equilibrium; phase separation; spinodal decomposition.

FIGURE 1

Binodals and spinodals of two model systems. (A) Van der Waals fluid, which satisfies the following equation of state: $\left(P+\frac{9}{8}{k}_{\mathrm{B}}{T}_{\mathrm{c}}\frac{{\rho }^2}{{\rho }_{\mathrm{c}}}\right)\left(1-\frac{\rho }{3{\rho }_{\mathrm{c}}}\right)=\rho k_{\mathrm{B}}T$ , where $P$ , $\rho$ , $T$ denote the pressure, density, and temperature, $T_{\mathrm{c}}$ and ${\rho }_{\mathrm{c}}$ denote the critical temperature and critical density, and $k_{\mathrm{B}}$ is the Boltzmann constant. (B) Symmetric polymer blend that follows the Flory-Huggins theory for the Helmholtz free energy: $\frac{F}{M{k}_{\mathrm{B}}T}=\frac{\varphi }{L}\ln \varphi +\frac{1-\varphi }{L}\ln \left(1-\varphi \right)+\chi \varphi \left(1-\varphi \right)$ , where $M$ is the total number of polymer chains, $L$ is the number of beads per chain, $\varphi$ is the mole fraction of one polymer species in the binary blend, and $\chi$ measures the energy gap between inter- and intra-species interactions. Inside the spinodal region, the system is unstable and phase separates by spinodal decomposition; between the binodal and spinodal, the system is metastable and phase separates by nucleation and growth. Note that there has been some controversy regarding the location of the spinodal in computer simulations (Binder et al., 2012; Díaz-Herrera et al., 2022).

FIGURE 2

Morphologies of the dense phases of LJ particles and LJ chains over a range of initial densities inside the spinodal. (A) Morphologies for an LJ particle system of N = 1,000 particles in a cubic box at T = 0.65. The dense phase appears as a sphere, cylinder, slab, hollow cylinder, and hollow sphere at ρ ₀ = 0.1, 0.2, 0.3, 0.6, and 0.7 respectively. The initial densities were changed by varying the box side lengths, which are shown here not to scale. At each density, the cubic box is cut by a plane (rendered as gray when the background is empty), and only the half behind the cut is displayed. (B) Corresponding results for an LJ chain system (100 10-bead chains) at T = 1.7 and ρ ₀ = 0.1, 0.2, 0.3, 0.5, and 0.6. (C) Boundaries between different morphologies for LJ particles in rectangular boxes with different L _z /L _x ratios. The density ranges for different morphologies are illustrated in the inset. L _x = 10 and T = 0.65. (D) Corresponding results for an LJ chain system at L _x = 13 and T = 1.7.

FIGURE 3

Slab formation at various L _z /L _x values. (A) LJ particle system at N = 4,000, T = 0.65, ρ ₀ = 0.3, and L _z /L _x = 1.25, 3.26, and 13.3. (B) LJ chain system at N = 4,000, T = 1.7, ρ ₀ = 0.25, and L _z /L _x = 1.16, 3.79, and 7.28. The z axis is along the vertical direction. For each system, the simulation boxes are drawn approximately to scale. Arrows indicate the fusion of multiple slabs into a single slab.

FIGURE 4

Slab formation at various L _z /L _x values. (A) HP chain system at N = 4,000, T = 1.05, ρ ₀ = 0.25, and L _z /L _x = 1.5, 3.0, and 5.0. (B) Patchy particle system at N = 500, T = 0.61, ρ ₀ = 0.36, and L _z /L _x = 1.5, 3.0, and 5.0. The z axis is along the vertical direction. For each system, the simulation boxes are drawn approximately to scale. Arrows indicate the fusion of multiple slabs into a single slab.

FIGURE 5

Binodals and interfacial tensions calculated from snapshots with a single slab. (A) Binodals of the LJ particle system at N = 4,000 and L _z /L _x = 1.25, 3.26, and 13.3. (B) Binodals of the LJ chain system at N = 4,000 and L _z /L _x = 1.16, 3.79, and 7.28. (C) Interfacial tension versus temperature for the LJ particle system at the three L _z/L _x ratios. (D) Interfacial tension versus temperature for LJ chain system at the three L _z/L _x ratios.

FIGURE 6

Application of SpiDec to a phase-separating peptide. (A) Structure of FFssFF. (B) The SpiDec simulation procedure. First, a random, loose configuration condensed into a single slab. Then the slab was solvated into an elongated box and the two phases reached equilibrium. The fourth snapshot shown was taken at 1.065 µs of a 2-µs simulation at 294 K. (C) Concentration profile (circles) at 294 K and fit to eq [6] (red curve). Concentrations were calculated by counting copy numbers in 2-Å slices along the z direction; the location of each copy was represented by the midpoint of the central disulfide bond. The conversion of concentrations to wt/wt used a molecular weight of 741.5 Da for the peptide and a density of 1 kg/L for water. (D) Dilute- and dense-phase concentrations at four temperatures, from 294 K to 360 K. The red curve shows a binodal to guide the eye.