VMD-L Mailing List

From: Giovanni Bellesia (gbellesia_at_chem.ucsb.edu)
Date: Wed Mar 05 2008 - 13:21:01 CST

Hi Axel,
thanks for your reply. I guess I have to give a few more details and
paste the short script I'm using.
The following is used to calculate the water oxygen-hydrogen radial
distribution fuction considering only the water
molecules within a distance of 10 angstroms from the protein.
I am interested in the number integral (third column in my output file)
since its value to 2.5 Angstroms gives an average of hydrogen bonds per
water molecule.
I ran this script for several distances and also for the entire system
(all the water molecules in my periodic box).
The value for the calculation over the entire system gives me a value
for the average number of Hbonds per molecule of 3.8 (@ 310K). The
experimental value @ 298K is 3.9 so I'm assuming that everything's fine
with my script at least in the case in which I'm considering the "whole"
water box.
The idea is that these numbers can give me a measure of how much my
small peptide is perturbating the water network.

set psf kffedimer-wb.psf
set dcd T1.dcd
mol load psf $psf dcd $dcd
set outfile1 [open gdr-OH-50.dat w]

        set sel1 [atomselect top "name OH2 and water and same residue as
within 5.0 of protein"]
        set sel2 [atomselect top "{name H1 or name H2} and water and
same residue as within 5.0 of protein"]
       
        set gr0 [measure gofr $sel1 $sel2 delta 0.1 rmax 10.0 usepbc 1
selupdate 1 first 0 last -1 step 1]
        set r [lindex $gr0 0]
        set gr [lindex $gr0 1]
        set igr [lindex $gr0 2]
        set isto [lindex $gr0 3]
        foreach j $r k $gr l $igr m $isto {
                puts $outfile1 [format "%.4f\t%.4f\t%.4f\t%.4f" $j $k $l $m]
        }
close $outfile1

exit

> On Wed, 5 Mar 2008, Giovanni Bellesia wrote:
>
> GB> Hi,
> GB> a quick question on the g(r) GUI plugin. I'm using it to calculate how the
> GB> presence of a small peptide perturbates
> GB> the water structure in its vicinity. I'm skipping the details of the
> GB> calculation here but I'm using the selection
> GB> "water and same residue as within XX from protein" with several values of XX
>
> hmmm... it is difficult to give a recommendation if you don't specify
> what you are looking for and what exactly you enter for _both_
> selections. e.g. if you enter the selection text from above for both
> entries you get something completely meaningless as you'll get the
> rdf from everything to everything. already using "water" is bad, as
> it includes the oxygens _and_ the hydrogens.
>
> GB> just to see how the perturbation decays.
>
> actually, if you are just interested into the the relative structure
> away from the protein, you don't need multiple selections. just pick
> a large value of XX and be happy.
>
> however, if you - as i suspect - want the relative structure changes
> orthogonal to the surface of the protein, you cannot use the gofr code
> at all, as it is computing a 3d isotrope distribution function.
>
> GB>
> GB> Now, I'm reading on the VMD website:
> GB>
> GB> "Note, that the normalization of g(r) has little meaning unless you enable
> GB> processing of the periodic boundary condition
> GB> and have a fixed set of atoms in the selection.
> GB> The number integrals are computed directly and thus provide accurate
> GB> coordination numbers."
>
> GB> Since I'm interested in the number integrals, I'd like to know if
> GB> the last line actually means that the number integrals are correct
> GB> even if the number of atoms varies from one frame to another
>
> yes. it is computed from directly summing up the histrograms before
> normalization. the most important step to get meaningful results is
> to be very careful in the formulation of the selection texts.
> the plugin cannot really help you there...
>
> cheers,
> axel.
>
>
> GB>
> GB> Thanks
> GB> Giovanni
> GB>
> GB>
>
>