Contents

DYNAMO Molecular Dynamics and Analysis System

INTRODUCTION

DYNAMO is a system of software tools and scripts for calculating and evaluating molecular structures. DYNAMO includes a cartesian-coordinate simulated annealing engine, and facilities for NMR homology search to assemble collections of molecular fragments which are consistent with NMR observables. The tools of DYNAMO are accessed via scripts written in the TCL/TK scripting language.

DYNAMO can use the following kinds of experimental restraints during a simulated annealing protocol:

Some Examples of DYNAMO Applications include:

Molecular Dynamics with DYNAMO

The DYNAMO molecular dynamics engine is used to generate structures of proteins and nucleic acids which are consistent with both standard covalent geometry and also any user-supplied NMR-derived structural constraints such as interproton distances. As such, DYNAMO does not perform a realistic physical simulation of molecular dynamics; rather it is used as an optimization tool to generate structures for best agreement with measured data. Structure determination from NMR data is particularly challenging because the number of structural constraints available from NMR spectroscopy is small compared to the number of degrees of freedom in a biological macromolecule.

DYNAMO is able to use a variety of NMR-derived constraints: NOE distances, torsion angles, J couplings (in the form of Karplus expressions), and dipolar couplings. In addition, protein backbone chemical shifts can be translated into backbone torsion angle constraints through the use of the program TALOS.

In addition to experimental constraints, DYNAMO also incorporates structural restraints that are known on the basis of the covalent geometry of the molecule. These are expressed in terms of Van der Waals radius parameters for the atoms, expected bond distances and angles, and torsions which describe any planarity or stereochemistry in the molecular system.

DYNAMO uses molecular dynamics based simulated annealing to find sets of atomic coordinates that are consistent with the experimental and covalent constraints. Molecular dynamics is a computational technique that uses numerical integration to simulate Newton's equations of motion. Simulated annealing is a computational technique for function optimization. In simulated annealing, a system's temperature begins at a high value and is slowly reduced. Thus, early in a simulated annealing procedure, the system has sufficient kinetic energy to jump over large potential energy barriers. But as the system cools, its movement is gradually restricted to lower and lower energy structures. A sufficiently high starting temperature ensures that all the relevant conformational space is sampled, and a sufficiently slow cooling rate ensures that a very low energy structure is finally chosen.

The DYNAMO energy function

Molecular dynamics uses numerical integration to simulate a system's movement by Newton's equations of motion. The forces on each atom are the partial derivatives of the potential energy function. DYNAMO's potential energy function includes terms that describe the current coordinates' agreement with both experimental and a priori constraints:

   Etot = Ebond + Eangle + Eimpr + Evdw + Enoe + Ejcoup + Etorsion + Edipo + Epcs Ergyr ...

where

Etot is the total energy.

Ebond describes the structure's agreement with expected covalent bond lengths:

   Ebond = kbond Σ (expected_distance - actual_distance)2

Eangle describes the structure's agreement with expected covalent bond angles:

   Eangle = kangle Σ (expected_angle - actual_angle)2

Eimpr describes the structure's agreement with expected "improper" torsion angles. Improper torsion angles are constraints on the geometry of groups of four atoms that are used to maintain correct chirality of chiral centers and the planarity of planar groups:

   Eimpr = kimpr Σ (expected_impr - actual_impr)2

Evdw describes the structure's steric overlap, using a simple repulsive potential:

   Evdw = kvdw Σ ((size_scale * expected_min_distance)2 - actual_distance2)2

Enoe describes the structure's agreement with observed NOE distances. It uses a "square well" potential that is zero if the actual distance is within a given range:

   Enoe = knoe Σ delta_dist2

        if actual_dist > observed_dist_upperbound:
           delta_dist = actual_dist - observed_dist_upperbound

        if actual_dist < observed_dist_lowerbound:
           delta_dist = actual_dist - observed_dist_lowerbound

        otherwise:
           delta_dist = 0

Ejcoup describes the structure's agreement with observed J-couplings. Karplus parameters can be specified for each type of J-coupling.

   Ejcoup = kjcoup Σ (observed_J - calculated_J)2

      calculated_J = A cos2 (phi + phase) + B cos (phi + phase) + C

Etorsion describes the structure's agreement with experimental torsion angle constraints. Like the NOE term, it uses a square well potential:

   Etorsion = ktorsion Σ delta_torsion2

        if actual_torsion > observed_torsion_upperbound:
           delta_torsion = actual_torsion - observed_torsion_upperbound 
   
        if actual_torsion < observed_torsion_lowerbound: 
           delta_torsion = actual_torsion - observed_torsion_lowerbound
   
        otherwise:   
           delta_torsion = 0

