Potential of mean force calculations of pore formation

Potentials of mean force (PMFs) of pore formation were computed for pure-lipid systems, as well as for the POPG systems with RSide, R4, R4cap, and R8. The PMFs were computed along a recently developed reaction coordinate, here referred to as ‘chain coordinate’ ξch (Hub and Awasthi 2017; Awasthi and Hub 2019). In brief, ξch is defined with a membrane-spanning cylinder that is decomposed into slices with a thickness of typically 1Å. Then, ξch is approximately given by the fraction of slices that are occupied by polar heavy atoms, in this work taken as oxygen atoms of water and lipid phosphate groups. As such, ξch is unit-less and probes the degree of connectivity of a polar transmembrane defect: a value of ξch≈0.25 indicates the flat membrane, with only few cylinder slices at the two headgroup regions filled by polar atoms. In contrast, a large value of ξch≳0.85 indicates a continuous, membrane-spanning polar defect. Critically, the membrane-spanning cylinder does not take a fixed lateral position in the membrane plane. Instead, the lateral position of the cylinder is dynamically defined so it can ‘follow’ the polar defect while it travels along the membrane plane. This property excludes that the system moves along the reaction coordinate by shifting the defect laterally out of the cylinder, which would lead to undesired hysteresis effects during PMF calculations (Awasthi and Hub 2016). The reaction coordinate was implemented into an in-house modification of GROMACS 2016, which is available from the authors upon request.

PMFs of pore formation were computed using umbrella sampling (Torrie and Valleau 1974), similar to previous work (Ting et al. 2018; Awasthi et al. 2019). The initial frames for umbrella windows were taken from slow-growth pulling simulations along ξch over 46 ns. The cylinder had a radius of 1.2 nm and contained 30 slices with a thickness of 0.1 nm. The parameter ζ, specifying the fraction to which a slice is filled upon the addition of the first polar atom, was set to 0.75. We used 24 umbrella windows restrained to the following, non-uniformly distributed ξch positions: 0.065–0.625 in steps of 0.08 using a force constant of 5000 kJ/mol, and 0.7– in steps of 0.02 using a force constant of 10,000 kJ/mol. For systems with arginine or arginine side chains, each window was simulated for 1 μs, where the first 0.5 μs were omitted for equilibration. For pure-lipid systems, each window was simulated for 150 ns, where the first 20 ns were omitted for equilibration. Here, the temperature was controlled using a stochastic dynamics integration scheme (τ=0.5 ps) (Gunsteren and Berendsen 1988). All other parameters were chosen as described above. The PMFs were computed with the weighted histogram analysis method (WHAM) (Kumar et al. 1992), as implemented in the gmx wham module of GROMACS (Hub et al. 2010). Statistical errors were estimated using 50 rounds of Bayesian bootstrapping of histograms (Hub et al. 2010). Before computing the average and standard deviation of the PMFs obtained from bootstrapped histograms, the PMFs were defined to zero at the PMF minimum (ξch≈0.25).