Dear Amber users,
Using cpptraj, I am trying to analyze two 10ns MDs in one principal component analysis (PCA). Both MDs are closely related to each other: it is the same protein but in two different conformations (conf1 and conf2) with same number of residues. I was interested to include the structures of each conformations just after equilibration as reference structures for the PCA. I autoimaged and removed water and ions from my trajectories and reference structures and generated the corresponding stripped topologies conf1.prmtop, conf2.prmtop for each conformation. Before calculating the mass-weighted matrix, I aligned all frames to an average structure generated from both trajectories (createAverageStructure_cpptraj file).
The problem is that I get two very different results in my evecs file when switching the loading order between the reference structures and trajectories and I have no idea why.
By loading the reference structure first (PCA_cpptraj_1), I get the following 5 first eigenvalues:
0.85663
3.57670
5.85045
7.69400
8.09997
...(and so on, increasing linearly)
And by loading the MD trajectories first (PCA_cpptraj_2):
26.45372
2640022.22255
2642822.92226
2769903.48459
2844226.43589
...(and so on, continuing with values ~2e6, 3e6)
If I try not to load any reference structures, I get almost the same results as when loading the reference structure first (PCA_cpptraj_1); so I checked the reference structures and they are both fine. I checked also the structures after the RMSD on the average structure with outtraj and in both cases the structures are correctly aligned.
If necessary, I can provide a 12MB archive containing the input files and shortened MDs trajectories (20 first frames of each) which is sufficient to reproduce this behavior.
I am using cpptraj V14.09 and calling cpptraj via $AMBERHOME/bin/cpptraj -i file_cpptraj
Thanks for your help,
Best regards,
G.R.
Input Files used:
createAverageStructure_cpptraj:
------------------------------------------------------------
parm conf1.prmtop
parm conf2.prmtop?
trajin conf1_MDtraj.mdcrd parmindex 0
trajin conf2_MDtraj.mdcrd parmindex 1?
rms RMSDfirst first out RMSDfirst.dat .CA
average averageCA.pdb .CA
------------------------------------------------------------
PCA_cpptraj_1 (loading reference structures first):
------------------------------------------------------------
parm conf1.prmtop
parm conf2.prmtop
trajin conf1_afterEquil.rst7 parmindex 0
trajin conf2_afterEquil.rst7 parmindex 1
trajin conf1_MDtraj.mdcrd parmindex 0 # 500 frames
trajin conf2_MDtraj.mdcrd parmindex 1 # 500 frames
parm averageCA.prmtop
reference averageCA.pdb parmindex 2
rms RMSDaverage reference out RMSDaverage.dat .CA
matrix mwcovar name mwcvmat out mwcvmat.dat .CA
diagmatrix mwcvmat out evecs.dat vecs 200 name mymodes
run
readdata evecs.dat
projection modes evecs.dat out project.dat beg 1 end 50 .CA
------------------------------------------------------------
PCA_cpptraj_2 (loading MD first):
------------------------------------------------------------
parm conf1.prmtop
parm conf2.prmtop
trajin conf1_MDtraj.mdcrd parmindex 0 # 500 frames
trajin conf2_MDtraj.mdcrd parmindex 1 # 500 frames
trajin conf1_afterEquil.rst7 parmindex 0
trajin conf2_afterEquil.rst7 parmindex 1
parm averageCA.prmtop
reference averageCA.pdb parmindex 2
rms RMSDaverage reference out RMSDaverage.dat .CA
matrix mwcovar name mwcvmat out mwcvmat.dat .CA
diagmatrix mwcvmat out evecs.dat vecs 200 name mymodes
run
readdata evecs.dat
projection modes evecs.dat out project.dat beg 1 end 50 .CA
------------------------------------------------------------
_______________________________________________
AMBER mailing list
AMBER.ambermd.org
http://lists.ambermd.org/mailman/listinfo/amber
Received on Tue Nov 25 2014 - 10:00:02 PST