Edipo describes the structure's agreement with experimentally observed dipolar couplings. Data from one or more different alignment tensors can be used simultaneously. The expected dipolar couplings are calculated from the structure as described the section on Options for Dipolar Coupling Data below.

   Epca = kpcs Σ (observed_pcs - calculated_pcs)2

Epcs describes the structure's agreement with experimentally observed pseudo-contact shifts between a given atom and a paramagnetic center. Data from one or more paramagnetic centers can be used simultaneously.

   Edipo = kdipo Σ (observed_dipo - calculated_dipo)2

Ergyr describes the structure's agreement with its expected overall size. NMR structures are typically expanded relative to crystal structures, because there are so many forces from the vdw term pushing atoms apart, with relatively few forces from the NOE term pulling them back together. The radius of gyration is a measure of a structure's overall size--specifically, it is the RMS distance of a group of atoms to their centroid. It can be predicted reasonably accurately for globular proteins just from the number of residues:

   expected_Rgyr = 2.2 Nresidues0.38

      Ergyr = krgyr Σ (expected_Rgyr - actual_Rgyr)2
Steps in a DYNAMO Structure Calculation

DYNAMO organizes a structure calculation project in a "GMC" (Generalized Molecular Coordinate) directory. This directory starts with a collection of tables which describe the covalent geometry of the molecular system. These tables are created by DYNAMO based on sequence information. After the GMC directory has been created, additional tables can be added containing the experimental restraints. The typical steps are as follows:

  1. Specify the residue sequences of the segments in the molecule or complex, which should be saved in the file gmc_record.txt. This can be done in a few ways:


  2. Based on the sequence information and options in the file gmc_record.txt, DYNAMO automatically creates the definitions of covalent geometry. In DYNAMO, these are stored as a series of NMRPipe-format tables with specific names in a General Molecular Coordinate (GMC) directory:

    File Name Contents
    atoms.tab Specification of All Atoms
    bonds.tab Specification of Bonds
    angles.tab Covalent Geometry Angles
    impropers.tab Covalent Geometry Torsions
    vdwex.tab Van der Waals Exclusions

    During creation of the GMC, DYNAMO also creates an initial structure in the GMC directory, called "random.pdb" with random (!) atomic coordinates.



  3. After GMC creation, one or more user-supplied experimental restraint tables, also with specific names, are added to the GMC directory. As a convenience, some example TCL scripts are provided to help prepare these tables from existing restraints in XPLOR format (Note: these scripts are only examples, they may not work with every variation in restraint format, and they may require that input atom names match the DYNAMO ones exacty).

    The DYNAMO user-supplied restraint tables can include:

    File Name Contents
    noes.tab NOE Distances
    torsions.tab Torsion Restraints
    radgyr.tab Radius of Gyration Restraints
    jcoup.tab J Couplings
    ac.tab Atom Coordinate Restraints
    ds.tab Distance Symmetries
    pcsObs.tab Pseudo-Contact Shifts
    dObsA.tab Dipolar Couplings (A, B, C etc.)

  4. A DYNAMO simulated annealing script is used to read the GMC directory contents and execute one or more simulated annealing runs.


  5. Subsequent DYNAMO scripts are visualize and evaluate the structure calculation results.

  6. Additional structure calculations may be performed, using previous results as a starting point, or with additional or adjusted restraints.

DYNAMO Restraint Table Format

The constraint files within a GMC directory are all variations on the nmrPipe table format. Generally, each line will specify a restrain involving one or more atoms. A given atom is specified according to four parameters: atom name (ATOMNAME) residue ID (RESIDRESNAME) and segment name (SEGNAME). In the case of restraints which involve two or more atoms, the atoms within a restraint are referred to as atom I, atom J, etc. For example, in the table of covalent bond length restraints ("bonds.tab") each entry describes a bond distance between two atoms I and J, with an expected bond distance (D) and a force constant for the bond (FC):

VARS SEGNAME_I RESNAME_I RESID_I ATOMNAME_I SEGNAME_J RESNAME_J RESID_J ATOMNAME_J D FC
FORMAT %8s %3s %4d %4s %8s %3s %4d %4s %8.3f %8.3f

    UBIQ MET    1    C     UBIQ MET    1    O    1.231 1000.000
    UBIQ MET    1   CA     UBIQ MET    1    C    1.525 1000.000
    UBIQ MET    1   CA     UBIQ MET    1   CB    1.530 1000.000
    UBIQ MET    1   CA     UBIQ MET    1   HA    1.080 1000.000
    UBIQ MET    1   CB     UBIQ MET    1   CG    1.530 1000.000

