I didn’t use this script, but I have obtained a more reliable method:
There are three steps:
1) Calculation of the PCA using g_covar
Then calculate the inner product of the eigenvectors:
2) g_anaeig -s npt.tpr -v eigenvec_traj1.trr -eig eigenval_traj1.xvg -v2 eigenvec_traj2.trr -eig2 eigenval_traj2.xvg -over overlap_1_2.xvg -inpr inprod_1_2.xpm -first 1 -last 10
3) apply the script on the xpm file produced by g_anaeig
./xpm2normal.csh inprod_1_2.xpm
You will find more details on Gromacs (Howto) for the PCA calculations.
If you used another program you can change your trajectory to gromacs format XTC or TRR using catdcd tool.
Here is the xpm2normal.csh script:
#! /bin/csh -ef
#
# Convert an Inner product matrix generated with Gromacs
# into a matrix with numbers
# Isabella Daidone Feb 2002
# La matrice deve essere generata con g_anaeig mattendo: -v "vec[set1]" -v2 "vec[set ref.]"
# In questo modo si ottengon in ascissa i vec[set1] e in ordinata i vec[set ref.].
###################################################################################
if ($#argv < 1) then
echo "ERROR: Missing command line arguments!"
echo "USAGE: extrInPr_averSUM_sqSUM.com Mat1.xpm "
echo
exit -1
endif
cat << _eof_ | gawk -v Mat1File=$1 -f -
BEGIN {
FileComp = "InPr_averSUM_sqSUM.log"
#
# Defaults variables
#
# FileComp: name of the file for the output
#
#
printf " Matrix of Square Roots of Scalar Products Mij " > FileComp
printf "\n" > FileComp
printf " i Vector for the Tunf, j Vector for the Native" > FileComp
printf "\n" > FileComp
printf "\n" > FileComp
ini = 1
while (getline < Mat1File) {
if ( \$1 == "static" ) {
getline < Mat1File
ifi = substr(\$1,2)
# ifi = \$2
}
if (substr(\$3,1,1) == "#") {
str= \$6
gsub(/["]/,"",str)
val[substr(\$1,2)] = str
labC++
}
if (substr(\$1,1,1) == "\"" && NF == 1) {
x=1
sum=0
sumsum=0
# printf " %2d",ifi-y > FileComp
for(i=ini+1;i<=ifi+1;i++) {
k2 = substr(\$0,i,1)
map[x,y] = val[k2]
sum = sum + (map[x,y])*(map[x,y])
sum1= sum1 + (map[x,y])*(map[x,y])
printf " %10.6f",map[x,y] > FileComp
x++
}
sqsum= sqrt(sum)
printf "\n" > FileComp
y++
z=ifi-y
sum=0
}
}
sumaver=sum1/ifi
sqsumaver=sqrt(sumaver)
printf "\n" > FileComp
printf "sqsumaver: %10.6f", sqsumaver > FileComp
}
_eof_
exit -1
Yip Yew Mun (Mr) | PhD Research Scholar | Division of Chemistry & Biological Chemistry
School of Physical & Mathematical Sciences | Nanyang Technological University | Singapore 639798
Tel: (+65) 97967803 | Email: yipy0010.e.ntu.edu.sg | GMT+8h
On 18 May, 2014, at 7:31 am, Hector A. Baldoni <hbaldoni.unsl.edu.ar> wrote:
> Hello,
>
> Try the attached file, but tell us if it works !
>
> Greeting,
> Hector.
>
>
> #---Perl script for computing R.M.S.I.P. value---------------
> # each file "eig1.xvg" and "eig2.xvg" contain 10 eigenvectors
> $F1 = "eig1.xvg";
> $F2 = "eig2.xvg";
> open(F1) or die("Could not open file.");
> $N1 = 0;
> foreach $line(<F1>)
> {
> $A = $line;
> open(F2) or die("Could not open file");
> foreach $line(<F2>)
> {
> $B = $line;
> $N = ($A*$B)*($A*$B);
> print " $N \n";
> $N1 = $N1 + $N ;
> #print " $N1 \n ";
> }
> close(F2)
> }
> print "The Value of N1 is $N1 \n";
> $N1 = $N1/1;
> $N1 = sqrt ($N1);
> print "$N1";
> close(F1) ;
> #-------End of code-----------------------------------
>
>
>
>
>> Hi, I really really would appreciate any help given.
>>
>> I understand that in the archives, someone tried asking for help on a
>> script to calculate RMSIP after doing PCA using ptraj in 2008. However,
>> there was completely no response to it. While the formula for RMSIP is
>> well known in literature, I’m unable to generate a bash script to process
>> the eigenvectors from PCA to calculate RMSIP. Therefore, I’m hoping if
>> there is anyone who is willing to share his/her script on calculating
>> RMSIP.
>>
>> Your help is very very much appreciated.
>>
>> Regards
>> Yip Yew Mun
>> Graduate, PhD
>> Chemistry and Biological Chemistry
>> School of Physical and Mathematical Sciences
>> Nanyang Technological University
>>
>>
>> _______________________________________________
>> AMBER mailing list
>> AMBER.ambermd.org
>> http://lists.ambermd.org/mailman/listinfo/amber
>>
>
>
> --------------------------------------
> Dr. Hector A. Baldoni
> Area de Quimica General e Inorganica
> Universidad Nacional de San Luis
> Chacabuco 917 (D5700BWS)
> San Luis - Argentina
> hbaldoni at unsl dot edu dot ar
> Tel.:+54-(0)266-4520300 ext. 6157
> --------------------------------------
>
>
> _______________________________________________
> 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
Received on Tue May 20 2014 - 20:00:02 PDT
