Re: [AMBER] Different PCA results when loading trajectories in a different order

From: Daniel Roe <daniel.r.roe.gmail.com>
Date: Tue, 25 Nov 2014 11:08:33 -0700

PS - Also, if this is the same protein but different conformations you
should only need 1 topology (in fact if this is the case then
conf1.prmtop and conf2.prmtop should only differ by their VERSION
line). The only time you need a different topology file for two
trajectories is if the atoms/residues do not match between the two
(such as your PDB that only contains CA atoms).

On Tue, Nov 25, 2014 at 11:03 AM, Daniel Roe <daniel.r.roe.gmail.com> wrote:
> Hi,
>
> Out of curiosity do you see the same thing happen if you use
> non-mass-weighted covariance matrices (i.e. covar instead of mwcovar)?
> Note that typically in PCA just plain covariance matrices are used.
> Mass-weighted covariance matrices are really only required if you want
> to perform quasi-harmonic analysis.
>
> -Dan
>
> On Tue, Nov 25, 2014 at 10:49 AM, Guillaume Roellinger
> <roellinger.bio.mx> wrote:
>> 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
>
>
>
> --
> -------------------------
> Daniel R. Roe, PhD
> Department of Medicinal Chemistry
> University of Utah
> 30 South 2000 East, Room 307
> Salt Lake City, UT 84112-5820
> http://home.chpc.utah.edu/~cheatham/
> (801) 587-9652
> (801) 585-6208 (Fax)



-- 
-------------------------
Daniel R. Roe, PhD
Department of Medicinal Chemistry
University of Utah
30 South 2000 East, Room 307
Salt Lake City, UT 84112-5820
http://home.chpc.utah.edu/~cheatham/
(801) 587-9652
(801) 585-6208 (Fax)
_______________________________________________
AMBER mailing list
AMBER.ambermd.org
http://lists.ambermd.org/mailman/listinfo/amber
Received on Tue Nov 25 2014 - 10:30:03 PST
Custom Search