<br><br>----- Original Message -----<br>From: ms <devicerandom@gmail.com><br>Date: Monday, October 18, 2010 8:23<br>Subject: [gmx-users] The trouble with dihedral restraints: frozen peptide backbones<br>To: Discussion list for GROMACS users <gmx-users@gromacs.org><br><br>> Hi,<br>> <br>> It's a bit long but bear with me if you can. I'm getting quite <br>> mad with this. Thanks :)<br>> <br>> Now. I am using Gromacs (4.0.7 currently) to create a custom <br>> coarse-grained, no-solvent MD peptide model -one that contains <br>> the backbone but has only C-alphas, no side chains.<br>> <br>> To enforce chirality in such a toy, I thought that a simple <br>> (naive?) but functional idea could be that of enforcing a not-<br>> too-hard and wide-bottomed dihedral restrain on the phi angle, <br>> like that:<br>> <br>> [ dihedral_restraints ]<br>> ; ai aj ak <br>> al type label phi dphi kfac power<br>> ; phi C'(n-1) - N - CA - C'<br>> 3 <br>> 5 7 8 1 1 -90 75 2 1<br><br>Try mdrun -debug on a short run. What does it have to say about the parameter values for this dihedral restraint? Likewise gmxdump -s on the .tpr (though you will have to search around for F_DIHRES among other things, probably).<br><br>Mark<br><br>> <br>> so that the system is free to move in a large phi-angle space <br>> without problems, but cannot go in the forbidden, symmetrical <br>> glycine-like space. In testing, I set the force constant small <br>> (10) so that the chain, even if with some difficulty, should be <br>> able to bend and move by spending some energy in the "forbidden" <br>> zones, if needed to achieve structure, etc.<br>> <br>> Now, the point is that *when* it works, it works quite well. But <br>> sometimes -with "sometimes" being from 0% to 90% of times, <br>> depending on details of what I am simulating- it totally <br>> doesn't. That is: sometimes, starting from an extended <br>> structure, the chain simply freezes. It is not static: it <br>> wiggles a bit here and there, but it stays mostly extended. <br>> Energies plots (g_energy output xvg's) are on average flat: that <br>> is, they show the right noise for all components, but no <br>> structure: nothing happens, for as much as 50 ns.<br>> <br>> Sometimes instead it works, and then you see that everything <br>> works alright; the peptide structures, potential energy goes <br>> down and down, etc.etc.<br>> <br>> There are no differences in minimization energies etc. from a <br>> simulation that freezes and one that works. Moving the restraint <br>> center doesn't help much. I *know* it's the restraints because <br>> if I keep everything the same but I lower the force constant or <br>> remove restraints altogether, the freezing bug happens no more -<br>> but then I have little or no more chirality, of course. I have <br>> the following further facts:<br>> - It tends to happen more often for longer chains: 10 or 20-<br>> residue long chains work often OK, but a 40-long chain is <br>> difficult to get "right".<br>> - It *never* happens, to my knowledge, if I put more than 1 <br>> chain to cluster in the box.<br>> - The Ramachandran plot of the extended frozen thing is "right", <br>> meaning that it populates a reasonable (mostly beta) angle <br>> space. During trajectory angles move quite a bit but stay in the <br>> right (top-left) zone of the plot.<br>> - "Good"-behaving simulation often start with a short "frozen" <br>> stretch: then the dihedral restraint energy suddenly jumps to <br>> some above noise positive value, to accomodate some bending I <br>> guess, and all goes downhill from there.<br>> - In frozen simulations, there's almost no LJ energy: the chain <br>> does not self-interact.<br>> <br>> All of this, added to the apparent "randomness" of the <br>> pathology, makes me think that the chain simply is stuck in a <br>> situation where it has a too high barrier to properly bend -a <br>> barrier which is not a "soft" one as the harmonic potential of <br>> the restrain would make me think. This would explain why putting <br>> more peptides and letting them bind together avoids the problem: <br>> the potential energy of binding compensates for the bending when <br>> you have more than one thing; but sometimes a single chain <br>> cannot find itself to interact with, because it should bend <br>> correctly to do so.<br>> <br>> And then, I read on the Gromacs website *this*<br>> http://www.gromacs.org/Documentation/How-tos/Dihedral_Restraints<br>> " 2. The manual is a bit unclear about whether this type of <br>> dihedral restraint is stable for use near 180 degrees. Chris <br>> Neale has found that everything appears to behave normally and <br>> as expected over the entire range of dihedral angles including <br>> 180 degrees. However, one must avoid the situation in which the <br>> actual dihedral is close to 180deg away from the restrained dihedral."<br>> <br>> And I know that somehow the potential is discontinuous and/or <br>> anyway not intended to be well-behaved everywhere:<br>> http://www.mail-archive.com/gmx-users@gromacs.org/msg12353.html<br>> <br>> "> On Tue, Feb 26, 2008 at 11:37 AM, David Mobley <[EMAIL <br>> PROTECTED]> wrote:<br>> > > The other way of putting what Mark said is that phi is <br>> only meaningful on some range (-pi to pi, or 0 to 2pi, depending <br>> on how you define it)<br>> > > and so what you require is that the potential be <br>> harmonic for the<br>> > > region in which phi is meaningful. You don't care what <br>> happens outside that. Or, in this case, you handle the issue by <br>> mapping phi values outside the allowed range back into <br>> that allowed range."<br>> <br>> (I don't get what does it mean to "map back" phi values, by the way).<br>> <br>> I am no MD expert, so maybe my hypothesis is totally off, but I <br>> am reasonably sure that either I do not understand correctly how <br>> dihedral restraints work, or I am suffering by the discontinuity <br>> of the way the potential is implemented (Or both). I would like <br>> to know:<br>> - How can I try to understand in full what is going on<br>> - If my diagnosis sounds right, how can I start to correct the <br>> problem.- If the problem is my naive approach to chirality, any <br>> suggestion is welcome.<br>> <br>> Thanks A LOT if you've been patient enough to read until here! <br>> And thanks even *more* if you can help!<br>> <br>> Yours truly,<br>> Massimo<br>> -- <br>> gmx-users mailing list gmx-users@gromacs.org<br>> http://lists.gromacs.org/mailman/listinfo/gmx-users<br>> Please search the archive at <br>> http://www.gromacs.org/Support/Mailing_Lists/Search before posting!<br>> Please don't post (un)subscribe requests to the list. Use the <br>> www interface or send it to gmx-users-request@gromacs.org.<br>> Can't post? Read http://www.gromacs.org/Support/Mailing_Lists