VMD-L Mailing List

From: Andrew Dalke (dalke_at_ks.uiuc.edu)
Date: Fri Mar 14 1997 - 16:40:47 CST

Continuing from last week we will finish describing some of the
details of how VMD makes its images. For the people who recently
joined the list, VMDTech is a weekly article we write that describes
some aspect of VMD in more detail. They are currently archived as
part of the general vmd-l archive at
   http://www.ks.uiuc.edu/Research/vmd/mailing_list/vmd-l/

If you have a suggestion for an issue or a contribution, please
contact us at vmd_at_ks.uiuc.edu.

BTW, you might be interested to know there are now 78 people on the
VMD mailing list.

                                                Andrew Dalke
                                                dalke_at_ks.uiuc.edu

 = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = =
                        Drawing Methods, Part II

  Last week we covered all of the full-atom and full-bond drawing
representations (Lines, Bonds, CPK, Points, VDW, Dotted and Licorice)
and Solvent, which is almost the same as Points. We also describe one
reduced atom representations, Trace, which draws the C-alpha trace.
This week we will finish off the other representations available in
VMD: Tube, Ribbons, Cartoon, HBonds, Surf and Alpha.

 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
                                Tube

  Last week we described how the Trace method finds the sequential
C-alpha coordinates that make up a protein. In that representation,
cylinders are drawn between successive residues to indicate that they
are bonded to each other. It is possible, however, to use information
about the previous residue before and the next residue after the two
C-alphas of the connection in order to produce a smoother curve.

  Shown here is a picture of the system to be dicussed. The
connection is to be drawn between the 2nd and 3rd C atoms, as
indicated with the double \\ bond, along with the two residues to
either side.

           C2 C4
          / \\ / \
   . / \\ / \
    \ / ^ \\ / ^ .
     \ / /|\ \\ / /|\
      C1 | C3 |
           |_____| |
      ^ Used for Trace
     /|\ |
      |________________|
              Also used for Tube

With two points you can fit a line, but with four points you can fit a
3rd-order function. Curves of order N-1 that are used to fit N data
points are called spline curves. There are many classes of splines,
and several of them have been used to fit protein backbondes. For the
most part these use 4 points, so they are cubic polynomials,
parameterized by some variable t:

        -> -> -> -> 2 -> 3
        X(t) = a0 + a1*t + a2*t + a3*t 0 <= t <= 1

All that remains is to find the coefficients, and that's where the
interesting things lie :) Of course, if a2 and a3 are 0, this
describes a line, which is the C-alpha trace. However, consider the
continuity of the connection from C1 to C2 with the connection from C2
to C3. The two endpoints meet, but the curve isn't smooth across C2
(I would say that it isn't C(1), but there's starting to be a few too
many C's hanging around.)

  The key point is that we want the curve between the residues to be
connected and smooth, so (for example) the curve from C1 to C2 meets
the curve from C2 to C3, and their derivatives are also continuous at
the point C2. If you want more details on how splines work, I strongly
advise you to read one of the many books on computer graphics; I
learned al of this from Foley and van Dam.

  There are many ways of constructing the coefficients based on the
different bondary conditions. The standard B-spline goes through the
1st and 4th points and uses the 2nd and 3rd points to determine an
intermediate slope, but the curve does not go through the middle two
points; instead, it is 2nd order continuous at the endpoints. For
some cases, like where you need smooth curves, this is useful. This
method (or "basis") has been used in several molecular modeling
programs, however, most people want a curve that goes through all of
the C-alphas so we needed another basis. (Actually, some programs try
to get the curve closer to correct by scaling a few parameters, but
this is highly dependent on the secondary structure type becuase the
factor for helicies is different than the one for beta sheets.)

  The Catmull-Rom basis goes through all four points, but by adding
this constraint we have to get rid of the requirement for second order
continuity. (A basis, to be technical, is the matrix used to
transform the input coordinates into a set of coefficients.) Also,
earlier I said "the 2nd and 3rd points [...] determine an intermediate
slope." For a Catmull-Rom basis, the 1st and 3rd points determine the
slope of the curve through the 2nd point, and the 2nd and 4th points
determine the slope of the curve through the 3rd points.

  However, this also isn't all that good looking. As it turns out,
this set of contraints causes the curve to bend too much around
helicies. It is enough that I can look at the output of some programs
and say "oh, they use an unmodified Catmull-Rom basis." The trick is
there is an additional scaling factor available which determines the
magnitude of the slope between the two points. By scaling the slope
up a bit, the curve looks much nicer.

  If you look at the code for DrawMolItem.C you can actually see the
different basises we tried out. Someday if you have the time you
might even want to try them out :)

  Now that the curve is parameterized, VMD breaks the curve down into
