[gmx-users] Unstable Minimizations

Tsjerk Wassenaar tsjerkw at gmail.com
Mon Jan 11 09:14:28 CET 2010


Hi Jack, Justin,

Is this related to last weeks post? I suggested that you might have
seen long bond warnings, which in this case certainly show up. Lots of
them. But one stands out: there's a break in the chain of 3.3 nm!
During minimization an attempt is made to bring these ends together,
which causes failure. Probably pdb2gmx missed that there should be two
chains rather than one.

When designing a workflow, or adapting one for a new situation, you
have to start stepwise. Well, trying an existing workflow is okay, but
if it fails, you'll have to trace steps and check for anything that is
out of the ordinary. You're stop is the step with pdb2gmx. Something
is not working there. Make sure to get to know what pdb2gmx requires
in order to process a structure and make sure that you make the input
comply with these requirements. The main requirements are:

1. every chain needs to have a unique chain identifier
2. there should be no missing atoms (other than hydrogens)
3. there should be no missing residues , except possibly at the termini.
4. there should be no residues (ligands) that are not in the building
block database for the force field you're using
5. atom names in the coordinate file should be consistent with those
in the building block database
6. there should be no atoms in residues that are not listed in the
building block entry, except possibly for hydrogen atoms, which can be
stripped using the -ignh flag

Like I replied last week, your problem probably relates to 1. Some
programs (CNS/X-PLOR) use the segment identifier field at the end of
the line, rather than the chain identifier more in the middle. You can
fix that with sed or python (the latter being more readable ;)):

sed -e '{/^\(ATOM  \|HETATM\)/
s/^\(.\{21\}\).\(.\{50\}\)\(.\)/\1\3\2\3/}' in.pdb > out.pdb

python -c "for i in [ l for l in open('in.pdb') if
l.startswith('ATOM') ]: print i[:21]+i[71]+i[22:]," > out.pdb

I'm not exactly sure that 71 (72) in the python line is the right
field for the segid; you'd better check the PDB specification for
that. The sed line should be correct as is.

Cheers,

Tsjerk
Cheers,

