NMRlipids VI: First results from polarizable Charmm-Drude force field

After analyzing the simulation data for the POPC lipids at various NaCl and CaCl₂ concentrations using the Drude polarizable model, I have compiled all the results in the manuscript. You can access all the figures, scripts, and the source files used in the analysis in the GitHub repository. Here I’d like to summarize the main results we have obtained so far.

1) The headgroup and glycerol backbone order parameters predicted by Drude force field are very different than the ones with Charmm36 force field.

The Fig. 1 in the manuscript show that Drude model predicts a forking at the beta and alpha order parameters, which are not predicted by the Charmm36 force field. 

Figure 1: Headgroup and glycerol backbone order parameters of POPC without ions from Charmm-Drude and Charmm36 simulations

We further see that the response of these order parameters to changing Na⁺/Ca²⁺ concentrations also deviate significantly from the Charmm36 results, as shown in Fig. 2.

Figure 2: Changes in the headgroup order parameters upon addition of ions in Charmm-Drude and Charmm36 simulations.

 

2) The dihedral distributions around the lipid headgroup and glycerol backbone have completely different profiles compared to the Charmm36 force field.

In Figs. 3 and 4, we show the dihedral angle distributions around the lipid headgroup and glycerol backbone at different ion concentrations. We observe the largest differences, as a function of the ion concentration, at the O_g3-P-O_α-C_α and the g₃-O_g3-P-O_α pairs. For the remaining of the distributions, we don’t observe large changes.

When compared to the Charmm36 results [Fig.2 of NMRlipids IVb manuscript], we see that the Drude model yields completely different dihedral angle distributions.

3) The density profiles for the ions suggest that the sodium binding is too strong.

Following Samuli’s comment: "Based on rough comparison of calcium ion density profiles to the results in Fig.5 in http://dx.doi.org/10.1021/acs.jpcp.7b12510, it seems that the binding affinity of calcium could be similar to the model which gives the correct order parameter responses. However, based on both order parameter changes and density profiles, the sodium binding to membranes seems to be clearly stronger than calcium, which is not line with the experiments."

Overall, our calculations show that the Charmm Drude polarizable model does not improve on the results of the Charmm36 force field. In fact, we see that the Drude polarizable model yields worse results than Charmm36, as evident from the forking of the beta and alpha order parameters in Drude model. This suggests that inclusion of the polarizability does not automatically improve the headgroup conformational ensembles and ion binding affinities to the membranes. Based on the current evidence, we are not planning to run any more simulations using the Drude model (for the POPC lipids) and switch our focus to the other polarizable force fields.

From now on, we are planning to concentrate on obtaining simulation results for the AMOEBA and other polarizable force fields. Still, we think it could be very helpful to have data from Charmm-Drude simulations (also older versions) for cross-checking our results. To be consistent with the Drude simulations, our first focus will be the pure POPC membranes at different CaCl₂ and NaCl concentrations. It seems like we will be obtaining some trajectories soon. However we kindly invite anyone who might have some trajectories with the AMOEBA and other force fields to share them with us.

Any comment on the current results or data contributions, as always, will be greatly appreciated.

 

Batuhan Kav

NMRlipids IVb: Toward submission of the manuscript with PE and PG results

After the previous post on the NMRlipids IVb manuscript, I have made further analysis on the structural differences between PC, PE, PG and PS headgroups using the development version of the NMRlipids databank and CHARMM36 simulations, which give the headgroup and glycerol backbone order parameters closest to experiments for all lipids. Based on these analyses, I have reformulated the NMRlipids IVb manuscript to focus on the headgroup conformations of different lipids.

Main points of the manuscript are now:

1) Differences in NMR order parameters between different headgroups can be explained by the changes in dihedral angle distributions, suggesting that similar conformations are accessed by all headgroups, but with different probabilities.

CHARMM36 simulations approximately reproduce the main differences of headgroup order parameters between PC, PE, PG and PS lipids (forking of the ⍺-carbon in PS and positive β-carbon value in PG in Figs. 1 A and B). These differences can be explained by the differences in the heavy atom dihedral distributions (Fig. 1 C).

Figure 1: Headgroup and glycerol backbone order parameters of POPC, POPE, POPG and POPS from (A) experiments and (B) CHARMM36 simulations. (C) Heavy atom dihedral angle distributions from CHARMM36 simulations.  

2) The changes in headgroup order parameters upon addition of charged molecules can be explained by the changes in dihedral angle distributions close to the phosphate region.

CHARMM36 simulations reproduce the experimental changes in headgroup order parameters upon addition of cationic surfactants which are related to the substantial tilt of P-N dipole (Fig. 2 A) and changes in dihedral distributions close to the phosphate region (Fig. 2 B). 

Figure 2: Response of POPC headgroup A) order parameters, P-N vector and B) dihedral angles to the addition of cationic surfactant from CHARMM36 simulations compared with experimental data.

3) Wide range of headgroup conformations are observed also in lipids that are directly bound to proteins, suggesting that the specific binding of lipids to proteins is dominated by the intermolecular lipid-protein interactions rather than the differences in conformational restrictions of lipids.

Dihedral angles calculated from protein bound lipid structures from the protein data bank (PDB) have wide distributions independently on the headgroup chemistry (Fig. 3).

Figure 3: Dihedral distributions of protein bound lipids analysed from the protein data bank (PDB).

Comparison of headgroup and glycerol backbone order parameters between different force fields and experiments, and results from lipid mixtures are now in the supplementary information because conclusions from these are essentially the same as in NMRlipids I and IV publications. None of the force fields correctly capture the conformational ensembles of lipid headgroups, but CHARMM36 is closest, and simulations containing charged lipids are complicated by the overestimated binding affinity of counterions. 

I believe that the manuscript is approaching the submission stage. Because MD simulations from the NMRlipids project are able to explain the differences in experimental order parameters between different lipids, thereby giving an interpretation of the lipid conformational ensembles, my current plan is to submit the manuscript to the ACS Central Science.

Most critical things to do before submission are related to the methods, see todo points in the manuscript and especially in the supplementary information. There are also some issues still open in the GitHub. All kind of comments on the manuscript are welcomed.