6 parts and draws each as a cylinder. As with Bonds, the radius of
the cylinder and the number of sides to use in the approximation can
be modified by the user. If either are too small, the segments are
instead drawn as lines. The cylinders are colored according to the
C-alpha of the residue being drawn.

 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
                                Ribbon

  The path of the ribbon is the same as the path of the tube, but we
need some way to determine the orientation of the cross section of the
ribbon. Everyone I know of uses the oxygen of the backbone, and so
does VMD. As one solution you can use the direction of the C=0 bond
to determine the perpendicular, but this causes some toe-in for
helicies -- the bottom of the ribbon point too much into the helix and
the top point out. You would expect that the ribbons for a helix are
roughly parallel to the helix axis, so the parameters need (again) to
be tweaked. (BTW, it is also possible to look at the direction of the
ribbons in some images and say "oh, they orient along the oxygen".)

  There are many possible ways to modify adjust the ribbon normal. I
tried a lot, but none worked very well. In frustion I looked at the
ribbon code from the Raster3D distribution and (with permission from
the author, Ethan Merritt) use that. This one keeps a running sum.
The new normal is computed as a sum of the current oxygen direction
and a scaled down contribution from the previous normal, for a nice
result. As a historical viewpoint, the Raster3D code itself came from
FRODO, which was one of the first molecular visualization programs.

  With the path and one normal computed it is a relatively simple
