PS - For a concrete example of performing PCA with cpptraj see
http://pubs.acs.org/doi/abs/10.1021/jp4125099 and supporting material
therein.
-Dan
On Tue, Jul 29, 2014 at 8:16 AM, Daniel Roe <daniel.r.roe.gmail.com> wrote:
> Hi,
>
> You're missing some elements from your script, and have a few things out
> of order. By default in cpptraj when you specify an action or an analysis
> it goes into either the action list or analysis list respectively; it is
> not immediately executed. I will go through your input and try to explain
> what's happening:
>
> On Fri, Jul 25, 2014 at 11:56 AM, Nicole Ippolito <
> ippolito.nicole89.gmail.com> wrote:
>
>> trajin ~/BRAF/MD/V600E505H/sor/V600E505H_sor_aMD_0_300ns.crd.gz
>> createcrd crd1
>> run
>>
>
> Here you are essentially loading a trajectory into the coords data set
> named 'crd1'. The 'run' statement means "process all currently loaded
> trajectories with any initialized actions", and 'createcrd' is the only
> action set up here.
>
>
>> rms first .CA
>> reference V600E505H_sor_aMD_avg_pdb.pdb
>>
>
> These statements seem out of order. They are placed after the run
> statement, and so they are not processed until the first 'run' is
> completed. Presumably your intent was to remove overall rotation and
> translation from your coordinates with the 'rms' statement. If that is the
> case then it should come before the createcrd statement. Also, since you
> also have a later 'rms' command that fits to a reference, this makes me
> think this first 'rms' was for creating an average structure to use as a
> reference, and the 'average' command is missing.
>
> The 'reference' command is actually what kills the run: if you check your
> output you will see the message "Error: reference: Could not set up
> trajectory". This indicates that "V600E505H_sor_aMD_avg_pdb.pdb" does not
> exist; perhaps you meant to create an average?
>
>
>> crdaction crd1 rms reference .CA
>> crdaction crd1 matrix out matrix.dat name mymatrix covar .CA
>>
> runanalysis diagmatrix mymatrix out evecs.dat vecs 20 name mymodes
>>
>
> These statements seem fine as long as a reference is successfully loaded
> and your 'crd1' coordinates is properly loaded.
>
>
>> crdaction crd1 projection modes evecs.dat start 1 stop 3 .CA myproj
>>
>
> This will not work, since you named your modes internally as 'mymodes' but
> you're trying to use the file name 'evecs.dat'. Change 'evecs.dat' to
> 'mymodes'.
>
>
>> hist myproj:1 bins 200 out PC.hist.agr
>> hist myproj:2 bins 200 out PC.hist.agr
>>
>
> This seems like it should work as well.
>
> So to summarize, your final script should look something like this:
>
> trajin ~/BRAF/MD/V600E505H/sor/V600E505H_sor_aMD_0_300ns.crd.gz
> # Fit to first frame for averaging
> rms first .CA
> # Average the trajectory
> average V600E505H_sor_aMD_avg_pdb.pdb pdb
> # Save the rms-fit coords
> createcrd crd1
> # First run
> run
> # Load the generated average as reference
> reference V600E505H_sor_aMD_avg_pdb.pdb
> # Fit saved coords to the reference
> crdaction crd1 rms reference .CA
> # Calculate covariance matrix from fit coords
> crdaction crd1 matrix out matrix.dat name mymatrix covar .CA
> # Diagonalize matrix to get eigenvectors/eigenvalues
> runanalysis diagmatrix mymatrix out evecs.dat vecs 20 name mymodes
> # Project fit coords along eigenvectors 1-3
> crdaction crd1 projection modes evecs.dat start 1 stop 3 .CA myproj
> # Set up histograms for PC projections
> hist myproj:1 bins 200 out PC.hist.agr
> hist myproj:2 bins 200 out PC.hist.agr
> hist myproj:3 bins 200 out PC.hist.agr
> # Run analyses
> runanalysis
>
> Note that the final 'runanalysis' is not strictly necessary as cpptraj
> will try to run any actions or analyses that are set up before quitting
> when reading a script, but I include it for completeness.
>
> Hope this helps,
>
> -Dan
>
>
>
>
>>
>>
>> CPPTRAJ: Trajectory Analysis. V14.01
>> ___ ___ ___ ___
>> | \/ | \/ | \/ |
>> _|_/\_|_/\_|_/\_|_
>> Reading 'parm' as Amber Topology
>> INPUT: Reading Input from STDIN
>> [trajin ~/BRAF/MD/V600E505H/sor/V600E505H_sor_aMD_0_300ns.crd.gz ]
>> Reading
>> '/home/aubndi/BRAF/MD/V600E505H/sor/V600E505H_sor_aMD_0_300ns.crd.gz' as
>> Amber Trajectory
>> [createcrd crd1 ]
>> CREATECRD: Saving coordinates from Top parm to "crd1"
>> [run]
>> ---------- RUN BEGIN -------------------------------------------------
>>
>> PARAMETER FILES:
>> 0: 'parm', 36610 atoms, 11002 res, box: Orthogonal, 10731 mol, 10722
>> solvent, 60000 frames
>>
>> INPUT TRAJECTORIES:
>> 0: 'V600E505H_sor_aMD_0_300ns.crd.gz' is an AMBER trajectory, Parm parm
>> (Orthogonal box) (reading 60000 of 60000)
>> Coordinate processing will occur on 60000 frames.
>> TIME: Run Initialization took 0.0000 seconds.
>>
>> BEGIN TRAJECTORY PROCESSING:
>> .....................................................
>> ACTION SETUP FOR PARM 'parm' (1 actions):
>> 0: [createcrd crd1 ]
>> ----- V600E505H_sor_aMD_0_300ns.crd.gz (1-60000, 1) -----
>> 0% 10% 20% 30% 40% 50% 60% 70% 80% 90% Error: reference: Could not set up
>> trajectory.
>> 1 errors encountered reading input.
>> 100% Complete.
>>
>> Read 60000 frames and processed 60000 frames.
>> TIME: Trajectory processing: 1743.3814 s
>> TIME: Avg. throughput= 34.4159 frames / second.
>>
>> ACTION OUTPUT:
>>
>> DATASETS:
>> 1 data set:
>> crd1 "crd1" (coordinates), size is 60000 (25138 MB) Box Coords,
>> 36610 atoms
>> ---------- RUN END ---------------------------------------------------
>> [rms first .CA]
>> RMSD: (.CA), reference is first frame (.CA), with fitting.
>> [reference V600E505H_sor_aMD_avg_pdb.pdb]
>> TIME: Total execution time: 1743.8619 seconds.
>> _______________________________________________
>> 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 201
> 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 201
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 Jul 29 2014 - 08:00:02 PDT