The first two lines identify the contents of each data column. The two atoms involved are identified by their segment name, residue name, residue number, and atom name. Other information relevant to the constraint (the equilibrium bond length in Ångstroms and the force constant in kcal/mol/Å2) are defined in the last two columns.

Similar formats are used to define the other covalent geometry restraints in atoms.tab, angles.tab, and impropers.tab. Since these are generated automatically by the GMC Editor, they rarely need to be examined or edited.

The experimental restraints tables follow a similar format. For example, each line in a torsion angle restraint table contains specifications for the four atoms decsribing the torsion, I, J, K and L. And, a given torsion restraint is specified in terms of allowed lower and upper bounds in degrees (ANGLE_LO and ANGLE_HI). Each restraint also includes a unique INDEX number to identify the restraint:

VARS   INDEX  SEGNAME_I RESID_I RESNAME_I ATOMNAME_I \
              SEGNAME_J RESID_J RESNAME_J ATOMNAME_J \
              SEGNAME_K RESID_K RESNAME_K ATOMNAME_K \
              SEGNAME_L RESID_L RESNAME_L ATOMNAME_L \
              ANGLE_LO  ANGLE_HI FC

FORMAT %4d %4d %4s %4s %4s %4d %4s %4s %4s %4d %4s %4s %4s %4d %4s %4s %4s %8.3f %8.3f %8.3f

   1   UBIQ  1 MET C    UBIQ   2 GLN N    UBIQ   2 GLN CA   UBIQ   2 GLN C   -138.00  -66.00    1.0
   2   UBIQ  2 GLN C    UBIQ   3 ILE N    UBIQ   3 ILE CA   UBIQ   3 ILE C   -173.00 -109.00    1.0
   3   UBIQ  3 ILE C    UBIQ   4 PHE N    UBIQ   4 PHE CA   UBIQ   4 PHE C   -149.00  -93.00    1.0
   4   UBIQ  4 PHE C    UBIQ   5 VAL N    UBIQ   5 VAL CA   UBIQ   5 VAL C   -163.00  -75.00    1.0
   5   UBIQ  5 VAL C    UBIQ   6 LYS N    UBIQ   6 LYS CA   UBIQ   6 LYS C   -133.00  -85.00    1.0
   6   UBIQ  6 LYS C    UBIQ   7 THR N    UBIQ   7 THR CA   UBIQ   7 THR C   -142.00  -58.00    1.0
   7   UBIQ  7 THR C    UBIQ   8 LEU N    UBIQ   8 LEU CA   UBIQ   8 LEU C    -82.00  -42.00    1.0

Sixteen columns are needed to specify the segment name, residue number, residue name, and atom name of each of the four atoms that define the torsion angle. Finally, two columns give the lower and upper bounds on the angle (in degrees), followed by a force constant (in kcal/mol/rad2).

The J-couplings restraints table has a format similar to the torsion angle restraints, but in addition, each restraint includes karplus parameters:

VARS INDEX SEGNAME_I RESID_I RESNAME_I ATOMNAME_I \
           SEGNAME_J RESID_J RESNAME_J ATOMNAME_J \
           SEGNAME_K RESID_K RESNAME_K ATOMNAME_K \
           SEGNAME_L RESID_L RESNAME_L ATOMNAME_L A B C PHASE OBSJ FC 

FORMAT %3d %5d %4s %4s %4s %5d %4s %4s %4s %5d %4s %4s %4s %5d %4s %4s %4s %9.5f %9.5f %9.5f %9.5f %9.3f %.2f 

   1    UBIQ  1 MET  C    UBIQ   2 GLN  N    UBIQ   2 GLN  CA   UBIQ   2 GLN  C     1.74000  -0.57000   0.25000  240.00000    2.30 1.0 
   2    UBIQ  3 ILE  C    UBIQ   4 PHE  N    UBIQ   4 PHE  CA   UBIQ   4 PHE  C     1.74000  -0.57000   0.25000  240.00000    1.30 1.0 
   3    UBIQ  4 PHE  C    UBIQ   5 VAL  N    UBIQ   5 VAL  CA   UBIQ   5 VAL  C     1.74000  -0.57000   0.25000  240.00000    0.80 1.0 
   4    UBIQ  5 VAL  C    UBIQ   6 LYS  N    UBIQ   6 LYS  CA   UBIQ   6 LYS  C     1.74000  -0.57000   0.25000  240.00000    1.60 1.0 

Once again, the first column defines an "index" number for each constraint. The next sixteen columns select the four atoms involved in the torsion angle. Three more columns define the Karplus parameters A, B, and C for that J-coupling (where Jcalc = A cos2 theta + B cos theta + C). The next column allows a phase correction (in degrees) to be applied to the Karplus calculation. This phase correction is necessary in cases where the atoms selected to define the torsion angle (in the above example, C..N..Ca..C atoms are used to define a backbone phi angle) are different from the atoms involved in creating the observed coupling (in this case, HN..Ha). Finally, the observed coupling constant (in Hz) is defined, together with a force constant (in kcal/mol/Hz2).