matter of making the ribbon. In addition, we add a cylinder along
each edge to give it more of a definition. The ribbon width is
user-adjustable and the cylinders are adjustable just like the regular
Tube is, except that if the radius or resolution is too small, nothing
is drawn (since the line of the tube is identical to the edge of the
ribbon which is already being drawn). Also, as with the Tube, each
ribbon section is colored according to the C-alpha to which it is
associated.

 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
                        Cartoon (nee' Structure)

  The Cartoon representation (called Structure in the 1.1 release)
draws solid ribbons for the beta sheets, cylinders for the helicies
and tube for the coil. Of course, this requires the VMD have
knowledge of secondary structure, and for information on how that
works, see the earlier VMDTech titled "Secondary Structure".

                                Alpha Helix

  The method is pretty easy. For each set of N contiguous helical
residues, compute the best-fit line through the C-alphas. Draw N
cylinders along the line, equally spaced, one for each residue. Color
each cylinder by the C-alpha of the residues associated with it. The
radius and resolution of the helix can be modified by the user, and
the default is set that the radius of the cylinder is about the same
as the radius of the helix. The only real complication that arises is
when the helix is drawn with a transparent color you can see the ends
of the cylinders in the middle of the helix. We're hoping to modify
the code so those extra disks aren't drawn.

                                Beta Sheet

  The structure of a beta sheet looks somewhat like this:

   ------------------------------------
         | CA | CA CA
         | / \ | / \ / \
         |/ \| / \ / \
   \ * * / \ /
    \ /| |\ / \ /
     \ / | | \ / \ /
      CA | | CA CA
   ------------------------------------

Where the --- and | lines will be used to indicate how this gets
turned into a solid ribbon. As with the Ribbon section, we need to
compute a path and a normal for the ribbon. The path is pretty
simple, it goes along the points halfway between successive C-alpha
atoms. These points are marked with an * in the above picture.

The normal is found by computing the cross product of two vectors.
The normal for the leftmost * is the cross product of the vector
between the 1st and 2nd CA atoms, and the vector between the 2nd and
3rd CA atoms. There is also a bit of extrapolation (trickery) to get
the normals for the points at the start and end of the ribbon.

Given the three vectors that define the sheet it is a pretty simple
matter of computing the corners for the ribbon that goes through each
section of the ribbon, and from there to draw the boxes connecting
each section. However, as usual, there are few more minor details
because beta sheets are not in the ideal conformation shown above;
they can twist quite a bit. If the normal vectors for the two ends of
the box are inverted, the later one is flipped. If the dot product of
the two normals exceedes some threshold, an intermediate coordinate is
computed to make the ribbon a bit smoother.

The ribbon width can be modified by the user. As always, each section
is colored by the color of the associated C-alpha atom.

                                Coil

  The residues that are neither helix nor coil, and those marked as
beta sheets by members of strands of length 2 or less, are drawn as a
Tube. In fact, we even call the draw_tube() function, but we do use
one trick. The Tube goes through the C-alpha atoms, but think about the
helix. You would like the Tube to go to the center of the end of the
helix cylinder, but in reality the protein wraps around the helix, so
it would actually go to the edge of the end. A simular situation
exists for the beta sheets.

  To fix this, we temporarily alter the coordinates of the C-alpha
atoms at the end of the helix cylinders and move them to the center of
the end faces, draw the Tube, and return the atoms to their original
coordinates. One way to see this is to draw things as Tube and also
draw things as Cartoon. You'll see the difference at the ends of the
helicies and sheets.

 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
                                HBonds

  The HBond representation will draw a dotted line between two atoms
if there is a possible hydrogen bond between them. A possible
hydrogen bond is defined as the following poem (this text is getting a
bit long :)

  Given an atom D with a hydrogen H bonded to it,
  And an atom A with no hydrogen bonded to it,
  A hydrogen bond exists between A and H iff,
  The distance ||D-A|| < dist and
  The angle D-H-A < ang,
  Where ang and dist are user-defined.

  Only the selected atoms are searched, so both the donor and acceptor
must be selected for the bond to be drawn. Also, you'll note that the
above doesn't check the atom type of the donor or acceptor; the only
criterion is if it already has or doesn't have a hydrogen. One
downfall, BTW, of the current implementation is that it does an
n-squared search of the selected atoms so you probably don't want to
show all the HBonds of a very large structure.

 There are several reasons why you might not see H-bonds, but the most
likely is that this algorithm didn't find any. The default angle and
distance criterion in VMD are too small, so you might want to try
increasing the angle value from 20 to 30 degrees and the distance
value from 3 to 4. Looking at 9pti with HBonds shows nothing with the
default and several with those parameters increased.

  The HBonds are drawn as dotted lines of a given width. The default
is 1 but you should probably increase that to 2. On most SGIs you
can't make it any wider than that, as described in the man page for
linewidth). We should also add a way to use cylinders in addition to
(instead of?) lines so they can be made more obvious by increasing
the radius. The bond is colored by the color associated with the
acceptor.

 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
                                Surf

  VMD does not directly compute the molecular surface. Instead, we
use Amitabh Varshney's SURF program (for more information please see
ftp://ftp.cs.unc.edu/pub/projects/GRIP/SURF/surf.tar.Z, which is also
in the VMD source distribution in lib/surf .) To use it, we write out
the atom coordinates and radii (remember from last week that you can
modify the radii yourself!) then call the appropriate SURF executable.
When it finishes, it writes a set of triangles to a file which VMD
reads and displays.

  The advantage for us of using someone else program is it was much
easier to incorporate his code than redeveloping it ourselves. (This
is the same reason we use STRIDE to compute the secondary structure.)
There are some disadvantages, and the biggest of these is it takes a
long time to write all the triangles to a file then read them back in.
Our estimates are that about half the time is spend in disk
reads/writes. Thus, we are looking for other means for computing the
surface. The most likely of these is Michael Sanner's MSMS code, from
Olsen's lab at Scripps.

  The only user changeable parameter to this option is the probe
radius. Every time you change the radius or move a coordinate, the
surface is recomputed. How are the colors determined? When the list
of triangles that make up the surface is computed, an atom index is
printed out along with the triangle data. VMD maps the atom index
back to the appropriate atom and finds and uses the appropriate atom
coloring.

 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
                                Alpha

  You probably don't even have this option on your version of VMD.
This is primarily a test interface to the Alpha Shapes code of Herbert
Edelsbrunner. We provide this code for the 1.2 version, but most of
the code is not redistributable so we have to check into it. I can't
even really tell you about how it works, so you might want to see the
Alpha Shapes web page at http://alpha.ncsa.uiuc.edu/alpha/index.html
instead.

 = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = =