normal DNA Junction Blowing Up

  • jonasori
  • jonasori's Avatar Topic Author
  • Offline
  • Fresh Boarder
More
7 years 9 months ago #5670 by jonasori
DNA Junction Blowing Up was created by jonasori
Hello everyone,

I'm trying to run a MARTINI simulation of a stacked DNA junction and am running into some problems. Fair warning, I'm quite new to doing MD sims, but I have been chasing this problem for a couple days now. I am working from the Martini DNA tutorial , modifying each command as I see appropriate.

The problem itself is, essentially, that the simulation is blowing up when I try to run a short-ish energy minimization (1000 steps at dt = 0.001). The first trouble comes in the first grompp (at step 6), when it dies by "Fatal Error: Too many warnings." A selection of these warnings are as follows:

Warning: atom name 102 in cg-junction_original.top and bw.gro does not match (BB2 - BB1)
Warning: atom name 103 in cg-junction_original.top and bw.gro does not match (BB3 - BB2)
Warning: atom name 104 in cg-junction_original.top and bw.gro does not match (SC1 - BB3)
Warning: atom name 105 in cg-junction_original.top and bw.gro does not match (SC2 - SC1)
Warning: atom name 106 in cg-junction_original.top and bw.gro does not match (SC3 - SC2)
Warning: atom name 107 in cg-junction_original.top and bw.gro does not match (BB1 - SC3)
Warning: atom name 108 in cg-junction_original.top and bw.gro does not match (BB2 - BB1)
Warning: atom name 109 in cg-junction_original.top and bw.gro does not match (BB3 - BB2)
Warning: atom name 110 in cg-junction_original.top and bw.gro does not match (SC1 - BB3)
Warning: atom name 111 in cg-junction_original.top and bw.gro does not match (SC2 - SC1)
Warning: atom name 112 in cg-junction_original.top and bw.gro does not match (SC3 - SC2)
Warning: atom name 113 in cg-junction_original.top and bw.gro does not match (BB1 - SC3)
Warning: atom name 114 in cg-junction_original.top and bw.gro does not match (BB2 - SC4)
Warning: atom name 115 in cg-junction_original.top and bw.gro does not match (BB3 - BB1)
Warning: atom name 116 in cg-junction_original.top and bw.gro does not match (SC1 - BB2)
Warning: atom name 117 in cg-junction_original.top and bw.gro does not match (SC2 - BB3)
Warning: atom name 118 in cg-junction_original.top and bw.gro does not match (SC3 - SC1)
Warning: atom name 119 in cg-junction_original.top and bw.gro does not match (BB1 - SC2)
Warning: atom name 120 in cg-junction_original.top and bw.gro does not match (BB2 - SC3)
Warning: atom name 121 in cg-junction_original.top and bw.gro does not match (BB3 - SC4)
(more than 20 non-matching atom names)

WARNING 1 : 3468 non-matching atom names atom names from cg-junction_original.top will be used atom names from bw.gro will be ignored I know that there is a -maxwarn option, but I don't want to use it. When I do use it (to allow the grompp to finish), I just run into the same problem on the next step (the mdrun) when Fmax refuses to get under ~XXe03. My configuration: My .top file That is the problem. Here is what I have tried/my thought process: The easy, obvious pattern in the problem is that the atoms being compared are offset by 1; that is to say, BB2 is trying to match with BB1, but then BB3 is trying to match with BB2, and so if I constructed an atom to align those, the problem should be gone. However, that just pushes the problem elsewhere. I'm also suspicious of correcting these files by hand, since the two files are involved are bw.gro (a box of water, constructed automatically with editconf) and cg-junction.top, which is calling a file that is generated by the provided martini-dna.py (called NucleicA+Nucleic_B.itp). Since both of these files involved are automatically generated, it seems unlikely that there are simple computational errors like a missing atom (row). So tracking down the errors by hand seems both tedious and counter-productive. Furthermore, even when I do go in by hand and try to fix it, it doesn't get immediately fixed. I'm not really sure what else could be going on here. I have cleaned my input .pdb file well, and visual inspection in VMD shows that both the original file and the coarse-grained version of it both look as they should. When I created the box of water, I did not put any -maxsol in, which should allow the box to fill as much as possible. For the cg-junction.top file, I checked how many atoms were added in the above step, subtracted the 407 atoms generated by the Nucleic_A ... .itp files, took off another 132 atoms to make ions out of (132 since I have two intersecting arms of two 34-base strands, so ((34-1)*2)*2 = 33*4 = 132), and ~10% of the remaining beads to place as antifreeze (WF). I put those numbers into my cg-junction.top file, as they do in the tutorial's step 5. So my cg-junction_original.top file (which grompp doesn't like) looks like this: [ system ] ; name Martini system from j3s-cleaned.gro [ molecules ] ; name number Nucleic_A+Nucleic_B 1 Nucleic_C 1 Nucleic_D 1 W 23955 WF 2660 NA 132 If you've made it this far, I deeply appreciate your willingness to stick around and continue reading. So my obvious question is how do I stop this system from blowing up? However, I'm also realizing I don't understand much of what's really going on in the martinize-dna.py file, and by extension, what the columns in its output files (Nucleic_X.itp) stand for. If there is anything else you want to know about my files, just let me know. I put down what I thought was relevant, but as a beginner, my sense of relevance is probably pretty far off. Any help is greatly appreciated. Many thanks. Jonas[file cg-junction_original.top, line 21]:
3468 non-matching atom names
atom names from cg-junction_original.top will be used
atom names from bw.gro will be ignored