The format of the dipolar couplings restraints file is:

VARS   SEGNAME_I RESID_I RESNAME_I ATOMNAME_I SEGNAME_J RESID_J RESNAME_J ATOMNAME_J D      DD    W
FORMAT %5d     %6s       %6s        %5d     %6s       %6s    %9.3f   %9.3f %.2f

UBIQ    1    MET      C    UBIQ   2    GLN     HN       3.993     0.333 3.00
UBIQ    2    GLN      C    UBIQ   3    ILE     HN      -5.646     0.333 3.00
UBIQ    3    ILE      C    UBIQ   4    PHE     HN       1.041     0.333 3.00
UBIQ    4    PHE      C    UBIQ   5    VAL     HN       0.835     0.333 3.00

The first eight columns select the two atoms involved in the dipolar coupling. The next column gives the observed value of the coupling (in Hz). The next column gives an indication of the possible range of values for this coupling, relative to the HN..N coupling. The final column defines a force constant (in kcal/mol/Hz2).

The NOE table format is somewhat more complex than the other experimental restraints tables:

VARS INDEX GROUP RESID_I RESNAME_I ATOMNAME_I SEGNAME_I RESID_J RESNAME_J ATOMNAME_J SEGNAME_J D_LO D_HI FC W S
FORMAT %4d %3d  %5d %6s %6s %4s  %5d %6s %6s %4s  %9.3f %9.3f %.2f %.2f %.2f

   1  1    UBIQ MET  1  HA      UBIQ MET  1  HG1     1.800   4.459 1.0 1.0 1.0
   1  2    UBIQ MET  1  HA      UBIQ MET  1  HG2     1.800   4.459 1.0 1.0 1.0
   2  1    UBIQ MET  1  HA      UBIQ MET  1  HG1     1.800   4.060 1.0 1.0 1.0
   2  2    UBIQ MET  1  HA      UBIQ MET  1  HG2     1.800   4.060 1.0 1.0 1.0
   3  1    UBIQ MET  1  HB2     UBIQ MET  1  HE1     1.800   3.565 1.0 1.0 1.0
   3  1    UBIQ MET  1  HB2     UBIQ MET  1  HE2     1.800   3.565 1.0 1.0 1.0
   3  1    UBIQ MET  1  HB2     UBIQ MET  1  HE3     1.800   3.565 1.0 1.0 1.0
   4  1    UBIQ MET  1  HE1     UBIQ MET  1  HG1     1.800   3.155 1.0 1.0 1.0
   4  1    UBIQ MET  1  HE2     UBIQ MET  1  HG1     1.800   3.155 1.0 1.0 1.0
   4  1    UBIQ MET  1  HE3     UBIQ MET  1  HG1     1.800   3.155 1.0 1.0 1.0
   4  2    UBIQ MET  1  HE1     UBIQ MET  1  HG2     1.800   3.155 1.0 1.0 1.0
   4  2    UBIQ MET  1  HE2     UBIQ MET  1  HG2     1.800   3.155 1.0 1.0 1.0
   4  2    UBIQ MET  1  HE3     UBIQ MET  1  HG2     1.800   3.155 1.0 1.0 1.0

As in the J-coupling restraint file, the first column is an index number to identify a particular constraint. Columns 3-10 identify two atoms whose distance will be used in calculating the NOE energy. The next two columns define that restraint's lower and upper distance bounds (in Ångstroms), and the last column sets a force constant for that restraint (in kcal/mol/Å2).

But there are two differences in the NOE restraints file's format. The first is that more than one line of the file can be given the same index number. Having more than one line with the same index number tells dynamo that these lines are all part of a single NOE constraint. The second difference in the NOE restraints file format is the addition of the GROUP column. Like the INDEX column, the group column contains an integer. It is used in combination with the index column to define different possible assignments for a particular NOE restraint. Lines in the NOE restraints file which have the same index number and the same group number are treated as part of a single possible assignment for an NOE. The interatomic distances from each line's atoms are averaged together using an R-6 sum average. Lines (or groups of lines) with the same index number but different group numbers are treated as separate possible assignments for a given NOE.

The DYNAMO Annealing Stages

      Keyword    Meaning

      init       Initialization
      high       High-Temperature
      coolStart  Cooling
      coolEnd    Done
      cool       Sets parameters for both "coolStart" and "coolEnd"
      all        Sets parameters for all stages to one given value.
There are three parameter classes:
      -sa   Annealing Schedule Parameters
      -fc   Force Constant Scaling
      -size Size or Radius Scaling
