From: Andrew Dalke (
Date: Wed Dec 09 2009 - 17:20:43 CST

On Dec 9, 2009, at 11:30 PM, Bernardo Sosa Padilla Araujo wrote:
> I was wondering if any of you had a script or at least any idea about how to write one, that can count the number of water molecules around a specific part of an enzyme.
> Any suggestions would be greatly appreciated.

In general you can do this on the interactive window using the 'atomselect' command

For example,

  set sel [atomselect top "water within 5 of resid 28"]

will find all the atoms in all the water molecules within 5A of residue 28 and put the result into the variable "sel", which is the variable name I most often use because it's short for "selection'.

You can then get information about how many atoms are in the selection like this:

  $sel num

However, bear in mind that the selection includes hydrogens, if they are present in your waters. (Most structure files omit hydrogens, but not all.) If that happens then you'll want to restrict the selection to just the oxygens

  set sel [atomselect top "element O and water within 5 of resid 28"]

Now, VMD knows that a residue is a water if the residue name is H2O, HHO, OHH, or HOH. This may not be true for all input structures. If you have something else then you'll have to know more about how your file defines what a water is.

Also, your definition of "specific part" might not be a residue. You can select a single atom, or multiple residues, or anything in between by using a different atom selection than "resid 28". Or, your active site might be a region, like inside a bounded box, or within a certain distance of a point in space. For those you can construct queries based on the x, y, z coordinates of each atom.

The easiest way to test out a selection is with the selection window. Once you've figured it out then you can use it in your script.


(Felt like answering this one for old time's sake. I wrote the selection language in VMD a long time ago and I'm still rather proud about it.)