I know that there is a -maxwarn option, but I don't want to use it. When I do use it (to allow the grompp to finish), I just run into the same problem on the next step (the mdrun) when Fmax refuses to get under ~XXe03.

My configuration:

My .top file



That is the problem. Here is what I have tried/my thought process:

The easy, obvious pattern in the problem is that the atoms being compared are offset by 1; that is to say, BB2 is trying to match with BB1, but then BB3 is trying to match with BB2, and so if I constructed an atom to align those, the problem should be gone. However, that just pushes the problem elsewhere. I'm also suspicious of correcting these files by hand, since the two files are involved are bw.gro (a box of water, constructed automatically with editconf) and cg-junction.top, which is calling a file that is generated by the provided martini-dna.py (called NucleicA+Nucleic_B.itp). Since both of these files involved are automatically generated, it seems unlikely that there are simple computational errors like a missing atom (row). So tracking down the errors by hand seems both tedious and counter-productive. Furthermore, even when I do go in by hand and try to fix it, it doesn't get immediately fixed.

I'm not really sure what else could be going on here. I have cleaned my input .pdb file well, and visual inspection in VMD shows that both the original file and the coarse-grained version of it both look as they should. When I created the box of water, I did not put any -maxsol in, which should allow the box to fill as much as possible.

For the cg-junction.top file, I checked how many atoms were added in the above step, subtracted the 407 atoms generated by the Nucleic_A ... .itp files, took off another 132 atoms to make ions out of (132 since I have two intersecting arms of two 34-base strands, so ((34-1)*2)*2 = 33*4 = 132), and ~10% of the remaining beads to place as antifreeze (WF). I put those numbers into my cg-junction.top file, as they do in the tutorial's step 5. So my cg-junction_original.top file (which grompp doesn't like) looks like this:

[ system ]
; name
Martini system from j3s-cleaned.gro

[ molecules ]
; name number
Nucleic_A+Nucleic_B 1
Nucleic_C 1
Nucleic_D 1
W 23955
WF 2660
NA 132



If you've made it this far, I deeply appreciate your willingness to stick around and continue reading.

So my obvious question is how do I stop this system from blowing up? However, I'm also realizing I don't understand much of what's really going on in the martinize-dna.py file, and by extension, what the columns in its output files (Nucleic_X.itp) stand for.

If there is anything else you want to know about my files, just let me know. I put down what I thought was relevant, but as a beginner, my sense of relevance is probably pretty far off.

Any help is greatly appreciated. Many thanks.


Jonas

Please Log in or Create an account to join the conversation.

More
7 years 9 months ago #5684 by jaakko
Replied by jaakko on topic DNA Junction Blowing Up
Hi Jonas,

