Hi, Kamali,
Glad to hear leap is now behaving itself.
On Thu, Jul 5, 2012 at 10:41 AM, Kamali Sripathi <ksripath.umich.edu> wrote:
> 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.
Make sure the ions don't replace any crystal waters, if you have them...
> 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?
Excellent question that depends entirely your goals and hypotheses. Is
it your goal to figure out what the ions are doing, or do you just
want to include explicit ions to model experimental conditions? Note
that PME is pretty powerful..
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?
Again, depends on what you want to learn. 550ps is a bit on the short
side compared to our usual 1000ps equilibration, which is still on the
short side of ion equilibration times.
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.
I've read that equilibration times for ions can range from ~50 ns to
100s of ns. This is different from experimental structure
equilibration. This just means that you need 50-100s ns of straight MD
for the ion grids to converge. The Amber Masters may have more learned
recommendations to make, but I would recommend equilibrating a couple
systems to generate independent runs (using the prmtop/crd you just
built): for example your first step in equilibration is 10k steps of
minimization. The second step would be to introduce some temperature -
since this is the first step of MD, you can set ig differently
(explicitly, for each run set this number to a different value, see
Amber manual...) and therefore generate "independent trajectories".
But as a disclosure, I don't know your system. Each system is
different and you will have to read a lot of papers to figure out what
you need to do. You'll most likely have to experiment. (I'm a DNA guy,
so Cheatham, Case, Orozco, Sponer, Lavery, Beveridge and oh so many
more 8^)
"If we knew what we were looking for, it wouldn't be called
_re_search" - Einstein...
>
> 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
--
-- - -
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
Received on Thu Jul 05 2012 - 08:30:03 PDT