Parameters are set via triplets of "parameterName stageName value", for example:
dynSimulateAnnealing -graph -print 10 -nocenter \
    -sa    stepCount   init       100 \
           stepCount   high       100 \
           stepCount   cool      3000 \
           temperature init       500 \
           temperature high       500 \
           temperature coolStart  500 \
           temperature coolEnd      0 \
    -fc    ac          init       8.0 \
           ac          high       8.0 \

Simulated Annealing Parameter Names and Default Values

      Parameter Name     Type Scaling init   high   coolStart coolEnd

      stepCount          sa   none       500   2000   12000        0
      temperature        sa   other   4000.0 4000.0  4000.0        0
      temperatureStep    sa   none       0.0    0.0    25.0     25.0
      temperatureControl sa   none       1.0   10.0    10.0     10.0
      timeStep           sa   none       3.0    5.0     5.0      5.0

      bond               fc   power      1.0    1.0     1.0      1.0
      dist               fc   power      1.0    1.0     1.0      1.0
      angle              fc   power      0.5    0.5     0.5      1.0
      improper           fc   power      0.1    0.1     0.5      1.0
      torsion            fc   power     10.0   10.0    10.0    200.0
      j                  fc   power      0.0    0.0     0.1      1.0
      dt                 fc   power      0.0    0.0     0.0      0.0
      ds                 fc   power      0.0    0.0     0.0      0.0
      noe                fc   power      2.0    2.0     2.0     30.0
      ac                 fc   power      0.0    0.1     0.1      3.0
      dc                 fc   power      0.0    0.0     0.005    0.5
      dcA                fc   power      0.0    0.0     0.005    0.5
      dcB                fc   power      0.0    0.0     0.005    0.5
      dcC                fc   power      0.0    0.0     0.005    0.5
      pcs                fc   power      0.0    0.0     1.0      1.0
      cs                 fc   power      0.0    0.0     1.0      1.0
      radGyr             fc   power      0.0  100.0    10.0    100.0
      vdw                fc   power      0.0    0.0     0.004    4.0

      radGyr             size linear     1.0    1.3     1.3     0.90
      vdw                size power      0.0    0.0     0.9     0.81

DYNAMO Arguments for Reading GMC Directory

dynReadGMC
   -gmc gmcDirName    Name of DYNAMO Input GMC Directory.
   -pdb pdbFileName   Name of Initial PDB Structure, Created by DYNAMO.

   -user              Read User-Supplied Restraint Files (Default).
   -nouser            Read Only DYNAMO Covalent Geometry Files.

   -vdw               Read and Use DYNAMO Van der Waals Table (Default).
   -novdw             No Van der Waals Terms Used.

DYNAMO Arguments for Executing a Simulated Annealing Schedule

dynSimulateAnnealing
   -sa   [saArgTriplets]    Simulated Annealing Protocol Arguments.
   -fc   [fcArgTriplets]    Force Constant Scaling Arguments.
   -size [sizeArgTriplets]  Size Scaling Arguments.

   -center                  Remove Center of Mass Motion at each step (Default).
   -nocenter                Do not remove motion.

   -graph                   Display Energy Curve Graphs.
   -nograph                 No Display of Energy Curves.

   -print  printCount       How often to print status.
   -rasmol drawCount        How often to draw current structure via rasmol.

   -tensorProc tclProcName  TCL Procedure for DC or PCS tensor calculation;
                            default proc is in dynamo/tcl/dynEvalTensor.tcl

The Ubiquitin NOE/DC Simulated Annealing Structure Calculation Demo

The files for this demo can be found in the "demo/ubiq" directory of the dynamo installation. The directory "orig" contains the experimental restraint tables for use with DYNAMO, along with some examples of table format conversion. The experimental tables used in this demo are:

noes.tab NOE Distances
jcoup.tab J-Coupling Values
torsions.tab TALOS phi/psi restraints
dObsA.tab Dipolar Couplings in Alignment Medium A
dObsB.tab Dipolar Couplings in Alignment Medium B

The script "all.com" goes through a complete example of NOE structure calculation and subsequent refinement including dipolar couplings.

The files for this demo can be found in the "dynamo/demos/ubiq" directory:

