From: James Starlight (jmsstarlight_at_gmail.com)
Date: Thu Aug 31 2017 - 05:49:39 CDT
Dear VMD/NAMD users!
I am studying the role of the internal water channels in protein function.
I am using a simple TCL script for the calculation of the radial pair
distribution function between water and some residue buried within the
active pocket of the protein - which is in contact with several water
molecules during MD.
Here is my script:
set dir_output ./output/result_water
set tittle [molinfo top get name]
set outfile1 [open $dir_output/gdr_$tittle.dat w]
set sel1 [atomselect top "resid 58 and name CA"]
set sel2 [atomselect top "water and oxygen"]
#set sel2 [atomselect top "water and oxygen within 15 of
(resid 58 and name CA)"]
set gr0 [measure gofr $sel1 $sel2 delta 0.1 rmax 10.0 usepbc 0
selupdate 1 first 0 last -1 step 1]
set r [lindex $gr0 0]
set gr [lindex $gr0 1]
set igr [lindex $gr0 2]
puts $outfile1 "@ title \"Radial Distribution Function
between between $sel1 and $sel2\""
puts $outfile1 "@ xaxis label \"r\""
puts $outfile1 "@ yaxis label \"G\(r\)\""
puts $outfile1 "@TYPE xyz"
foreach j $r k $gr l $igr {
puts $outfile1 [format "%.4f\t%.4f\t%.4f" $j $k $l]
}
close $outfile1
exit
here my questions:
1 - how the selection 2 (water) should be properly defined?
#set sel2 [atomselect top "water and oxygen"]
#or
#set sel2 [atomselect top "water and oxygen within 15 of
(resid 58 and name CA)"]
meaning that in the second case I am considering of only water within
some distance from the active site.
2 - I have tried to use just simplified version with
#set sel2 [atomselect top "water and oxygen"]
but in that case all values within the second column of the results
(the Gr) are always zero althought the values in the third column
(integral Gr) are correct. Meaning that the method is not sensitive
for my case. Here the example of the end of the results:
0.0500 0.0000 0.0000
0.1500 0.0000 0.0000
0.2500 0.0000 0.0000
0.3500 0.0000 0.0000
0.4500 0.0000 0.0000
.......
8.0500 0.0000 10.2085
8.1500 0.0000 10.4270
8.2500 0.0000 10.6279
8.3500 0.0000 10.8017
8.4500 0.0000 10.9664
8.5500 0.0000 11.1299
8.6500 0.0000 11.2811
8.7500 0.0000 11.4331
8.8500 0.0000 11.5808
8.9500 0.0000 11.7312
9.0500 0.0000 11.8822
9.1500 0.0000 12.0517
9.2500 0.0000 12.2319
9.3500 0.0000 12.4157
9.4500 0.0000 12.6212
9.5500 0.0000 12.8356
Finally if I do the same just in VMD GUI the values in the second
column are correct but they are very small like 0.0000000000002 etc.
Should i correct some script somehow to introduce normalization factor
for G(r) ?
Thanks !
James
This archive was generated by hypermail 2.1.6 : Mon Dec 31 2018 - 23:20:32 CST