Use MDAnalysis package of Python, It will make it much easier.

First make a atom group collecting only the water molecules within 5 angstrom of the two desired residues (they can be from anywhere in your structure file).

If you want the water to be five angstrom from both then modify the condition by taking water which are within 5 of first residue and within 5 of second residue.

Next job is to calculate their retention time within the 5 angstrom cutoff region, using a function that returns 1, while inside that cutoff else returns 0. Add this return values of this function and multiply by the timestep of the trajectory, which will give the retention time.

For making atom groups, selecting atoms within cutoff in MDAnalysis,

refers to https://www.mdanalysis.org/MDAnalysisTutorial/basics.html

I doubt that you can get meaningful results for your problem using MD.if you think that you got something meaningful after trying,please send me a message about it.

You can run such a script below in TCL in VMD:

mol new *.psf

mol addfile *.dcd

set frames [molinfo top get numframes]

set fp [open "waters.txt" w]

set ligand "resname MRV"

for {set i 0} {$i < $frames} {incr i} {

                puts "Frame: $i"

                set a [atomselect top "(water within 3 of ($ligand) and oxygen)" frame $i]

                set num [$a num]

                puts $fp "$i $num"

                $a delete

}

Dear Aniruddha Seal and Ramin Ekhteiari Salmas ..Thanks a lot for your valuable suggestion and input. I am flexible with VMD, hence, I would start with VMD.

Thanks and Regards

Budheswar

Budheswar Dehury this is the procedure for calculating the number of water molecules within 5 A of two catalytci residues belonging to two different chains of a lipid bilayer:

1-you have to create a group including both residues (you can call it "residues")

2-run this command on gromacs:

gmx select -f .trajectory.xtc -s topol.tpr -oi data.dat -select 'group "SOL" and same residue as within 0.5 of group "residues"' -n index.ndx

3-After running this command you will obtain a .dat file containing all the information you need:

1st column: time

2nd column: how many water molecules do you have contacting the group residues in that frame.

rest of the columns: the residue_id of those water molecules

so, if you want to calculate the retention period of one particular water molecule you only have to count how many times it appears in the file and divide by the number of frames...on the other hand, if you want to calculate the retention time of all of them, it's time to do a little bit of programming :) python should be a good an easy option...

hope this helps!

It's a simple problem to handle in VMD. Create a representation that shows only waters within 5A of one residue or the other, or both by using "and". Then you can save this as a pdb and count the number of water molecules in the pdb file.

For retention period, you'd need to save this for multiple frames and probably need to write a script that can count each frame and see how many of the frames each molecule is present in.

Gromacs may also have a tool like this that I've used in the past. Check out the g_dist or gmx dist tool (depending on the version of Gromacs, they're the same basically). It may be able to perform this for you IIRC from my use in the past.

http://manual.gromacs.org/current/user-guide/cmdline.html#commands-by-topic

Mario E. Valdés-Tresanco and Tim Cholko ..Thank you so much for your kind response and suggestion. I will try and give it a go.

Thanks and regards

Budheswar

Mario E. Valdés-Tresanco is it possible to write the output to pdb files show that I can make movie ..I want to create some movies for the same.

Thanking You

Budheswar

Budheswar Dehury Yes, sure! just add the flag -on and an index file will be written as output. This index file contains as many groups as frames you analyzed. Each one of these groups contains the water molecules within 5A of your residues in that frame. Then, in order to get your .pdb files just go over your trajectory through a loop and generate a water trajectory. You can easily generate a protein trajectory by using trjconv in GROMACS...

See attached a python script to do so...

then, if you want to generate a movie, you can load both files in pymol for example...

pymol water_traj.pdb protein_traj.pdb

let me know if everything goes well...

cheers!

Mario E. Valdés-Tresanco ,

Thanks a lot for your attention and valuable time in contributing a script to do the job.

However I am getting the following error code:

import-im6.q16: not authorized `os' @ error/constitute.c/WriteImage/1037.

./script_water.py: line 21: syntax error near unexpected token `"""gmx select -f md1.xtc -s md.tpr -oi data.dat -select 'group "TIP3" and same residue as within 0.5 of group "residues"' -n index.ndx -on index_water.ndx"""'

./script_water.py: line 21: `os.system("""gmx select -f md1.xtc -s md.tpr -oi data.dat -select 'group "TIP3" and same residue as within 0.5 of group "residues"' -n index.ndx -on index_water.ndx""")'

I have attached the modified script as per my requirement.

Thanks and regards

Budheswar

Budheswar Dehury try using python for executing the code:

python script_water.py

or

python3 script_water.py

cheers!

Mario E. Valdés-Tresanco

Thanks a lot for your kind suggestion and time. Now its working. I am wondering, whether do you have any idea to compute how many water molecules go through the water channel and the pore lining residues in a MD system of channels. In VMD I clearly observe how waters pass through the panel. I want to calculate the pore radius, and pore lining residue and how many water molecules pass through the channel across the lipid bilayer over the course of time.

Do you have any suggestion or input.

Looking forward to hear from you.

Thanks and regards

Budheswar

Budheswar Dehury yes, I have an idea how to do that, and I would like to help you, just I am super busy right now... I'll get back to you as soon as I have something...anyway if you find a script or program to do that just let me know...

cheers!

Hi,

What you are trying to accomplish is a part of the Desmond MD module from the Schrodinger LLC software suite. It is not how many water molecules are there but if they are stably bound in the place. Look up my paper in JMGM from 2019 and you will see that the solution is delivered to you on a plate.

Best,

Wes