editGMC Creates DYNAMO Molecular Decscription by interactive specification of the sequence
init.com An alternative to interactively specifying the sequence; this script uses tools to extract sequence from existing file (in this case, a TALOS chemical shift file), and to create an initial extended structure. Output: ubiq.gmc GMC Directory and contents. ext.pdb Initial Extended structure.
ac.com Like "init.com", but creates a DYNAMO structure which conforms to an input PDB of the X-ray structure. Output: ubiq.gmc GMC Directory and contents. ac.pdb Conforms to X-ray structure. ext.pdb Initial Extended structure. ac.tab Restraints for conforming to PDB input. init.pdb Initial structure based on PDB input.
sa.tcl Computes 54 structures via high-temperature annealing, with DC terms off. Output: ubiq.com/dyn_*.pdb
avg.com Computes average structure of "sa.tcl" annealing results. Output: avg.pdb
pdbSelect.tcl Used by "avg.com" to select the lowest-energy structures produced by "sa.tcl".
refineAvg.tcl Refines average structure to restore proper geometry. Output: refinedAvg.pdb
refineDC.tcl Low temperature annealing to refine average structure by including dipolar couplings. Output: dc.pdb
ov.tcl Computes Backbone RMSD between calculated structure and reference structure, and also overlays both molecules in a single PDB output. Output: overlay.pdb
saDC.tcl Computes a series of structures via high-temperature annealing, with DC terms on. Output: ubiq.com/dynDC_*.pdb
simDC.tcl Simulate Dipolar Couplings for Given Structure
simCS.tcl Simulate Chemical Shifts for Given Structure

Some General Scripts for DYNAMO

Script Name Purpose
noeConvert.tcl Convert XPLOR NOE Table to DYNAMO Format
jConvert.tcl Convert XPLOR NOE Table to DYNAMO Format
dcConvert.tcl Convert XPLOR NOE Table to DYNAMO Format
torsionConvert.tcl Convert XPLOR NOE Table to DYNAMO Format
talos2dyn.tcl Convert TALOS Phi/Psi Table to DYNAMO Format
name2torsion.tcl Convert Named Torsions to DYNAMO Format
pdb2dyn.tcl Extract DYNAMO Residue Creation Info from PDB
pdb2ac.tcl Extract DYNAMO Atom Coord Restraints from PDB
pdb2torsion.tcl Extract DYNAMO Torsions Restraints from PDB
pdb2gmc.tcl Extract DYNAMO Sequence Information from PDB
seq2gmc.tcl Extract DYNAMO Sequence Information from Table
dynAvg.tcl Compute and Average PDB Structure
dynCenter.tcl Center a PDB Structure at the Origin
dynBasicExt.tcl Create an Extended Structure
dynEval.tcl Evaluate PDB According to Restraints
ov.tcl Overlay PDB files, Report Coord and Torsion RMS
addPDBNoise.tcl Add Random Structural Noise to PDB
addTabNoise.tcl Add Random Noise to a Table
dcNoise.tcl Dipolar Coupling Noise Analysis
pdbSelect.tcl Select DYNAMO PDB Based on Eneregy, Etc
resetPhiPsi.tcl Reset Protein Backbone Angles
ss.tcl Analyze Secondary Structure and H-Bond Info
mapPDB.tcl Map DYNAMO Table Values onto PDB
mfr.tcl Perform Molecular Fragment Search
mfr2init.tcl Create Initial Structure from MFR Angles
mfr2dyn.tcl Create Torsion Restraints from MFR Angles
dynAngles.tcl Display Backbone and Sidechain Angles
showCS.tcl Show Chemical Shift Table
showDC.tcl Show Dipolar Coupling Table
rotDC.tcl Add Dipolar Coupling Tensor Info to PDB
scrollRama.tcl Display Backbone Angle Trajectories
showTab.tcl Show X/Y Table Graphs

Creating New Residue Types

DYNAMO residue type parameters are located in a series of parameter directories defined by the environment variable DYNAMO_PARAMS. DYNAMO provides for three general classes of residues, with parameters in these default locations:


   dynamo/params/protein   Amino Acids
   dynamo/params/dna_rna   Nucleic Acids
   dynamo/params/other     Other Molecules

A given residue is defined by a TCL procedure in one of the above directories. The TCL procedure defines the atoms, bonds, angles, and torsions that describe the residue.

The master list of all the residues is stored in a table which defines the naming details of each residue, and the name of the TCL procedure which creates that residue:


   dynamo/params/params.tab

So, in order to create a new DYNAMO residue, there are two basic steps:

  1. Create a TCL procedure to define the new residue, and place it in the appropriate parameter directory (protein dna_rna or other).
  2. Edit the "params.tab" file to include an entry for the new residue.

TCL Residue Definition

As an example of how residues are defined in DYNAMO, consider the TCL procedure "add_ala" which defines an alanine residue, in the file "params/protein/ala":

