Hi Kevin,
Thank you very much for your long and detailed response! I've been trying
out some suggestions from you and from others in this thread. I've found
that what works best for me (for future reference of others with similar
problems) is
1. Solvating my molecule
2. Using Dr. Cheatham's molarity.perl script to calculate the number of
ions needed for a certain concentration
3. Adding net neutralizing ions (Cl in my case)
4. Adding ions one molecule (of KCl in my case) at a time. This results in
non-crystallized KCl around my protein.
5. Using the "randomizeions" command in ptraj to result in a more random
distribution of the ions throughout the water box.
Because of the size of my system, I've had to write a script containing all
of the above leap commands and one containing the ptraj command, and submit
these as two separate cluster jobs. However, others with smaller systems
probably won't have to do that.
Robin's AddIonsRand command in Leap does seem faster; however, I haven't
installed the latest version of AmberTools yet, so I haven't tried it out.
I did have a question, though: is it possible to have ions that are too
well-distributed throughout the water box? I ask this mainly because I'm
concerned, with the use of the randomizeions command, that I will move the
neutralizing chloride ions too far away from the most electrostatically
favorable locations. Conversely, I know that I can specify residue numbers
in ptraj, and simply leave out the residue numbers for the neutralizing
chloride ions, but I feel that this might result in artificially
long-residence chloride ions near certain areas of the protein. Is it bad
practice to randomize ions using ptraj, and then hope that equilibration is
long enough (550 ps in my case) to allow the ions to find or return to the
most electrostatically favorable locations? I haven't found any threads in
the AMBER mailing list archive specifically addressing this question, but
please do let me know if any exist, and also if there are publications that
would be particularly helpful.
Thank you all very much, and have a great weekend,
Kamali
2012/6/22 Kevin Hauser <84hauser.gmail.com>
> Hi, Kamali,
>
>
> I haven't tried Dan's or Robin's method but will likely do so soon enough.
>
> I'm not sure why you're using addions2 or why you're adding ions before
> solvating the system. My responses, inserted below, assume all you care
> about is that you end up with 0.15 M salt and that the initial
> configuration of ions is evenly distributed around your system.
>
> On Fri, Jun 22, 2012 at 2:51 PM, Kamali Sripathi <ksripath.umich.edu>
> wrote:
>
> > Thanks a lot, Robin; I'll give that a try.
> >
> > Dr. Cheatham, I actually had a follow-up question for you. I've been
> using
> > your molarity.perl script to calculate the number of ions I'd need for a
> > certain box size. My problem is that I solvate after adding ions, so that
> > the box size is always a bit off. For example, I calculated by hand
> (before
> > I had found the script) that I needed 111 molecules of KCl to reach a
> > concentration of 0.15 M. Just to double-check, I plugged in my box
> > dimensions after adding these 111 molecules and solvated my system. The
> > script then told me that I'd need 128 molecules. I started over, and
> added
> > 128 molecules to my system and solvated. As another double-check, I
> plugged
> > in the parameters of this box into the molarity.perl script, which told
> me
> > that I now needed 131 molecules. I know the best way to do this would be
> to
> > solvate my system and use addions2 to add molecules after solvation.
> > However (as I think I mentioned), I tried this, and my system is so big
> > that Leap crashed.
>
>
> Did you try adding ions one pair at a time (as Tom suggested)? Did you try
> using addions instead of addions2 (see AmberTools15 manual page 52) or
> Robin's AddIonsRand?
>
>
> > Do you have any suggestions for how I might get my
> > desired salt concentration?
> >
>
> In case a future Google search is made that this thread comes up on:
>
> *Step 1.* Run leap with the final solvent size you want. Do not add any
> ions yet. Just solvent.
> *Step 2.* Use *tail* (see bash, please) to get the box information (last
> line of *crd file), e.g.,
> 96.1984 97.1984 98.1984 109.4712190 109.4712190 109.4712190
>
> *Step 3.* Copy and paste the box information into Professor Tom's script
> appropriately, e.g.,
> ./Tom.pl *0.15* 96.1984 97.1984 98.1984 109.4712190 109.4712190 109.4712190
>
> Resulting in,
>
> MOLARITY = 0.150
> Box size = 96.198 97.198 98.198 109.471 109.471 109.471
> Volume = 706821.110
>
> 63.847 molecules are necessary to make a molarity of 0.15M
>
> *Step 4. Bash a bit,*
>
> #!/bin/bash
> for i in $(seq 1 64); do
> echo addions mol K+ 1 Cl- 1
> done
>
> Resulting in 64 lines of,
> addions mol K+ 1 Cl- 1
>
>
> *Step 4.1* Add the 64 lines of "addions mol K+ 1 Cl-" to the leap
> script *_*after
> solvate line_**, e.g.
>
> #!/bin/bash
> #leap stuff here
> mol = loadpdb ./solute.pdb
> ...
> *solvateOCT mol TIP3PBOX 84*
> *addions mol K+ 0*
> *addions mol K+ 1 Cl- 1*
> ...
> ...
>
> *Step 4.2. Run your leap script like so,*
> ./leap_script.sh .& leap_script.out-2012-06.22.log &
> then you can
> grep "Done adding ions" leap_script.out-2012-06.22.log|wc
> to see how far along leap is, if you're impatient and in need of instant
> gratification.
>
> *Step 5. Randomize your ions. *This step ought to address your second issue
> of ions hanging around some domains and not others.
>
> ptraj prmtop << EOF
> trajin inpcrd
> trajout inpcrd.rand restart
>
> prnlev 3
>
> randomizeions :K+,Cl- around :$X by Y overlap Z noimage
> EOF
>
> Where X= the resnames of all the residues in the solute, to suit; Y and Z
> can be found in the literature. Professor Tom has some very detailed
> Methods sections 8^)
>
> Things to consider: adding net neutralizing ions; box size after
> equilibration; whether you really want ions to be randomized or you want
> them to go where the electrostatics are good. For this you may want to look
> into Dan's work. Also, Robin's method is probably not only easier but
> faster than what I've outlined here...
>
>
> HTH,
> Kevin
>
> --
> -- - -
> HK
>
>
> -----------------------------------------------------------------------------
> Kevin Hauser
> The Louis and Beatrice Laufer
> Center for Physical and Quantitative Biology
> at Stony Brook University
>
> National Institutes of Health,
> Chemical Biology Interface Training Program Fellow
>
> The Department of Chemistry
> Stony Brook University
> Stony Brook, New York 11794
>
> Phone: (561) 635.1848
> Email: 84hauser.gmail.com
>
> -----------------------------------------------------------------------------
> _______________________________________________
> AMBER mailing list
> AMBER.ambermd.org
> http://lists.ambermd.org/mailman/listinfo/amber
>
>
>
_______________________________________________
AMBER mailing list
AMBER.ambermd.org
http://lists.ambermd.org/mailman/listinfo/amber
Received on Thu Jul 05 2012 - 08:00:03 PDT
