Converting Bond-Angle-Torsion to Anchored Cartesian coordinates
This tutorial describes the use of a set of MATLAB scripts for converting given BAT coordinates of a molecule to anchored cartesian coordinates. See Figure 1 of J. Chem. Phys. 127, (024107), 2007 for a description of BAT coordinates.
Contents
Overview
The various steps of the method are demonstrated using all-atom Cyclohexane as a test system. It contains 18 atoms - 6 carbons and 12 hydrogens which correspond to 3*18-6 = 48 internal BAT coordinates - 17 bonds, 16 angles and 15 torsions.
Click here to download a zip file containing files required for this tutorial.
In following sections, file names and matlab variables are in this font.
Connectivity Matrix
The set of BAT coordinates used to define a conformation are specified by the connectivity matrices. load cyclohexaneConn.mat to import the connectivity data structure.
The connectivity is defined in terms of the atom numbers in the pdb file (see cyclohexaneOrig.pdb). Thus atom 1 is c1, atom 2 is c2 .. atom18 is h18.
Bond connectivity is contained in matrix bondconn which has 17 rows - one for each bond. The first two numbers in each row are the atom numbers corresponding to the bond. The third column is irrelevant for this tutorial. The fourth column gives the corresponding column number in matrix bonds which contains values of this bond in different conformations.
In addition to bondconn, function genConfBATConn requires another data structure b which contains information about the columns of the bondconn matrix. For this tutorial only b.SampleIdCol need to be set which here is 4 indicating the 4th column of bondconn matrix
Angle connectivity matrix angleconn and data structure a are analogous to bonds. Now first three columns of angleconn contain the atom numbers forming the angle. Again, only cols 1:3 and 5 of angleconn and a.SampleIdCol are relevant here.
Torsions are defined by torsionconn and t. Cols [1:4 6] of torsionconn and t.SampleIdCol required.
The connectivity matrices here were constructed 'by hand'.
BAT input data
load BATinput.mat to import BAT coordinates of 10 conformations in the pdb trajectory file cyclohexaneOrig.pdb. bonds is 10x17 matrix with each row containing the 17 bond lengths of corresponding pdb structure in cyclohexaneOrig.pdb Similarly, angles and torsions matrices contain the angle and torsion values in radians.
Conversion command
[conform] = genConfBATConn(bonds,angles,torsions,bondconn,angleconn,torsionconn,b,a,t);
conform contains the converted cartesian coordinateis for each atom. Each row contains one conformation with entries as
For the current case it contains 10 rows - one for each conformation.
Create pdb trajectory of reconstructed frames
The cartesian coordinates for the 10 conformations in conform can be used to write a pdb trajectory file as described here.
load atomNameCyclohexane.mat for atom name, and load molNameCyclohexane.mat for molecule name to use in pdb file.
coord2pdbWname(conform,'cyclohexaneRecons',-1,atomName,molName) then creates the pdb trajectory file reconsCyclohexane.pdb.
Note that due to way BAT coordinate system is set up, for each conformation atom h7 is placed on the origin, atom c1 is on the x-axis and atom c2 is in the x-y plane.
In other words, we from all have removed six rigid body translational and rotational degrees of freedom.
VMD check
To verify if the conversion was done correctly compare the original trajectory cyclohexaneOrig.pdb with reconstructed cyclohexaneRecons.pdb one in VMD 1.8.6
To superimpose the two trajectories:
1. goto VMD Main -> Extensions -> Analysis -> RMSD Trajectory Tool
2. replace 'protein' by 'all'
3. click 'Add All'
4. click 'Align'
If everything went well, corresponding frames of the two trajectories should overlap exactly!