proc add_ala { seg resID } \
{
   add_Atom $seg ALA $resID C    12.0 0.0903 3.2072   0.48 -0.471  0.723 -1.504
   add_Atom $seg ALA $resID CA   12.0 0.0903 3.2072   0.22 -0.376  0.324 -0.028
   add_Atom $seg ALA $resID CB   12.0 0.0903 3.2072  -0.30  0.908  0.867  0.597
   add_Atom $seg ALA $resID HA    1.0 0.0045 2.6157   0.10 -1.225  0.735  0.499
   add_Atom $seg ALA $resID HB1   1.0 0.0045 2.6157   0.10  1.720  0.809 -0.116
   add_Atom $seg ALA $resID HB2   1.0 0.0045 2.6157   0.10  1.157  0.278  1.467
   add_Atom $seg ALA $resID HB3   1.0 0.0045 2.6157   0.10  0.763  1.895  0.893
   add_Atom $seg ALA $resID HN    1.0 0.0498 1.4254   0.26  0.391 -0.265  1.965
   add_Atom $seg ALA $resID N    14.0 0.1592 2.7618  -0.10 -0.397 -1.128  0.109
   add_Atom $seg ALA $resID O    16.0 0.2342 2.6406  -0.48 -0.877 -0.080 -2.345

   add_Bond_Intrares $seg ALA $resID C    O     1.231 1000.000
   add_Bond_Intrares $seg ALA $resID CA   C     1.525 1000.000
   add_Bond_Intrares $seg ALA $resID CA   CB    1.521 1000.000
   add_Bond_Intrares $seg ALA $resID CA   HA    1.080 1000.000
   add_Bond_Intrares $seg ALA $resID CB   HB1   1.080 1000.000
   add_Bond_Intrares $seg ALA $resID CB   HB2   1.080 1000.000
   add_Bond_Intrares $seg ALA $resID CB   HB3   1.080 1000.000
   add_Bond_Intrares $seg ALA $resID N    CA    1.458 1000.000
   add_Bond_Intrares $seg ALA $resID N    HN     0.98 1000.000

   add_Angle_Intrares $seg ALA $resID CA   C    O    120.800 500.000
   add_Angle_Intrares $seg ALA $resID CA   CB   HB1  109.500 500.000
   add_Angle_Intrares $seg ALA $resID CA   CB   HB2  109.500 500.000
   add_Angle_Intrares $seg ALA $resID CA   CB   HB3  109.500 500.000
   add_Angle_Intrares $seg ALA $resID CB   CA   C    110.500 500.000
   add_Angle_Intrares $seg ALA $resID HA   CA   C    109.500 500.000
   add_Angle_Intrares $seg ALA $resID HA   CA   CB   109.500 500.000
   add_Angle_Intrares $seg ALA $resID HB1  CB   HB2  109.500 500.000
   add_Angle_Intrares $seg ALA $resID HB1  CB   HB3  109.500 500.000
   add_Angle_Intrares $seg ALA $resID HB2  CB   HB3  109.500 500.000
   add_Angle_Intrares $seg ALA $resID HN   N    CA   120.000 500.000
   add_Angle_Intrares $seg ALA $resID N    CA   C    111.200 500.000
   add_Angle_Intrares $seg ALA $resID N    CA   CB   110.400 500.000
   add_Angle_Intrares $seg ALA $resID N    CA   HA   109.500 500.000

   add_Improper_Intrares $seg ALA $resID HA   N    C    CB    65.977 500.000
   add_Improper_Intrares $seg ALA $resID HB1  HB2  CA   HB3  -66.514 500.000
}

The first section is a series of calls to add_Atom, which defines each atom in turn, supplying its segment name, residue name, residue number, atom name, standard mass, Lennard-Jones epsilon and sigma parameters, electric charge, and optional initial coordinates. (Note: the electric charge is not currently used by the DYNAMO energy function).

The second section defines all the covalent bonds. Each one is defined by the segment name, residue name, residue number, and the two atom names that are to be connected, together with the equilibrium distance (in Angstroms) and the bond's force constant (in kcal/mol-A^2).

The third section defines all the bond angle constraints. It works just like the bond length section, but three atom names are given, the equilibrium value is in degrees, and the force constant is in kcal/mol-radian^2.

The fourth section defines the improper torsion angle constraints. Impropers are torsion angles used to maintain chirality or planarity. They are defined just like the bond angles, except that four atom names are given.

Commonly, the TCL procedure will be created by copying a related existing one, and adjusting it manually. Alternatively, there is a crude script "pdb2dyn.tcl" which attempts to generate an initial TCL procedure based on a residue in an existing PDB file, and this can also be manually adjusted.

Adding the New Residue to the Master Table

In order for DYNAMO to be able to use the TCL procedure for a new residue, it must be listed in the master parameter table:


   dynamo/params/params.tab

The parameter table contains a header and comments to assist in adding new residue entries. In the case of alanine, the table contains this line with five items:


   ala protein add_ala ALA ala