the warnings tell you that you coordinate file (gro) and topology file (top and itp's) do not match. You should not ignore this warning with maxwarn flag but first fix this before trying to minimize the system. You said you had cleaned up your pdb file, did you also convert it to gro file before running martinize? To figure out what's wrong you need to know what you expect to get out of martinize, i.e., what the CG representation of your system should be like: how many strands are there in your junction, should all of those be within the same elastic network etc. (Based on the top file you show I'm guessing you have 4 strands and you would like to get them in one itp file, but now you seem to have three. Using flag all-stiff for the dnatype should help with that.)

I also saw your other question about what all the sections in the itp file are. Gromacs manual has an explanation of each section and goes column by column what each parameter means: www.gromacs.org/Documentation/Manual (Chapter 5.5 in version 4.6.7 at least).

- Jaakko

Please Log in or Create an account to join the conversation.

  • jonasori
  • jonasori's Avatar Topic Author
  • Offline
  • Fresh Boarder
More
7 years 9 months ago #5690 by jonasori
Replied by jonasori on topic DNA Junction Blowing Up
Hi Jakko,

Thanks so much for the quick reply.

Yup, I did convert it over to a gro file. Also, thank you for that link to the Gromacs manual. That helped a lot. I knew that there was some simple answer like that, but just couldn't find it.

Since submitting this post, I also found that all-stiff option, which did prove to be helpful; I now have one file called Nucleic_A+Nucleic_B+Nucleic_C+Nucleic_D.itp. However, I am still having the same non-matching atom name errors.

One thing that I have noticed is that the gro file contains 828 atoms in 134 groups, whereas the itp file contains only 407 atoms in 64 groups. I looked back at my run when I just used a normal double-strand, and noticed that the corresponding itp and gro files both had 154 atoms in 24 groups. Should those numbers match up in my junction run as well? Intuitively, it seems like both should be around the ~136 mark: if two 12-base strands gave me 24 groups in the tutorial, then four 34 base strands should give me 136. And if that is the case, where could those other 72 groups be going?

This leads me to still believe that there is some going wrong in the martinize-dna code; it must somehow be reading the input gro file incorrectly or something. That input file has 4306 atoms in 136 groups, which again seems to confirm that that 136 number is about right, but I'm again not sure where to go with it after that.


Again, any further help is greatly appreciated.

~Jonas

Please Log in or Create an account to join the conversation.

More
7 years 9 months ago #5691 by jaakko
Replied by jaakko on topic DNA Junction Blowing Up
Hi Jonas,

If your structure is 4 strands of 34 bases each then the itp file should have 4*34 residues as well. When you run martinize it writes in the output how many and what kind of residues it detected in the input file. Have you checked that? That should at least tell you if the problem is that it doesn't detect the residues at all or doesn't recognize them as being nucleic. Copy paste the output here if you still have problems with it.

- Jaakko

Please Log in or Create an account to join the conversation.

  • jonasori
  • jonasori's Avatar Topic Author
  • Offline
  • Fresh Boarder
More
7 years 9 months ago #5692 by jonasori
Replied by jonasori on topic DNA Junction Blowing Up
Hi Jaakko,

Wow, that is totally it - great catch. Here is the output:


INFO Chain termini will be charged
INFO Residues at chain brakes will not be charged
INFO The elnedyn22dna forcefield will be used.
INFO Local elastic bonds will be used for extended regions.
INFO Read input structure from file.
INFO Input structure is a GRO file. Chains will be labeled consecutively.
INFO Found 12 chains:
INFO 1: A (Unknown), 28 atoms in 1 residues.
INFO 2: A (Nucleic), 1020 atoms in 16 residues.
INFO 3: A (Unknown), 34 atoms in 17 residues.
INFO 4: B (Unknown), 28 atoms in 1 residues.
INFO 5: B (Nucleic), 1002 atoms in 16 residues.
INFO 6: B (Unknown), 34 atoms in 17 residues.
INFO 7: C (Unknown), 28 atoms in 1 residues.
INFO 8: C (Nucleic), 1021 atoms in 16 residues.
INFO 9: C (Unknown), 34 atoms in 17 residues.
INFO 10: D (Unknown), 28 atoms in 1 residues.
INFO 11: D (Nucleic), 1015 atoms in 16 residues.
INFO 12: D (Unknown), 34 atoms in 17 residues.
INFO Removing HETATM chain A consisting of 1 residues.
INFO Removing HETATM chain A consisting of 17 residues.
INFO Removing HETATM chain B consisting of 1 residues.
INFO Removing HETATM chain B consisting of 17 residues.
INFO Removing HETATM chain C consisting of 1 residues.
INFO Removing HETATM chain C consisting of 17 residues.
INFO Removing HETATM chain D consisting of 1 residues.
INFO Removing HETATM chain D consisting of 17 residues.
INFO All chains will be merged in a single moleculetype.
INFO Total size of the system: 64 residues.
WARNING No secondary structure or determination method speficied. Protein chains will be set to 'COIL'.
INFO Writing coarse grained structure.
INFO (Average) Secondary structure has been determined (see head of .itp-file).
INFO Created coarsegrained topology
INFO Written 1 ITP file
INFO Output contains 1 molecules:
INFO 1-> Nucleic_A+Nucleic_B+Nucleic_C+Nucleic_D (chains A B C D)
INFO Written topology files



So clearly it is only recognizing four of the chains as being nucleic, and in those four there are 16 residues each (which gives the 64 residues I've found). When I sum up the "Unknown" residues, it does in fact come out to 136, my golden number. So that does seem to be the problem.

However, I have no idea what that means or how I could fix it. It seems like that indicates a problem in the source file with naming; would that be reasonable?

Thanks again.

Please Log in or Create an account to join the conversation.

More
7 years 9 months ago #5693 by ifaustino
Replied by ifaustino on topic DNA Junction Blowing Up
Dear Jonas,

It looks like the terminal residues of your GRO file are named as X5 or X3, which martinize-dna does not like. Can you try to rename the terminal residues as "standard" ones?

Let us know what it comes out from it.

Thanks.
Ignacio

Please Log in or Create an account to join the conversation.

More
7 years 9 months ago #5694 by jaakko
Replied by jaakko on topic DNA Junction Blowing Up
Hi Jonas,

Basically what Ignacio already said. Considering how many errors there are in your gro file it's probably best if you add a link to it, based on the output there are other suspect things in the input file in addition to the termini.

Just to reiterate the tutorial, the gro file should not have anything else than DNA and the residue names should only use DA, DC, DG and DT.

- Jaakko

Please Log in or Create an account to join the conversation.

  • jonasori
  • jonasori's Avatar Topic Author
  • Offline
  • Fresh Boarder
More
7 years 9 months ago #5696 by jonasori
Replied by jonasori on topic DNA Junction Blowing Up
Hi guys,

Thanks for those ideas. I looked at the original file and, indeed, some of the atoms have DC5 or DG3 or whatever else. I will try removing those and running it again.


This is the file that I'm using. It was created by a neighboring lab at my school, and thus is formatted pretty differently than the traditional PDB files. Still, at its core it seems the same. Have a look and let me know what you guys think.



Thank you

~Jonas

Please Log in or Create an account to join the conversation.

More
7 years 9 months ago #5698 by jaakko
Replied by jaakko on topic DNA Junction Blowing Up
Hi Jonas,

So first of all there's 250,000 (!) lines of unrelated things in the file that you should remove. Leave only the lines that relate to DNA. Second of all that's a pdb file, not a gro file so you need to convert it with editconf before running martinize. Fixing these things together with the converting the DC5/DG3 etc. residue names to two letter names should produce a proper input file that can be converted to CG topology and structure.

- Jaakko

Please Log in or Create an account to join the conversation.

  • jonasori
  • jonasori's Avatar Topic Author
  • Offline
  • Fresh Boarder
More
7 years 9 months ago - 7 years 9 months ago #5699 by jonasori
Replied by jonasori on topic DNA Junction Blowing Up
Success! Getting rid of the 3s and 5s cleaned everything up beautifully into 4 groups of 64 beads/group, allowing the grompp to work out nicely. So it seems that I now have a good (or at least error-free) coarse-grained model.

@jaakko: yes - sorry; I must have misread what you suggested I post. I posted that pdb since that was the base file that I was working on before any other modifications. I have been taking that file, cleaning out the waters, K+ and Cl- and am left with a 4307 atom file, before turning it into a gro.


However, I have now got the structure right, the energy minimization is still not working ("Steepest Descents did not converge to Fmax < 10 in 1001 steps.") I am running this minimization over 10000 steps at dt = 0.00010. This timestep was shrunk down from the suggested 0.01, since it was suggested on here that a smaller timestep could help solve this. Any thoughts?



Again, thank you guys for the help. I know that these are very basic questions and that I'm pretty reliant on you all, but hopefully as I continue to learn, that will be less the case.


EDIT: it was suggested elsewhere that an Fmax of Xe03 is on the high side, but not unreasonable since Gromacs underestimates what the reasonable max force would be. Learning this, I ran the rest of the simulation and found that it worked out pretty well. Still, it seems a little counterintuitive to me.

Best,
Jonas
Last edit: 7 years 9 months ago by jonasori. Reason: Learned new stuff

Please Log in or Create an account to join the conversation.

More
7 years 9 months ago #5701 by jaakko
Replied by jaakko on topic DNA Junction Blowing Up
Hi Jonas,

A max force of around 10e3 is fairly common with Martini, even one or two orders of magnitude larger values equilibrate fine in many cases. I wouldn't worry about the minimization Fmax since the equilibration works fine.

- Jaakko

Please Log in or Create an account to join the conversation.

Time to create page: 0.118 seconds