tag:blogger.com,1999:blog-1283194685237019772.post5149020296844720567..comments2024-03-26T17:11:28.045+02:00Comments on The NMRlipids Project: NMRlipids IVb: Assembling the PE & PG resultsmarkus.http://www.blogger.com/profile/05973383391755193687noreply@blogger.comBlogger37125tag:blogger.com,1999:blog-1283194685237019772.post-42354463238438877802020-03-27T12:46:43.309+02:002020-03-27T12:46:43.309+02:00Hi Josef,
sorry I didn't notice your comment! ...Hi Josef,<br />sorry I didn't notice your comment! I don't know about this kernel density estimation. I'll have a look into it.Patrick Fuchshttps://www.blogger.com/profile/16910815195856162826noreply@blogger.comtag:blogger.com,1999:blog-1283194685237019772.post-52081924719347625312019-11-18T18:26:45.792+02:002019-11-18T18:26:45.792+02:00Hi -
do you think that using kernel density estim...Hi -<br /> do you think that using kernel density estimation (KDE) could be potentially useful to capture the uncertainty in the hydrogen direction?<br /><br />If yes, I would use the uncertainty of the hydrogen reconstruction as the width of the kernel function (usually gaussian) used to estimate the distribution of orientations (i.e. also mean and standard deviation). Josef Melcrhttps://www.blogger.com/profile/12422047696918877379noreply@blogger.comtag:blogger.com,1999:blog-1283194685237019772.post-35313877955414804992019-11-04T18:29:24.041+02:002019-11-04T18:29:24.041+02:00Hi all, here are my (late) comments on the two que...Hi all, here are my (late) comments on the two questions.<br /><br />1) About the difference between reconstructed Hs and all-atom (AA) Hs. To me there are two separate things:<br />- The fact that we cannot recalculate the exact same position of some Hs (oleoyl double bond, glycerol Hs, in general Hs close to a polar atom) because there are other subtleties in the FF (bonded parameters according to Thomas, see his comment of october 14) that we cannot catch with simple geometrical rules. Let's call it geometric error.<br />- The fact that for an AA force field, the "lipid backbone" (I mean the succession of heavy atoms on which the Hs are attached) will be incorrectly oriented on average because of the FF. Let's call it FF error.<br />Sure we have to accept a difference and indeed some AA models will be away from experimental values. But I'd like to discuss with you up to what value this difference is acceptable. From what I could see from buildH / CHARMM36 comparison on a full traj, the max absolute difference is slightly less than 0.02 (on one H of glycerol), and on average I found an error of 0.0034 on all Hs that are not CH3. To me, it is acceptable since it is on the order of the experimental uncertainty on the worst case, and almost always it is below. Does everyone agree?<br />If we place ourself in the context of UA simulations, we want to evaluate the FF error but we will always have a systematic geometric error. So that it is important to evaluate the latter accurately (see below).<br /><br />2) I agree about the statistical error. In our simulations of hundred(s) of ns, each lipid samples quite well its conformational space. So with the length (> 100 ns) & nb of lipids (>= 100) we use, the error is very small (most of the time <= 0.001). I think we can forget about it for simple lipid compositions like those we studied so far.<br />But, as suggested by Angel, it would be important to provide an error which reports on the accuracy of H reconstruction (the one I called reconstruction error). The only way is to compare H reconstruction to existing AA trajectories. Angel proposes to try using different reconstruction methods (different codes). I think it's indeed a good idea to make a thorough analysis. One easy way could be:<br />- to agree on a set of AA trajectories (also I never tested two different trajectories for a given FF, maybe that'd be worth a try?) ;<br />- then use buildH, Angel's script, g_protonate+awk script (some others could go in) on these trajectories ;<br />- report on a shared document (e.g. framacalc) the OP of the AA traj and for the different methods of H reconstruction (one for each column), one sheet per traj.<br />We'll likely find differences between Hs (e.g. higher error for glycerol Hs than for aliphatic CH2). So it may be doable to come up with different errors according to each H.<br />What do you think?Patrick Fuchshttps://www.blogger.com/profile/16910815195856162826noreply@blogger.comtag:blogger.com,1999:blog-1283194685237019772.post-74042136692817008032019-10-14T16:23:19.532+03:002019-10-14T16:23:19.532+03:00As to different protonation methods, I didn't ...As to different protonation methods, I didn't really see differences when looking at the tails (well, if they were doing what they claimed to and not something incorrect). <br /><br />The main differences arose in the double-bond due to the use of two fairly different approaches. One approach (e.g. as employed by protonate and g_lomepro) measures the explicit C-C=C angles and places the hydrogen atom bisecting this. The other method (e.g. as used by the fixed version of g_order) assumes the C-C=C angle is 120 degrees and builds this into the order parameter calculation. It turns out that neither of these place the hydrogen atoms where they are in the all-atom simulations.Thomas Piggothttps://www.blogger.com/profile/14030778642272576601noreply@blogger.comtag:blogger.com,1999:blog-1283194685237019772.post-6522379297756384402019-10-14T16:14:48.015+03:002019-10-14T16:14:48.015+03:00Sorry for the slow reply, I've just seen this ...Sorry for the slow reply, I've just seen this now.<br /><br />In the paper you mention, I go into all this in some detail. The main gist of it is that if you are trying to reproduce all-atom results using a united-atom methodology it can be very difficult. This is because in the all-atom simulations the explicit hydrogen atoms may not just be in some idealised geometry (as the united-atom analysis methods assume in two different main ways) but deviate from these due things like the bonded force field parameters. For CH2 groups in lipid tails, this deviation from an ideal tetrahedral geometry doesn't tend to be very large and so the united-atom analysis methods are in general very good. However, for the CH groups in the unsaturated double-bonds, this is a different case. Here the positions of the explicit hydrogen atoms in the all-atom force fields isn't where either of the two different united-atom methods predict.<br /><br />Cheers<br /><br />TomThomas Piggothttps://www.blogger.com/profile/14030778642272576601noreply@blogger.comtag:blogger.com,1999:blog-1283194685237019772.post-84134628895613145432019-10-07T03:13:09.384+03:002019-10-07T03:13:09.384+03:00Thank you very much for your answer and explanatio...Thank you very much for your answer and explanation. I do agree with your last comment, so for the representation of OP values obtained from reasonable trajectories (bilayers with several hundred lipids over several tens of ns) I guess we can forget the statistical error, since it is negligible. Other way of obtaining uncertainties could be applying different methods to rebuild the Hydrogens... as I stated in my previous post, different reference atoms can be employed to get the position of the same hydrogen. I guess this explains some small differences between calculations performed with different codes.Ángel Piñeirohttps://www.blogger.com/profile/09529986464564873340noreply@blogger.comtag:blogger.com,1999:blog-1283194685237019772.post-9176228208344324902019-09-30T21:15:11.131+03:002019-09-30T21:15:11.131+03:00I agree with the first point.
When calculating th...I agree with the first point.<br /><br />When calculating the standard deviation of the mean, the samples should be independent. However, successive frames in trajectory are probably not independent and block averaging methods should be used. In NMRlipids, we have instead assumed that each lipid samples its conformations independently from other lipids, which is the case in fluid phase by definition. Therefore, we first calculated the order parameters for individual lipids averaging over time, and then calculated the standard deviation of the mean over independent lipids. The STD is smaller in this case, because individual lipid values are already averaged over time. However, the error is similar as in the other error calculation method because the number of independent samples is smaller in the latter case. Different kinds of arguments can be made to justify different kind of statistical error calculation methods, but all approaches seem to give similar and very small errors from simulations.<br /><br />As demonstrated in the discussion here, uncertainty from other sources can be larger than the statistical error. Especially for the results from united atom simulations, I think that the statistical error underestimates the accuracy of the calculation.Samuli Ollilahttps://www.blogger.com/profile/06106569992787533569noreply@blogger.comtag:blogger.com,1999:blog-1283194685237019772.post-6952402161575466982019-09-26T22:43:56.068+03:002019-09-26T22:43:56.068+03:00Dear Samuli, Patrick, Josef and any other interest...Dear Samuli, Patrick, Josef and any other interested or involved in the calculation of order parameters. I would like to know your opinion, and to reach an agreement, on a couple of points:<br /><br />It seems that there are some differences in the order parameter values depending on the calculation method and on the force field. Rebeca found that the force fields for which the OP determined using the explicit and the rebuilt H atoms (based on optimized geometry) are significantly different, match with those for which the theoretical OP differ more from the experimental values. The comparison between OP values obtained from explicit and rebuilt H atoms is useful for the validation of the scripts that determine OP values but I think we should accept the relatively small differences observed between both methods since they seem to come from the force field more than from errors in the scripts, do you agree?<br /><br />On the uncertainties for the OP values obtained for a trajectory: in my scriptI take the OP values for a given C atom over frames and lipids all together, then the average value, the standard deviation of the sample (STD), and finally the standard deviation of the mean (STD/sqrt(len(OP)) are determined. As far as I understand, the uncertainty of the final OP value with a confidence level of ~68% should be directly 2*STDmean. This value is typically very small (ridiculously small) for a reasonable trajectory (several thousand frames and hundreds of lipids). In contrast, the STD of the sample, determined as I did, is quite large. The other two scripts do the calculation in a different way: the average of the OP for each C-atom for all the lipids within each frame is first determined and then the average over these average values over all the frames and the corresponding STD are provided. A much smaller value is obtained for the standard deviation of the sample in this way but I do not know how this is justified. I guess I am missing something basic…Ángel Piñeirohttps://www.blogger.com/profile/09529986464564873340noreply@blogger.comtag:blogger.com,1999:blog-1283194685237019772.post-23028858230294233712019-09-08T15:11:12.421+03:002019-09-08T15:11:12.421+03:00Hi Samuli, I'm glad that buildH is useful! I w...Hi Samuli, I'm glad that buildH is useful! I will have a look to the tiny changes you made and will update the files. In the coming months we plan to i) use threads, ii) add tests and iii) release it on PyPI (so that it's installable via pip).<br />For the automatic topology detection, I have now a first draft code working for AA or UA POPC. The code reads a graph of POPC (that I built by hand) using generic atom names. Then it reads a single POPC molecule (in pdb format) and builds a new graph based on atom connectivities (calculated with distances). By using an isomorphism between the two graphs the code automatically generates a dictionnary like {'generic_name_atom1': 'pdb_name_atom1', ...}. So now, it'll be trivial to write all the files needed by our programs (such as dic_lipids.py) from that dictionnary.<br />So I'm trying to clean up this draft and add a little documentation ASAP. Unfortunately, I have very little time now since my teaching duties just restarted. But I do my best to push something on github (hopefully) early october.<br />Patrick Fuchshttps://www.blogger.com/profile/16910815195856162826noreply@blogger.comtag:blogger.com,1999:blog-1283194685237019772.post-32369897383085947792019-08-30T11:32:10.471+03:002019-08-30T11:32:10.471+03:00I have now analysed the united atom simulations fo...I have now analysed the united atom simulations for NMRlipids IVb using buildH and added the results into the manuscript:<br />https://github.com/NMRLipids/NMRlipidsIVPEandPG<br /><br />I also forked the code into the NMRlipids organization (https://github.com/NMRLipids/buildH) and committed the files that I needed for the analysis. I did not make a pull request because I only did the modifications required for the NMRlipids IVb project, and I was not sure if you want to have such incomplete files in the main version. If you want, you can of course pull it yourself.<br /><br />Thanks again for the code, this is really useful for the project! Also the automatic topology detection sounds highly interesting. How close it would be to get that included?Samuli Ollilahttps://www.blogger.com/profile/06106569992787533569noreply@blogger.comtag:blogger.com,1999:blog-1283194685237019772.post-58294935596514837992019-08-23T13:42:58.458+03:002019-08-23T13:42:58.458+03:00Here, I update the paths for the different dynamic...Here, I update the paths for the different dynamics both for the 150 mM NaCl and the new ones without ionic concentration:<br /><br />POPE<br />150 mM NaCl<br /> • LIPID17 https://zenodo.org/record/2577305#.XR9fUyZ7nes <br /> • GROMOS https://zenodo.org/record/2574491#.XR9fUiZ7nes <br /> • CHARMM36 https://zenodo.org/record/2577454#.XR9fOiZ7nes <br /> • SLIPIDS https://zenodo.org/record/2578069#.XR9fUCZ7nes <br />0 nM NaCl<br /> • LIPID17 https://zenodo.org/record/3237701#.XR9azCZ7nes <br /> • GROMOS https://zenodo.org/record/3237754#.XR9axSZ7nes <br /> • CHARMM36 https://zenodo.org/record/3237461#.XR9fNyZ7nes <br /> • SLIPIDS https://zenodo.org/record/3231342#.XR9atiZ7nes <br />POPG<br />150 mM NaCl<br /> • LIPID17 https://zenodo.org/record/2573905#.XR9fVyZ7nes <br /> • GROMOS https://zenodo.org/record/3257649#.XR9fOSZ7nes <br /> • CHARMM36 https://zenodo.org/record/2573531#.XR9fViZ7nes <br /> • SLIPIDS https://zenodo.org/record/3269320#.XR9atSZ7nes <br />0 nM NaCl<br /> • LIPID17 https://zenodo.org/record/3247659#.XR9aziZ7nes <br /> • GROMOS https://zenodo.org/record/3266166#.XR9awiZ7nes <br /> • CHARMM36 https://zenodo.org/record/3237463#.XR9fNiZ7nes <br /> • SLIPIDS https://zenodo.org/record/3364460#.XV_B__x7nes<br />POPC<br />150 mM NaCl<br /> • LIPID17 https://zenodo.org/record/2574959#.XR9fOyZ7nes <br /> • GROMOS https://zenodo.org/record/2574691#.XR9fVSZ7nes <br /> • CHARMM36 https://zenodo.org/record/2628335#.XR9fPSZ7nes <br /> • SLIPIDS https://zenodo.org/record/2574689#.XR9fUSZ7nes <br />0 nM NaCl<br /> • LIPID17 https://zenodo.org/record/3237657#.XR9azCZ7nes <br /> • GROMOS https://zenodo.org/record/3247435#.XR9awyZ7nes <br /> • CHARMM36 https://zenodo.org/record/3247813#.XR9fNSZ7nes <br /> • SLIPIDS https://zenodo.org/record/3235552#.XR9atyZ7nes <br />POPG-POPE 1:3<br />150 mM NaCl<br /> • LIPID17 https://zenodo.org/record/2579061#.XR9fTyZ7nes <br /> • GROMOS https://zenodo.org/record/2579063#.XR9fSyZ7nes <br /> • CHARMM36 https://zenodo.org/record/2579108#.XR9fSiZ7nes <br /> • SLIPIDS https://zenodo.org/record/2579224#.XR9fSCZ7nes <br />0 nM NaCl<br /> • LIPID17 https://zenodo.org/record/3241269#.XR9ayCZ7nes <br /> • GROMOS https://zenodo.org/record/3266238#.XR9awCZ7nes <br /> • CHARMM36 https://zenodo.org/record/3249960#.XR9fNCZ7nes <br /> • SLIPIDS https://zenodo.org/record/3238258#.XR9auCZ7nes <br />POPG-POPE 3:1<br />150 mM NaCl<br /> • LIPID17 https://zenodo.org/record/2579344#.XR9fRyZ7nes <br /> • GROMOS https://zenodo.org/record/2580158#.XR9fQyZ7nes <br /> • CHARMM36 https://zenodo.org/record/2580153#.XR9fRCZ7nes <br /> • SLIPIDS https://zenodo.org/record/2580153#.XR9fRCZ7nes <br />0 nM NaCl<br /> • LIPID17 https://zenodo.org/record/3256490#.XR9aySZ7nes <br /> • GROMOS https://zenodo.org/record/3267234#.XR9avSZ7nes <br /> • CHARMM36 https://zenodo.org/record/3250127#.XR9atCZ7nes https://zenodo.org/record/3258112#.XR9azyZ7nes <br /> • SLIPIDS https://zenodo.org/record/3240403#.XR9auyZ7nes <br />POPC-POPG 7:3<br />150 mM NaCl<br /> • LIPID17 https://zenodo.org/record/2585523#.XR9fQSZ7nes <br /> • GROMOS https://zenodo.org/record/2582721#.XR9fPiZ7nes <br /> • CHARMM36 https://zenodo.org/record/2580902#.XR9fQCZ7nes <br /> • SLIPIDS https://zenodo.org/record/2581186#.XR9fPyZ7nes <br />0 nM NaCl<br /> • LIPID17 https://zenodo.org/record/3241243#.XR9ayiZ7nes <br /> • GROMOS https://zenodo.org/record/3266240#.XR9aviZ7nes <br /> • CHARMM36 https://zenodo.org/record/3248689#.XR9fMyZ7nes <br /> • SLIPIDS https://zenodo.org/record/3240156#.XR9auSZ7nes Antonio PEÓNhttps://www.blogger.com/profile/06166908064718301905noreply@blogger.comtag:blogger.com,1999:blog-1283194685237019772.post-41304189566000033442019-08-23T12:52:04.635+03:002019-08-23T12:52:04.635+03:00Running a similar orden parameters analysis for ta...Running a similar orden parameters analysis for tails as with heads for the different all atoms force fields (CHARMM36, SLIPID) using explicit and rebuilt H atoms (implicit, Ángel’s script). For both hydrogens in a tail the analyses are very similar to the implicit hydrogens. <br />https://drive.google.com/file/d/1viLnB9UAIN8Kpo7LKvOxWbiBG512TfjP/view?usp=sharing<br />Antonio PEÓNhttps://www.blogger.com/profile/06166908064718301905noreply@blogger.comtag:blogger.com,1999:blog-1283194685237019772.post-49445635152348845852019-08-22T09:04:35.371+03:002019-08-22T09:04:35.371+03:00Thanks a lot Patrick!
Based on your report it se...Thanks a lot Patrick! <br /><br />Based on your report it seems to me that there is no essential difference between results from different codes used for united atom reconstruction, but there are differences in hydrogen locations with atomistic models. <br /><br />I think that the solution in NMRlipids IVb is to use slightly larger error bars for order parameters calculated from united atom simulations.<br /><br />I will try your code and start the analysis of the available data asap.Samuli Ollilahttps://www.blogger.com/profile/06106569992787533569noreply@blogger.comtag:blogger.com,1999:blog-1283194685237019772.post-55790472034077564182019-08-22T08:55:56.640+03:002019-08-22T08:55:56.640+03:00I think that this a good point. The glycerol backb...I think that this a good point. The glycerol backbone dihedrals are fitted also in CHARMM36 parametrization to optimize the order parameters, but in my understanding they do not fit dihedrals involving hydrogens. If hydrogens are not in ideal conformation in all atom simulations, the question is then if the ideal or twisted geometry of hydrogens is more realistic?Samuli Ollilahttps://www.blogger.com/profile/06106569992787533569noreply@blogger.comtag:blogger.com,1999:blog-1283194685237019772.post-51523622866947235022019-08-17T18:08:08.802+03:002019-08-17T18:08:08.802+03:00Hi,
In my opinion, it is noticeable the difference...Hi,<br />In my opinion, it is noticeable the difference found in the order parameters calculated with the script from Ángel Piñeiro for some of the hydrogens in carbons g2 and g3 in CHARMM and SLIPIDs (notice the graphics uploaded by Antonio Peón: https://drive.google.com/open?id=10tT0YLT4OHajEbsLNDP1XBLm71sKXtoO). It does not seem to be a problem with the reconstruction of implicit hydrogens done with the script, since the script’s reconstruction matches much better with the explicit hydrogens found in other force-fields (CHARMM, in particular). I think it is a more general problem that has the origin in the different torsion angles considered for the angles involving g2 and g3 in both force-fields (CHARMM and SLIPIDS). This difference was also observed in your previous manuscript (Figure 2 in https://pubs.acs.org/doi/pdf/10.1021/acs.jpcb.5b04878?rand=p3onqvgq). There you can notice that again Slipids is the force-field fitting worse with experimental order parameters calculated for DPPC and POPC for these carbons. <br />Curiously, if you go to the original paper of Slipids (https://pubs.acs.org/doi/pdf/10.1021/jp212503e?rand=mo3qdgjs), they say they take the parameters for all covalent bonds and angles, LJ and torsional parameters for the lipid head group from the original CHARMM36 (C36) FF (https://pubs.acs.org/doi/pdf/10.1021/jp101759q?rand=19u4avp3), BUT they also say they derived new parameters “because of the known difficulties of reproducing the vdW dispersion interaction by ab initio methods, experimental heats of vaporization and densities were used during the fitting of the LJ parameters. First, the new charges were used with the original parameters from C36 (LJ and torsional) and the LJ parameters were then altered to agree with experiments experimental results. After this, the torsional potentials were fitted from ab initio computations for the model compound. After one round in the parametrization scheme, it was necessary to refit the LJ parameters again and the torsional potentials until self-consistency was obtained.” So, it seems that people from Slipids “manually” changed some torsion parameters, affecting the order parameters for atoms g2 and g3, as our results (and previous results) are showing. <br />As a conclusion I would say that inconsistences between OP coming from explicit and rebuilt H atoms are due to the mis-geometrical-optimization of the explicit H atoms in the target forcefield.<br />Rebecahttps://www.blogger.com/profile/07699158975198107956noreply@blogger.comtag:blogger.com,1999:blog-1283194685237019772.post-14021581642931281092019-08-14T21:19:39.911+03:002019-08-14T21:19:39.911+03:00Hello all, I'm very sorry for my late answer b...Hello all, I'm very sorry for my late answer but it took me longer than expected. Especially because I wanted to make a thorough comparison / validation of my program for rebuilding hydrogens / calculating OPs vs the programs from Josef and Angel. And also because I wanted to add docstrings and an understandable documentation.<br /><br />So my program is called "buildH" and the repository is on github: https://github.com/patrickfuchs/buildH. I didn't upload any script on nmrlipids github repo for now since more than one file is required, I thought it would be easier like that. If you use it, please do some regular "git pull" since we will work on it on august / september. There are still many little things to do, but the core of the program is working.<br />Importantly, I made a thorough validation of buildH which is also on github (https://github.com/patrickfuchs/buildH/tree/master/CHARMM36_POPC_validation, see report_buildH.ipynb or report_buildH.pdf). For easier comparison, buildH outputs the OPs in the format of the program from Josef *and* that from Angel. For the validation, I used the strategy like in the paper from Thomas Piggot: starting from an all-atom (POPC CHARMM36) traj, I remove the Hs, reconstruct them and calculate the OPs and compare them to the initial traj with Hs (OPs calculated with Josef's prog). You will see in the report that both Angel program and buildH agree very well (except for a few details). I also put there some links on output files with numbers in case they are useful. Don't hesitate if you have questions / comments.<br /><br />For the other script building automatically the required files for buildH, I didn't have time to complete it.<br /><br />In the near future, we plan to add some improvements to buildH (such as multithreading, unit tests, automatic topology detection using the other script, etc).Patrick Fuchshttps://www.blogger.com/profile/16910815195856162826noreply@blogger.comtag:blogger.com,1999:blog-1283194685237019772.post-84752877504700564512019-07-23T10:30:09.118+03:002019-07-23T10:30:09.118+03:00Hi Samuli. Yes, I meant that... but we should prob...Hi Samuli. Yes, I meant that... but we should probably do the same comparison we did with the heads. For the tails we just did a couple of plots.Ángel Piñeirohttps://www.blogger.com/profile/09529986464564873340noreply@blogger.comtag:blogger.com,1999:blog-1283194685237019772.post-80893043230495672722019-07-20T16:14:52.379+03:002019-07-20T16:14:52.379+03:00Thanks Antonio and Angel for the data. I think tha...Thanks Antonio and Angel for the data. I think that we should do similar analysis with the code by Fuchs et al. once it is available, and then decide how we calculate the united atom order parameters for the publication and determine the error bars. The numerical data from your plots would be useful for this. Samuli Ollilahttps://www.blogger.com/profile/06106569992787533569noreply@blogger.comtag:blogger.com,1999:blog-1283194685237019772.post-89493536884960903162019-07-20T16:07:41.455+03:002019-07-20T16:07:41.455+03:00"I am commenting here just the OP for the lip..."I am commenting here just the OP for the lipid heads since the results for the tails were almost indistinguishable." <br /><br />Let me know if I understood this correctly. Do you mean with this that the order parameters from regenerated and all-atom hydrogens are essentially identical for tails, but not for glycerol backbone and headgroup?<br />Samuli Ollilahttps://www.blogger.com/profile/06106569992787533569noreply@blogger.comtag:blogger.com,1999:blog-1283194685237019772.post-6231886730619848402019-07-17T14:03:43.028+03:002019-07-17T14:03:43.028+03:00Hello, I have compared the head orden parameters f...Hello, I have compared the head orden parameters for the different all atoms force fields (CHARMM36, SLIPID) using explicit and rebuilt H atoms. The analysis was performed over the last 100 ns of the trajectories (500 ns long with 500 lipids) with POPE, POPC and POPG bilayers in 150 mM NaCl. You can see the OP for both explicit Hs(dot) and rebuilt Hs ( line, calculated using the Ángel's code) giving very similar results (within the statistical uncertainties). The biggest difference is in hydrogen g2 and specifically in SLIPID. Link: https://drive.google.com/open?id=10tT0YLT4OHajEbsLNDP1XBLm71sKXtoOAntonio PEÓNhttps://www.blogger.com/profile/06166908064718301905noreply@blogger.comtag:blogger.com,1999:blog-1283194685237019772.post-2975099285823454042019-07-08T19:30:10.604+03:002019-07-08T19:30:10.604+03:00Sorry for taking so much time. I have done a compa...Sorry for taking so much time. I have done a comparison of the OP using the explicit and rebuilt H atoms for Charmm36 POPG and for slipids POPC using my code. I took two of the trajectories recently uploaded by Antonio Peón and Rebeca García (see the list just above). In order to directly see the impact of the H-rebuilding I determined the OP for a single lipid, for a single frame and also the average OP over 20 ns. You can see the results here: http://smmb.usc.es/docs/OPcomparisonAAvsUA.pdf. I have also generated a pdb file of POPC with the rebuilt H atoms of the lipid head over the original explicit H atoms (http://smmb.usc.es/docs/POPC_rebuiltHs.pdb). For the rebuilt ones I used X as atomname in the pdb, to facilitate the selection and comparison with a molecular viewer (note that the the bond distance for my rebuilt H atoms is exactly 1 A, since this is irrelevant for the OP calculation). I am commenting here just the OP for the lipid heads since the results for the tails were almost indistinguishable.<br /><br />Some comments:<br />As you can see in the first column of the plots (http://smmb.usc.es/docs/OPcomparisonAAvsUA.pdf) the impact of my rebuilt H atoms in the OP is significant for some positions and almost negligible for others. G1-G3 seem to be more sensitive than alpha and beta. <br />If the differences are just random noise, they are expected to be diluted when averaging over many lipids. Actually this is visible in the second and third columns (average over a single frame with 500 lipid molecules and average over 20 ns for the same system).<br />The results for charmm-POPG are reasonably good for all the atoms but for slipids-POPC the difference between explicit and rebuilt H atoms is significant for G3S and for G2R. The reason of the differences is not the method employed for the calculation since I used the same code. I do not know if the rest of the methods behave similarly. I will try to repeat the same analysis for other trajectories (lipids and forcefields).Ángel Piñeirohttps://www.blogger.com/profile/09529986464564873340noreply@blogger.comtag:blogger.com,1999:blog-1283194685237019772.post-31439729949728755122019-06-26T16:19:41.299+03:002019-06-26T16:19:41.299+03:00Hi,
we leave here the routes of the MD trajectorie...Hi,<br />we leave here the routes of the MD trajectories we have uploaded to the respository, indicating the membrane composition of each one and the different force-fields used. All of them were carried out using a concentration of 150 mM NaCl (we are now carrying out analogous MD simulations but using just the minimum amount of ions to neutralize the systems):<br /><br /><br />POPC-POPG 7:3<br /><br />Gromos-CKP<br /><br />https://zenodo.org/record/2582721#.XIjZg4XgrLA<br /><br />SLIPID<br /><br />https://zenodo.org/record/2581186#.XIjZhIXgrLA<br /><br />CHARMM36<br /><br />https://zenodo.org/record/2580902#.XIjZhYXgrLA<br /><br />LIPID17<br /><br />https://zenodo.org/record/2585523#.XIjZhoXgrLA<br /><br /><br />POPG-POPE 3:1<br /><br />Gromos-CKP<br /><br />https://zenodo.org/record/2580158#.XIjZh4XgrLA<br /><br />SLIPID<br /><br />https://zenodo.org/record/2579675#.XIjZiYXgrLA<br /><br />CHARMM36<br /><br />https://zenodo.org/record/2580153#.XIjZiIXgrLA<br /><br />LIPID17<br /><br />https://zenodo.org/record/2579344#.XIjZioXgrLA<br /><br /><br />POPG-POPE 1:3<br /><br />Gromos-Ckp<br /><br />https://zenodo.org/record/2579063#.XIjZjoXgrLA<br /><br />SLIPID<br /><br />https://zenodo.org/record/2579224#.XIjZi4XgrLA<br /><br />CHARMM36<br /><br />https://zenodo.org/record/2579108#.XIjZjYXgrLA<br /><br />LIPID17<br /><br />https://zenodo.org/record/2579061#.XIjZj4XgrLA<br /><br /><br />POPE<br /><br />Gromos-CKP<br /><br />https://zenodo.org/record/2574491#.XIjZk4XgrLA<br /><br />SLIPID<br /><br />https://zenodo.org/record/2578069#.XIjZkYXgrLA<br /><br />CHARMM36<br /><br />https://zenodo.org/record/2577454#.XIjZlIXgrLA<br /><br />LIPID17<br /><br />https://zenodo.org/record/2577305#.XIjZlYXgrLA<br /><br /><br />POPG<br /><br />Gromos-CKP<br /><br />https://zenodo.org/record/2572897#.XIjZnYXgrLA<br /><br />SLIPID<br /><br />https://zenodo.org/record/2562853#.XIjZnIXgrLA<br /><br />CHARMM36<br /><br />https://zenodo.org/record/2573531#.XIjZm4XgrLA<br /><br />LIPID17<br /><br />https://zenodo.org/record/2573905#.XIjZmoXgrLA<br /><br /><br />POPC<br /><br />Gromos-CKP<br /><br />https://zenodo.org/record/2574691#.XIjZloXgrLA<br /><br />SLIPID<br /><br />https://zenodo.org/record/2574689#.XIjZkoXgrLA<br /><br />CHARMM36<br /><br />https://zenodo.org/record/2575587#.XIjZlYXgrLA<br /><br />LIPID17<br /><br />https://zenodo.org/record/2574959#.XIjZl4XgrLA<br /><br />Rebecahttps://www.blogger.com/profile/07699158975198107956noreply@blogger.comtag:blogger.com,1999:blog-1283194685237019772.post-57404141985014496192019-06-26T08:17:44.016+03:002019-06-26T08:17:44.016+03:00I am looking now the figure 6 from the recent pape...I am looking now the figure 6 from the recent paper by Tom Piggot (http://dx.doi.org/10.1021/acs.jctc.7b00643). If I understand it correctly, they make a similar comparison as you have done (comparison between regenerated and original hydrogens in all atom trajectory). In their results, NMRlipids protonate approach does not give exactly the same values as the original trajectory. Could you comment or send a figure about the results from your code, are they closer to originals than NMRlipids protonate in the figure by Piggot et al.? I have not read the paper very carefully, but maybe Tom could have some insight why the results from different protonation methods are slightly different?<br /><br />It may be that we just have to put larger error bars to the united atom results which include the results from different hydrogen generation methods, although my feeling is that the results closer the original all atom trajectories would be better. Samuli Ollilahttps://www.blogger.com/profile/06106569992787533569noreply@blogger.comtag:blogger.com,1999:blog-1283194685237019772.post-48607450911102773762019-06-25T19:39:13.723+03:002019-06-25T19:39:13.723+03:00As far as I see the three codes (Josef, Patrick an...As far as I see the three codes (Josef, Patrick and mine) provide exactly the same OP values once the Hydrogens are rebuilt. There is a difference in the STEM values (my calculation is a bit more conservative than yours) but this can be corrected easily, if necessary. At the moment I do not see a reason to change my calculation. The only possible source of errors that I can imagine is in the rebuilding of Hydrogens. This is more tricky than the determination of OP. Since my code can ignore explicit H atoms (using the "-a ignH" option) I checked this for AA trajectories with explicit and rebuild Hydrogens. The results were not identical (as expected) but very similar. As explained before, I validated my rebuilding of H atoms by comparing the explicit Hydrogens in different AA trajectories with those rebuilt by using it.<br /><br />@Patrick, in case of doubts perhaps it is useful to compare the output of your code for UA with that coming from my code... although I understand that your code has already been well validated and it is robust. In the end I guess your code will be the best one since it will not need any input other than the trajectory.Ángel Piñeirohttps://www.blogger.com/profile/09529986464564873340noreply@blogger.comtag:blogger.com,1999:blog-1283194685237019772.post-46429068417292880072019-06-19T20:40:23.989+03:002019-06-19T20:40:23.989+03:00Hi Samuli, thanks for the precisions. So I ran my ...Hi Samuli, thanks for the precisions. So I ran my code on the full trajectory, and the mean absolute difference between the OP calculated with my script and the awk one is 0.002. It might be even less because it is possible that some hydrogens were mixed (in the database, for the case we have 2 or 3 hydrogens, I didn't check whether my code outputs these hydrogens in the same order, that's why it's useful to have the name of each H as in Josef's script). I can put the results here if you want, let me know. I do my best to make my code available as soon as I can.Patrick Fuchshttps://www.blogger.com/profile/16910815195856162826noreply@blogger.com