Friday, December 22, 2017

NMRlipids IV: Current status and reorganization of the manuscript

First of all, I want to thank again all the contributors for delivering significant amount of useful data about PE, PG and PS headgroups for the NMRlipids IV project.

In addition to simulation results, also experimental signs of headgroup and glycerol backbone order parameters were contributed for the POPS lipid. However, experimental signs for PE and PG are not yet known, which makes the comparison between simulations and experiments more ambiguous. Therefore, I have divided the manuscript into two parts. The first one contains results only from systems with PS lipids and the other contains results from systems with PE and PG lipids. This should also ease the management and completion of the projects. The manuscript about PS lipid systems should be more straightforward to finish and I have started to compile it towards a submittable version. Some of the current results and the most important open tasks are listed below.

Headgroup & glycerol backbone structures of PS lipid bilayers

Since the order parameter signs are known for the PS lipid headgroup and glycerol backbone, we can perform similar comparison between simulations and experiments as was done in NMRlipids I for the PC lipids. This is shown Fig. 1; the subjective quality assessment is also available. The tested models perform generally less well than the PC lipid models discussed in the NMRlipids I publication. Many relevant contributions have already been made to give structural insight to the differences, however, I think this deserves more attention (see the ToDo list below).
Figure 1: Headgroup and glycerol backbone order parameters of the PS lipids from different simulation models and experiments. 

Headgroups in mixtures of PS and PC lipids

Lipid mixtures are biologically relevant, and in fact experimental data for ion binding to PS lipid bilayers seems to be available only for PS/PC mixtures. Therefore, we also address the mutual interaction between PS and PC lipid headgroups. Figure 2 shows how the addition of POPS changes the order parameters of POPC (left column), as well as the changes of POPS order parameters with increasing amount of POPC (right column). It seems that the tested simulation models do not reproduce the interactions between PC and PS headgroups very well. However, the inaccuracies in counterion binding may also disturb the lipid headgroup structures when the amount on PS (and counterions) is increasing. I think that we need data with a few more different models before drawing general conclusions, and possibly also some data with different counterion concentrations (see the ToDo list below).
Figure 2: Headgroup order parameters from PC:PS mixtures from different simulation models and experiments. Left panel shows the PC headgroup order parameters and right panel shows the PS headgroup order parameters.

Interactions between cations and lipid bilayers containing PS

As already discussed in the opening post of NMRlipids IV project, molecular electrometer experiments show that the presence of negatively charged lipids enhance cation binding in lipid bilayers. The available experimental data for the lipid headgroup order parameters of POPC:POPS (5:1) mixture as a function of CaCl2 concentration are shown in Fig. 3 together with the simulations ran with the MacRog model. In line with the NMRlipids II publication, the PC headgroup order parameters' decrease with increasing CaCl2 concentration is overestimated in the MacRog simulations, indicating overestimated binding of Ca2+ to the lipid bilayer. It should, however, be noted that the point with the lowest CaCl2 concentration is in better agreement with the experiments. Potential explanation could be an overestimated screening effect by the overbound counterions, but this requires further analysis. Order parameters of the POPS headgroup rapidly increase or decrease even with low CaCl2 concentration and reach almost a plateau value above 100 mM. Also the order parameter changes of the POPS headgroup are overestimated in the MacRog simulations and a qualitative agreement for the alpha carbon is unclear. Simulation data of PC:PS mixtures with different CaCl2 concentrations is required also from other than MacRog model for more general conclusions (see the ToDo list below).
Figure 3: Headgroup order parameters from PC:PS (5:1) mixtures with different CaCl2 concentrations from the MacRog simulation model and experiments. Left column shows the PC headgroup order parameters and right column the POPS headgroup order parameters.

ToDo list

  1. Structural relevance of the observed order parameter differences between different simulation models and experiments are now analyzed using dihedral distributions (see Fig. 13 in the manuscript) and we have pictures of the sampled conformations of glycerol backbone and phosphate. I think that we should apply the tools contributed by Pavel Buslaev to calculate dihedral distributions and to visualize the structural sampling also for other contributed models. Among these data we should then select the parts giving the best representation for structural differences related to order parameters. To do this in practice, we need more simulation trajectories available in Zenodo.
  2. We need data for PC:PS mixtures without additional ions from few more simulation models. Some data from Berger model has already been delivered by Lukasz Cwiklik (see also data with calcium), but we still need simulations with counterions only. Simulations of PC:PS mixtures with the Slipids and CKP models may also be useful.
  3. Simulations of POPC:POPS mixtures with different CaCl2 concentrations are needed also from other force fields than MacRog. The above mentioned dataset with Berger force field complemented with simulations containing only counterions would be highly useful. I think that simulations with CHARMM36 is a must. Also simulations with the Slipids and CKP force fields may be worth of doing. Based on NMRlipids II, it is expected that cations will overbind to lipid bilayers in these simulations. However, the behavior of PS headgroup as a function of salt concentration in different force fields with respect to experiments is not yet known.