The information about each residue is:

  1. Name to put on the graphical interface button.
  2. Residue type (must be "protein", "DNA", "RNA", or "other").
  3. Unique name of TCL procedure that can generate one of this residue type.
  4. Residue name that will be used in the PDB file.
  5. DYNAMO's unique internal residue name used during GMC creation.

Once the TCL residue procedure is created, and a suitable entry for the residue is inserted into the parameter table, DYNAMO will be able to create instances of that residue via "gmcEdit".

DYNAMO Torsion Angle Names and Definitions

Some DYNAMO facilities associate particular names with molecular torsion angles, for example, the commonly-used protein backbone angles phi and psi. A list of torsion angle names and definitions used by DYNAMO follows.

Residue Type Torsion Name Atom 1 Atom 2 Atom 3 Atom 4
All Amino Acids  phi   C(I--1) N(I) CA(I) C(I)
All Amino Acids  psi   N(I) CA(I) C(I) N(I-1)
All Amino Acids  alpha   CA(I--1) CA(I) CA(I-1) CA(I-2)
All Amino Acids  omega   CA(I--1) C(I--1) N(I) CA(I)
 
ALA  chi1   N CA CB HB1
ARG  chi1   N CA CB CG
   chi2   CA CB CG CD
   chi3   CB CG CD NE
   chi4   CG CD NE CZ
   chi5   CD NE CZ NH1
   chi6   NE CZ NH1 HH11
   chi7   NE CZ NH2 HH21
ASN  chi1   N CA CB CG
   chi2   CA CB CG OD1
   chi3   CB CG ND2 HD21
ASP  chi1   N CA CB CG
   chi2   CA CB CG OD1
CYS  chi1   N CA CB SG
   chi2   CA CB SG HG
GLN  chi1   N CA CB CG
   chi2   CA CB CG CD
   chi3   CB CG CD OE1
   chi4   CG CD NE2 HE21
GLU  chi1   N CA CB CG
   chi2   CA CB CG CD
   chi3   CB CG CD OE1
HIS  chi1   N CA CB CG
   chi2   CA CB CG ND1
ILE  chi1   N CA CB CG1
   chi2   CA CB CG1 CD
   chi3   CB CG1 CD HD1
LEU  chi1   N CA CB CG
   chi2   CA CB CG CD1
   chi3   CB CG CD1 HD11
   chi4   CB CG CD2 HD21
LYS  chi1   N CA CB CG
   chi2   CA CB CG CD
   chi3   CB CG CD CE
   chi4   CG CD CE NZ
   chi5   CD CE NZ HZ1
Residue Type Torsion Name Atom 1 Atom 2 Atom 3 Atom 4
MET  chi1   N CA CB CG
   chi2   CA CB CG SD
   chi3   CB CG SD CE
   chi4   CG SD CE HE1
PHE  chi1   N CA CB CG
   chi2   CA CB CG CD1
PRO  chi1   N CA CB CG
   chi2   CA CB CG CD
   chi3   CB CG CD N
   chi3   CG CD N CA
SER  chi1   N CA CB OG
   chi2   CA CB OG HG
THR  chi1   N CA CB OG1
   chi2   CA CB OG1 HG1
   chi3   CA CB CG2 HG21
TRP  chi1   N CA CB CG
   chi2   CA CB CG CD1
TYR  chi1   N CA CB CG
   chi2   CA CB CG CD1
VAL  chi1   N CA CB CG1
   chi2   CA CB CG1 HG11
   chi3   CA CB CG2 HG21
 
All Nucleic Acids  alpha   O3'(I--1) P(I) O5'(I) C5'(I)
All Nucleic Acids  beta   P O5' C5' C4'
All Nucleic Acids  gamma   O5' C5' C4' C3'
All Nucleic Acids  delta   C5' C4' C3' O3'
All Nucleic Acids  epsilon   C4'(I) C3'(I) O3'(I) P(I-1)
All Nucleic Acids  zeta   C3'(I) O3'(I) P(I-1) O5'(I-1)
All Nucleic Acids  v0   C4' O4' C1' C2'
All Nucleic Acids  v1   O4' C1' C2' C3'
All Nucleic Acids  v2   C1' C2' C3' C4'
All Nucleic Acids  v3   C2' C3' C4' O4'
All Nucleic Acids  v4   C3' C4' O4' C1'
 
ADE  chi   C4 C1 N9 C4
GUA  chi   C4 C1 N9 C4
CYT  chi   C4 C1 N1 C2
THY  chi   C4 C1 N1 C2
URI  chi   C4 C1 N1 C2
 



[ Home ] [ NIH ] [ NIDDK ]
last updated: Jan 30 2012 / big fd