Tsjerk
On Mon, Jan 11, 2010 at 3:58 AM, Justin A. Lemkul <jalemkul at vt.edu> wrote:
>
>
> On 1/10/10 9:47 PM, Jack Shultz wrote:
>>
>> Thanks Justin,
>> I went back to the original pdb files. These were conformations of the
>> same protein derived from molecular dynamics simulations performed by
>> Andrey.
>> What I intially attempted was preping the structures using tleap, hoping
>> to paint in missing atoms for residues. Then use this to replace
>
> Well, it seems that you may be hoping for too much :)  Your log file shows a
> whole bunch of failures that look to be related to some early processing of
> your structure, and other warnings about close contacts detected in tleap.
>
> I think you may need to start with an actual intact structure, or else coax
> your preparation steps to make this happen.  I am not too familiar with
> tleap and sleap, do they magically fix missing atoms?
>
> -Justin
>
>> non-standard residues
>> sed s/PRO\ A\ \ \ 1/NPROA\ \ \ 1/g fzd2_md7-8_c6_cc.pdb | sed s/PRO\ B\
>> \ \ 1/NPROB\ \ \ 1/g | sed s/PHE\ A\ \ 99/CPHEA\ \ 99/g | sed s/PHE\ B\
>> \ 99/CPHEB\ \ 99/g | sed s/O\ \ \ CPHE/OC1\ CPHE/g | sed s/OXT\
>> CPHE/OC2\ CPHE/g | sed s/HIS\ /HID\ /g | sed s/LYS\ /LYP\ /g | sed
>> s/CYS\ /CYN\ /g > protein2.pdb
>> Then fixed the nterminal residue name. Finally replaced all CYS to CYS2
>> I went back and did the same thing except for tleap. It pdb2gmx seems to
>> process these files without needing the tleap step.
>> Still I see some of the same lincs errors.
>> rms 10.669050, max 173.182678 (between atoms 1857 and 1859)
>> rms 10.669803, max 173.177811 (between atoms 1857 and 1859)
>> rms 10.670179, max 173.175400 (between atoms 1857 and 1859)
>> rms 10.670368, max 173.174149 (between atoms 1857 and 1859)
>> rms 10.670460, max 173.173553 (between atoms 1857 and 1859)
>> rms 10.670508, max 173.173141 (between atoms 1857 and 1859)
>> rms 10.670531, max 173.173035 (between atoms 1857 and 1859)
>> rms 10.670543, max 173.172958 (between atoms 1857 and 1859)
>> rms 10.670549, max 173.172928 (between atoms 1857 and 1859)
>> rms 10.670552, max 173.172913 (between atoms 1857 and 1859)
>> rms 10.670554, max 173.172913 (between atoms 1857 and 1859)
>> rms 10.670554, max 173.172913 (between atoms 1857 and 1859)
>> ATOM   1857  CA  HIE   120      43.362  28.084  25.727  1.00  0.00
>> ATOM   1858  HA  HIE   120      43.677  27.135  25.748  1.00  0.00
>> ATOM   1859  CB  HIE   120      42.112  28.226  24.788  1.00  0.00
>> also this atom consistently has a very high Fmax
>> Step=    3, Dmax= 1.4e-02 nm, Epot=  1.45860e+10 Fmax= 2.82224e+12,
>> atom= 19392
>> Step=    4, Dmax= 7.2e-03 nm, Epot=  1.45396e+10 Fmax= 2.82207e+12,
>> atom= 19392
>> Step=    5, Dmax= 3.6e-03 nm, Epot=  1.45106e+10 Fmax= 2.82194e+12,
>> atom= 19392
>> Step=    6, Dmax= 1.8e-03 nm, Epot=  1.44953e+10 Fmax= 2.82181e+12,
>> atom= 19392
>> Step=    7, Dmax= 9.0e-04 nm, Epot=  1.44887e+10 Fmax= 2.82196e+12,
>> atom= 19392
>> Step=    8, Dmax= 4.5e-04 nm, Epot=  1.44850e+10 Fmax= 2.82196e+12,
>> atom= 19392
>> Step=    9, Dmax= 2.2e-04 nm, Epot=  1.44832e+10 Fmax= 2.82196e+12,
>> atom= 19392
>> Step=   10, Dmax= 1.1e-04 nm, Epot=  1.44822e+10 Fmax= 2.82196e+12,
>> atom= 19392
>> Step=   11, Dmax= 5.6e-05 nm, Epot=  1.44818e+10 Fmax= 2.82196e+12,
>> atom= 19392
>> Step=   12, Dmax= 2.8e-05 nm, Epot=  1.44815e+10 Fmax= 2.82196e+12,
>> atom= 19392
>> Step=   13, Dmax= 1.4e-05 nm, Epot=  1.44814e+10 Fmax= 2.82196e+12,
>> atom= 19392
>> Its not clear to me what we should do to correct this structures...maybe
>> Andrey has some input.
>> http://boinc.drugdiscoveryathome.com/em_restrained_rcs_mdrun2.txt
>> On Sun, Jan 10, 2010 at 5:37 PM, Justin A. Lemkul <jalemkul at vt.edu
>> <mailto:jalemkul at vt.edu>> wrote:
>>
>>
>>
>>    On 1/10/10 5:18 PM, Jack Shultz wrote:
>>
>>        I am trying to get this workflow opperational. However, my
>>        systems are
>>        getting unstable. I have preped two mdp files: 1) one for
>>        restrained 2)
>>        unrestrained. LINCS errors appear for restrained and
>>        unrestrained has
>>        infinite energy appearing.
>>        http://boinc.drugdiscoveryathome.com/_em_restrained_rcs_mdrun.txt_
>>
>>        <http://boinc.drugdiscoveryathome.com/em_restrained_rcs_mdrun.txt>
>>        __
>>
>>
>>    This log file shows several "long bond" warnings, which may be the
>>    root of your problem.  See here:
>>
>>
>>  http://www.gromacs.org/Documentation/Errors#Long_bonds_and.2for_missing_atoms
>>
>>    Since your minimization is failing immediately, there is something
>>    physically unreasonable about your structure, such that EM cannot
>>    resolve the problem. Note, too, that one of the long bond warnings
>>    pertained to atom 1668, which is the location of the first LINCS
>>    warning.  Coincidence?  Not likely.  Re-examine the starting
>>    structure and figure out if anything is missing or poorly
>>    reconstructed (e.g., from initially missing atoms).
>>
>>
>>        This is where I get the LINCS Warnings
>>        Step -1, time -0.001 (ps)  LINCS WARNING
>>        relative constraint deviation after LINCS:
>>        rms 0.461520, max 14.428611 (between atoms 1668 and 1669)
>>        bonds that rotated more than 30 degrees:
>>        atom 1 atom 2  angle  previous, current, constraint length
>>        Steepest Descents:
>>            Tolerance (Fmax)   =  1.00000e+04
>>            Number of steps    =          100
>>        Warning: 1-4 interaction between 1658 and 1672 at distance 2.655
>>        which
>>        is larger than the 1-4 table size 2.400 nm
>>        These are ignored for the rest of the simulation
>>        This usually means your system is exploding,
>>        if not, you should increase table-extension in your mdp file
>>        or with user tables increase the table size
>>        Step=    0, Dmax= 1.0e-02 nm, Epot=  1.09364e+09 Fmax= 2.21154e+11,
>>        atom= 3292
>>        Step 1, time 0.001 (ps)  LINCS WARNING
>>        relative constraint deviation after LINCS:
>>        rms 0.680509, max 22.293625 (between atoms 1668 and 1670)
>>        bonds that rotated more than 30 degrees:
>>        atom 1 atom 2  angle  previous, current, constraint length
>>        What is a reasonable increase in table-extension. Is this a
>>        mis-leading
>>        suggestion?
>>
>>
>>    You should not adjust the table-extension.  The other part of the
>>    error message is what you need to pay attention to ("your system is
>>    exploding").
>>
>>    -Justin
>>
>>        Here is the log from the unrestrained minimization.
>>        http://boinc.drugdiscoveryathome.com/_em_rcs_mdrun.txt_
>>
>>        <http://boinc.drugdiscoveryathome.com/em_rcs_mdrun.txt>
>>        Here is a zip archive containing the working directory for this
>>        minimization. Its about 428 kb
>>
>>  http://boinc.drugdiscoveryathome.com/rcs_ga_run_10_bt_Fzd2-MD7-MD8-7.zip_lig_24205_ChemDiv_5754-2873_ts_1263004110202172000.zip
>>
>>        --
>>        Jack
>>
>>        http://drugdiscoveryathome.com <http://drugdiscoveryathome.com/>
>>        http://hydrogenathome.org <http://hydrogenathome.org/>
>>
>>
>>    --
>>    ========================================
>>
>>    Justin A. Lemkul
>>    Ph.D. Candidate
>>    ICTAS Doctoral Scholar
>>    MILES-IGERT Trainee
>>    Department of Biochemistry
>>    Virginia Tech
>>    Blacksburg, VA
>>    jalemkul[at]vt.edu <http://vt.edu/> | (540) 231-9080
>>    http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin
>>
>>    ========================================
>>    --
>>    gmx-users mailing list gmx-users at gromacs.org
>>    <mailto:gmx-users at gromacs.org>
>>    http://lists.gromacs.org/mailman/listinfo/gmx-users
>>    Please search the archive at http://www.gromacs.org/search before
>>    posting!
>>    Please don't post (un)subscribe requests to the list. Use the www
>>    interface or send it to gmx-users-request at gromacs.org
>>    <mailto:gmx-users-request at gromacs.org>.
>>    Can't post? Read http://www.gromacs.org/mailing_lists/users.php
>>
>>
>>
>>
>> --
>> Jack
>>
>> http://drugdiscoveryathome.com
>> http://hydrogenathome.org
>
> --
> ========================================
>
> Justin A. Lemkul
> Ph.D. Candidate
> ICTAS Doctoral Scholar
> MILES-IGERT Trainee
> Department of Biochemistry
> Virginia Tech
> Blacksburg, VA
> jalemkul[at]vt.edu | (540) 231-9080
> http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin
>
> ========================================
> --
> gmx-users mailing list    gmx-users at gromacs.org
> http://lists.gromacs.org/mailman/listinfo/gmx-users
> Please search the archive at http://www.gromacs.org/search before posting!
> Please don't post (un)subscribe requests to the list. Use the www interface
> or send it to gmx-users-request at gromacs.org.
> Can't post? Read http://www.gromacs.org/mailing_lists/users.php
>



-- 
Tsjerk A. Wassenaar, Ph.D.

Computational Chemist
Medicinal Chemist
Neuropharmacologist



More information about the gromacs.org_gmx-users mailing list