27 comments:

  1. Dear Samuli,

    I have a few trajectories (and will generate the rest) for PC:PS (5:1 ratio) mixtures with the charmm36 force field, both neutral and with CaCl_2.

    Best,
    Jesper

    ReplyDelete
    Replies
    1. Excellent! Looking forward the results.

      Delete
    2. The results are now available here: https://www.dropbox.com/s/97w8l2fofo4e1qw/PCPS-results.txt?dl=0

      Best,
      Jesper

      Delete
    3. Thanks a lot!

      I added the results in current Figs. 8 and 9:
      https://github.com/NMRLipids/NMRlipidsIVotherHGs/blob/master/Figs/HGorderparametersPCvsPS-eps-converted-to.pdf

      https://github.com/NMRLipids/NMRlipidsIVotherHGs/blob/master/Figs/CHANGESwithCaClPS-eps-converted-to.pdf

      Unexpectedly, the results seem to indicate that the CaCl_2 binding to PS containing bilayer would be underestimated in CHARMM36. At the same time the headgroup order parameter changes in mixtures with different PC:PS ratios without additional CaCl_2 are not very well reproduced.

      I am suspecting now that these results may be related to the counterion binding. I am also realizing now that, obviously, in the presense of PS we always have counterions present and we have to take the ion binding into account also in the results without added CaCl_2.

      Could you deliver also the density profiles of lipids water and ions (also counterions) calculated from this data (and/or the trajectories for the further analyses)?
      Also, there seems to be the order parameter values for other beta carbon of PC missing from the results.

      Delete
    4. No problem!

      -The order parameter(s) for the other beta carbon of PC are now added to the results file (same URL as before).
      -The two "neutral" (w/zero added salt) trajectories are neutralized using potassium counterions, not sodium.
      -Density profiles of components in the trajectories are available here: https://www.dropbox.com/s/nnwee8z9q2rg9x7/dens_profiles_PCPS.tar.gz?dl=0

      Furthermore, I will look into getting the trajectories uploaded.
      Best,
      Jesper

      Delete
    5. Are the simulations with CaCl_2 ran then without any other salt, i.e. potassium or sodium?

      Delete
    6. I plotted the density profile of Ca2+ from these mixed PC/PS simulations together with previously reported PC/PG simulations and data from NMRlipids II for pure PC bilayer:
      https://github.com/NMRLipids/NMRlipidsIVotherHGs/blob/master/scratch/CAdensPSvsPCcharmm36100mM.pdf

      There is least binding in this PC/PS mixture and most in pure PC. Obviously, this should not be the case as PS and PG are negatively charged.

      The most obvious explanation would be the equilibration, which is known from NMRlipids II to be very slow:
      https://github.com/NMRLipids/lipid_ionINTERACTION/blob/master/Fig/bindingINlongRUNS-eps-converted-to.pdf

      However, also the pure PC data used is only 200ns and therefore not fully equilibrated. Maybe the equilibrium is even slower in the PS system? How long were these simulations actually?

      Delete
    7. Yes, in the simulations with CaCl_2 the concentration number refers to that of CaCl_2 and no other salt is added.

      Delete
    8. I think that it would be good to extend these simulations at least another 200ns or even longer if possible. Do you have possibility to do that or do you want to share the files such that someone can continue?

      Delete
    9. Yes, I can do this. Will report back once the results are in.

      Delete
    10. Jesper, can you check from your charmm36.itp file if it contains the lines

      [ nonbond_params ]
      ; i j func sigma epsilon
      CAL CLA 1 3.29454345968e-01 5.61342176000e-01
      CAL O2L 1 2.94352936474e-01 5.02080000000e-01

      These lines were not present in simulations ran during NMRlipidsII for PC lipids (e.g., http://doi.org/10.5281/zenodo.35159).

      It seems to that they have added some special interaction between oxygens and calcium at some point. I think now that this might be the reason for lower binding in PC:PS mixtures compared to pure PC, instead of equilibration issues.

      Delete
    11. Dear Samuli,
      Good catch! I can see that my most recent simulations with the POPC/POPS/CaCl2 mixtures have the interaction definitions that you mention:
      "[ nonbond_params ]
      ; i j func sigma epsilon
      CAL CLA 1 3.29454345968e-01 5.61342176000e-01
      CAL O2L 1 2.94352936474e-01 5.02080000000e-01
      CAL OCL 1 2.87314836600e-01 5.02080000000e-01"

      whereas my previous runs of POPC/POPG mixtures do not! Perhaps this is due to the new standard choice of force field in the charmm-gui (CHARMM36 vs. CHARMM36m)? In the brief paper by Huang et al. (Nature Methods volume 14, pages 71–73 (2017)), they describe changes in certain salt bridge interactions in their force field going from C36 to C36m. Apparently, some parameters were adjusted to better agree with experimental equilibrium association constants and osmotic pressure measurements.

      It may therefore be interesting to simulate the same systems for the original version of the C36 FF in order to get a direct comparison with previous results.

      Best,
      Jesper

      Delete
    12. I contacted CHARMM-GUI people and they told me that the additional interactions are the same as in this publication: https://doi.org/10.1016/j.bpj.2016.09.001
      The origin of the numbers is unpublished work by B. Roux and S. Rong.

      It seems to me that the calcium binding to lipid bilayers would be too weak with these parameters.

      I agree that we should run the simulations with the original version of the C36FF, i.e. without these additional parameters. Could you rerun your simulations after removing these lines?

      Delete
    13. Thank you for the clarification! Yes, I'll re-run the simulations using the original C36 FF parameters and report back these results.

      Delete
    14. Dear Jesper,

      after a discussion with Samuli, I run a simulation of POPC+CaCl with C36FF with this new set of parameters and compare the binding of calcium with the simulation without these parameters and with the experimental values.

      It seems that binding of calcium is too weak.

      Data are here

      Order parameter:
      https://github.com/NMRLipids/NMRlipidsIVotherHGs/blob/master/Figs/OP_CHARMM_CaCl_POPC_NBFix.pdf

      Density profiles:
      https://github.com/NMRLipids/NMRlipidsIVotherHGs/blob/master/Figs/density_profile_CHARMM_CaCl_POPC_NBFix.pdf

      Simulations data:
      https://zenodo.org/record/1198171#.WqqBF-ZG290

      Delete
    15. Hi Ricky, Thanks that's great to know! That NBfix appears to make a relatively large difference on the OPs indeed. It is very striking that essentially all calcium binds around the phosphate groups in C36, whereas most remains in the water in NBfix, so I agree with you that the binding looks too weak from that point of view.

      Best,
      Jesper

      Delete
    16. Thanks for tha data!

      I now merged also the second pull request by Ricky and now there are also the results for sodium in the plots behind above links. Notably, CHARMM simulations with NBfix in CaCl_2 gives similar binding and order parameter changes for sodium and calcium. Also, the order parameter changes with calcium are too small in simulations of POPC and POPC/POPS mixtures with CHARMM/NBfix parameters from CHARMM-GUI. Therefore, I think that we can conclude that these parameters gives too weak binding for Ca2+ in both POPC and POPC/POPS bilayers.

      Delete
  2. Hi,

    If it's going to be any use, in the following month I will generate different PC:PS mixtures, as well as POPC:POPS mixtures with different CaCl2 concentrations with Amber Lipid 17 force field.

    ReplyDelete
    Replies
    1. As discussed in another post (http://nmrlipids.blogspot.com/2017/03/nmrlipids-iv-headgroup-glycerol.html?showComment=1515574773714#c1375415424280576642),
      the Lipid17 do not give very good structure for PS lipids.

      However, the mixtures would be somewhat interesting to see if the response of PS to PC or ions is qualitatively correct despite of the incorrect initial structure. This was the case for PC response to dehydration and ions in NMRlipids I and II.

      Delete
    2. Is there by any chance a POPS Lipid17 topology in a Gromacs format (i.e. *.itp) somewhere?

      Delete
    3. Sorry, we ran our DOPS and POPS directly in Amber. However, I suppose one could create the Gromacs files from the Amber files using ACPYPE?

      Delete
    4. Hi Batuhan,

      when reading my previous response, I realise that the wording "somewhat interesting" is clearly a understatement. This data would be actually very interesting. Especially because I cannot currently see almost any logic on the responses of PS headgroups to ions. In addition to CaCl_2, I would be also interested to see the results for the excess amount of counterions (Na+ and/or K+) for POPC:POPS (5:1) mixture, see the discussion below.

      I have also a question about Lipid17: Do you know if there exist a Lipid17 model also for POPC? In other words, would POPC be described as Lipid14 or Lipid17 in POPC:POPS mixture?

      Delete
    5. Hi Samuli,

      To the best of my knowledge, POPC should be described as Lipid14. Lipid17 paper is still not out. I've checked the available forcefield files, and it seems like they did not make any changes to Lipid14 parameters for POPC. So I think POPC will be described as Lipid14, whereas POPS will be as Lipid17.

      I will start the simulations of POPC:POPS mixture with Na+ and K+ this week. I'll keep everyone posted about the outcomes.

      Delete
  3. When collecting the data and composing the manuscript, I have realized that we must focus also on the monovalent ion binding to PS containing bilayers. There are always counterions present in these systems and their binding presumably changes the PS conformation and order parameters. Therefore we cannot have a fully neutral refence system when PS lipids are present.

    I noticed that the sodium binding to PS lipids was studied by Roux and Neumann (https://doi.org/10.1016/0014-5793(86)81218-2) by measuring the changes of the alpha carbon order parameters in PS headgroup in DMPC:DMPS (3:1) mixture with excess amounts of NaCl. I prepared the corresponding simulations using CHARMM-GUI and plotted the results with experimental values: https://github.com/NMRLipids/NMRlipidsIVotherHGs/blob/master/Figs/PSresponseTONaCl-eps-converted-to.pdf

    It seems that the changes in simulations are too large. Interestingly, the changes in simulations seem similar as observed with CaCl_2 in experiments: https://github.com/NMRLipids/NMRlipidsIVotherHGs/blob/master/Figs/CHANGESwithCaClPS-eps-converted-to.pdf

    However, the results are not fully conclusive because the smallest buffer concentration in the experimental data was 100 mM. It is possible that the changes in experiments are already saturated before this, which is the case with CaCl_2. Therefore this dataset is actually not optimal for these studies.

    I realized that the response of headgroup order parameters to monovalent cations is also measured for the same POPC:POPS (5:1) mixture for which the response to CaCl_2 was measured (Figs. 7 and 8 in https://doi.org/10.1021/bi00482a019). Also the responses of beta carbon of POPS and order parameters of PC to added KCl concentration are measured in this publication for the POPC:POPS (5:1) mixture. In addition, the PS data for excess sodium is reported by starting from the system with only counterions.

    Therefore I think that we should use this data (https://doi.org/10.1021/bi00482a019) to compare the monovalent ion binding affinity between simulations and experiments. We already have the simulations of POPC:POPS (5:1) mixtures with K+ and Na+ counterions from CHARMM36 model (see the updated table I for the list of simulated systems). Therefore, CHARMM36 simulations of this mixture with excess concentration of roughly 500 mM and 1000 mM of NaCl and KCl using CHARMM36 would be highly useful. For MacRog we already have the system with K+ conterions, thus the systems with the larger concentrations and maybe also with sodium counterions would be useful. Corresponding simulations also with other models would be useful.

    ReplyDelete
    Replies
    1. I have now calculated the headgroup responses as a function of additional NaCl concentration from Berger simulation of POPC:POPS (4:1) mixture (data reported previously http://nmrlipids.blogspot.com/2017/03/nmrlipids-iv-headgroup-glycerol.html?showComment=1501687775751#c6314031175464531229).

      The results are shown together with the experimental numbers from https://doi.org/10.1021/bi00482a019 (extracted and committed to GitHub repository by J. Melcr) here:
      https://github.com/NMRLipids/NMRlipidsIVotherHGs/blob/master/Figs/CHANGESwithMONVALENTwithPS-eps-converted-to.pdf

      There is no experimental data available for the PC headgroup response to NaCl in mixtures with PS. However, comparison between the PC headgroup response to NaCl in simulations and the response to KCl and LiCl in experiments suggests that the binding of Na+ would be similar to Li+. Large amount of experimental data, however, suggests that the binding affinity of Na+ should be lower than for Li+. Therefore, the binding affinity of additional Na+ seems to be too high, as expected based on NMRlipids II.

      The response of PS beta carbon is also close to Li+ response in experiments, but alpha carbon is completely different than in any of the experiments (the values shoot out from the y-axis range).

      The results from other models would be highly useful for this plot.

      Delete
    2. I have now analyzed the MacRog simulations with the additional KCl concentrations, delivered by Matti Javanainen
      https://doi.org/10.5281/zenodo.1210255
      http://nmrlipids.blogspot.com/2017/03/nmrlipids-iv-headgroup-glycerol.html?showComment=1522563298136#c4868166750009073094

      The results are added into the figure:
      https://github.com/NMRLipids/NMRlipidsIVotherHGs/blob/master/Figs/CHANGESwithMONVALENTwithPS-eps-converted-to.pdf

      Rough conclusions could be:
      - Sodium in Berger binds like Lithium
      - Binding affinity of potassium in MacRog agrees with experiments within (quite large) error bars.
      - Response of PS to potassium is overestimated in MacRog, despite the binding affinity based on PC results is not.
      - Response of PS alpha carbon to monovalent cations in both models is strange.

      Delete

Please sign in before writing your comment.