Examples

In this section, we provide examples that illustrate how to use LEAP under a variety of conditions and for several different molecules. For each example, we first enumerate the principles that are demonstrated then print the log file. The log files have been edited in order to provide comments and explanations. The user may wish to recreate these examples and can do so by running the example log files found in the LEAP "demos" directory. (Refer to the "Log File" section of the "Concepts" section.)

A Simple Steroid: Cholesterol

In this example, the steroid cholesterol is created within xleap. The example illustrates:

• ·creating a LEAP UNIT using the Unit Editor
• ·using the Build and Relax Selection commands
• ·modeling a cyclic molecule
• ·creating LEAP OFF libraries.

This example is mainly a discussion of how to use the Unit Editor Build and Relax Selection commands to create a complex cyclic molecule. There are two points to mention that the user needs to know in order to create such a structure in the Unit Editor. First, since one builds molecules on a two-dimensional surface in xleap, structures that are minimized with the Relax Selection command may be planar unless they are first altered in the Manipulate mode or generated using the Build command. Second, the Build command frequently improperly generates fused cyclic structures. In the following example, the cholesterol molecule is drawn in the order: ring A, B, C, D, and alkyl chain. The LEAP structure was saved in an OFF library ("chol.lib") at each stage of refinement. The user is encouraged to look at the library entries if they have difficulty constructing the molecule. The user should note that this version of LEAP does a less than ideal job of generating a molecule such as cholesterol since vdw and electrostatic force field terms are not used in the Unit Editor minimizer.

> #
> # First, execute the commands in the distribution
> # "leaprc" file. Because of space considerations,
> # these commands have been deleted from this example.
> # Next, create a new UNIT and edit it in the Unit
> # Editor. Choose the "Draw" button and the drawing
> # element of carbon. Draw a six-membered ring (ring
> # A). Select the menu item "Add H & Build". Erase
> # two hydrogens and draw the ring A hydroxyl oxygen
> # and methyl carbon atoms. Select the menu item
> # "Add H & Build". Save the "ch1" UNIT in an OFF
> # library:
> #
> edit ch1
> saveOff ch1 chol.lib
> #
> # Make a copy of the "ch1" UNIT so we can leave
> # that UNIT as an illustration. Edit the copy,
> # "ch2", in the Unit Editor in order to add ring
> # B. First, remove the three hydrogens at the
> # ring fusion carbons. Choose the "Draw" button
> # and the drawing element of carbon. Draw a
> # six-membered ring (ring B). Select the ring B
> # carbons not associated with the ring A junction.
> # Select the menu items "Relax selected atoms" and
> # "Mark built". Save the UNIT in an OFF file:
> #
> ch2 = copy ch1
> edit ch2
> saveOff ch2 chol.lib
> #
> # Make a copy of the "ch2" UNIT so we can leave that
> # UNIT as an illustration. Edit the copy,
> # "ch3", in the Unit Editor in order to add ring
> # C. Choose the "Draw" button and the drawing
> # element of carbon. Draw a six-membered
> # ring (ring C). Select the ring C
> # carbons not associated with the ring B junction.
> # Select the menu items "Relax selected atoms" and
> # "Mark built". Save the UNIT in an OFF file:
> #
> ch3 = copy ch2
> edit ch3
> saveOff ch3 chol.lib
> #
> # Make a copy of the "ch3" UNIT so we can leave that
> # UNIT as an illustration. Edit the copy,
> # "ch4", in the Unit Editor in order to add ring
> # D. Choose the "Draw" button and the drawing
> # element of carbon. Draw a five-membered
> # ring (ring D). Select the ring D
> # carbons not associated with the ring C junction.
> # Select the menu items "Relax selected atoms" and
> # "Mark built". Save the UNIT in an OFF file:
> #
> ch4 = copy ch3
> edit ch4
> saveOff ch4 chol.lib
> #
> # Add the alkyl group to ring D. Select the
> # menu items "Relax selected atoms" and
> # "Mark built". Save the UNIT in an OFF file:
> #
> ch5 = copy ch4
> edit ch5
> saveOff ch5 chol.lib
> #
> # Add all remaining hydrogens to the molecule.
> # Select the menu items "Relax selected atoms"
> # and "Mark built". Save the UNIT in an OFF file:
> #
> ch6 = copy ch5
> edit ch6
> saveOff ch6 chol.lib
> #
> # Quit the xLEaP program:
> #
> quit

An Ion-Molecule Complex: 18-Crown-6 and Potassium Cation

In this example, the crown ether/ion complex 18-crown-6/K⁺ is created within xleap and then solvated with methanol. The example illustrates:

• ·creating a LEAP UNIT starting from AMBER PREP input files
• ·modeling a cyclic molecule
• ·loading a PDB file into a LEAP UNIT
• ·building a box of solvent molecules
• ·solvating a molecule
• ·converting an AMBER parameter file to a LEAP PARMSET
• ·creating LEAP OFF libraries
• ·creating topology and coordinate files for AMBER or LEAP.

  > #
  > # First, execute the commands in the distribution
  > # "leaprc" file. Because of space considerations,
  > # these commands have been deleted from this example.
  > # Next, 3 UNITs will be created from AMBER PREP
  > # files:
  > #
  > loadAmberPrep cra.in
  > loadAmberPrep crb.in
  > loadAmberPrep crc.in
  > #
  > # A PDB file is loaded into the UNIT "cr":
  > #
  > cr = loadPdb crown.pdb
  > #
  > # The "cr" UNIT is edited in the xLEaP X-Windows
  > # Unit Editor and a bond is added between RESIDUEs
  > # CRA(1) and CRB(6). This bond could also have been
  > # created with the command line editor using the command:
  > # "crossLink cr.1 lastend cr.6 firstend".
  > # The same editor is then used to create
  > # a potassium ion and assign atom properties to
  > # that ion. The ion is placed at the center of the
  > # crown ether. Finally, a methanol UNIT ("meoh") is created
  > # in the Unit Editor. Atom Properties are defined
  > # for the methanol UNIT and the ATOMs in this UNIT
  > # are minimized using the Unit Editor "Relax Selection"
  > # menu item:
  > #
  > edit ion
  > edit cr
  > edit meoh
  > #
  > # We will then create methanol solvent box.
  > # In order to do so, it is first
  > # necessary to align the principle axes of the
  > # "meoh" UNIT along the Cartesian coordinate axes
  > # with its geometric center at ( 0, 0, 0 ).
  > # The bounding box is defined to be a 4 angstrom
  > # cube:
  > #
  > alignAxes meoh
  > set meoh box 4
  > set meoh restype solvent
  > #
  > # Parameters will be needed for molecular mechanics.
  > # The parameters had previously been used in the
  > # AMBER modules and therefore an AMBER parameter
  > # file is used to load the 18-crown-6 parameters.
  > # We now will solvate the crown with methanol
  > # molecules. The original crown/ion UNIT will be
  > # copied so the "cr" UNIT will not be altered by solvent:
  > #
  > cr_parms = loadAmberParams cr_parms.mod
  > cr-tem = copy cr
  > solvateBox cr-tem meoh 10
  > #
  > # The parameters, and all of the UNITs that have
  > # been created within xLEaP, are then saved as
  > # OFF library files. Finally, topology and coordinate
  > # files are created for use in AMBER or SPASMS:
  > #
  > saveOff cr_parms crown.lib
  > saveOff meoh crown.lib
  > saveOff cr crown.lib
  > saveOff ion crown.lib
  > saveOff cr-tem crown.lib
  > saveAmberParm cr-tem cr-meoh.top cr-meoh.crd
  > #
  > # Quit the xLEaP program:
  > #
  > quit
  

Free Energy Perturbation: Guanine To Adenine

One type of calculation frequently performed in molecular mechanics is that of Free Energy Perturbation. During Free Energy Perturbation, the potential energy function of one molecular species is slowly perturbed to that of a second molecular species through a perturbation parameter called $lambda$. This calculation allows one to determine the relative Gibbs or Helmholtz free energy differences between the two species. In this example, LEAP is used to produce AMBER/SPASMS input files for the perturbation of Guanine to Adenine, in vacuo and in water, using the 1991 force field residues. This perturbation involves the a) replacement of GUA ( (C_{(6)} = O_{(6)}) ) into ADE ( (C_{(6)} - NH_{2}) ), b) replacement of GUA ( (N_{(1)} - H) ) into ADE ( (N_{1}) ), and c) replacement of GUA ( (C_{(2)} - N_{(2)}H_{(2)}) ) into ADE ( (C_{2} - H) ).

This example illustrates:

• ·creating Free Energy Perturbation input files within LEAP
• ·loading a AMBER PREP RESIDUE into LEAP
• ·converting an AMBER parameter file to a LEAP PARMSET
• ·creating LEAP OFF libraries
• ·creating topology and coordinate files for AMBER or LEAP.

  > #
  > # First, execute the commands in the distribution
  > # "leaprc" file. Because of space considerations,
  > # these commands have been deleted from this example.
  > # An AMBER PREP residue for Guanine (GUB) is loaded
  > # into LEaP:
  > #
  > loadAmberPrep guan.in
  > #
  > # An AMBER format parameter set file is loaded and placed
  > # in "parmatgc.dat" file. The resulting PARMSET is saved
  > # as an OFF file:
  > #
  > newparm = loadAmberParams parmatgc.dat
  > saveOff newparm newparm.lib
  > #
  > # The GUB UNIT is edited in the Unit Editor to
  > # prepare the perturbation properties. The user
  > # should examine the GUB UNIT in "gub.lib",
  > # using the Unit Editor and Atom Properties
  > # Editor, to see the changes that are made for
  > # perturbation. The UNIT is saved as an OFF file
  > # and coordinate and topology files are prepared
  > # for AMBER and SPASMS :
  > #
  > edit GUB
  > saveOff GUB gub_leap.lib
  > saveAmberParm GUB gub_vacuo.top gub_vacuo.crd
  > #
  > # The GUB UNIT is then solvated with a box of
  > # TIP3P waters.
  > # The GUB UNIT is saved in an OFF library file.
  > # Topology and coordinate files for AMBER or
  > # SPASMS perturbation calculations (TIP3P
  > # solvated) are created:
  > #
  > loadOff water.lib
  > solvateBox GUB WATERBOX216 10
  > saveAmberParm GUB gub_water.top gub_water.crd
  > #
  > # Quit the program
  > #
  > quit
  

Creating Polynucleotides: B DNA

A double helix of dna is built using LEAP in this example. Two examples of adding counterions to the dna are illustrated. In one example, both counterions and tip3p water are added. The second example is an in vacuo model and the counterions have a large vdw radius in order to simulate a water shell surrounding them.

The Cartesian coordinates for the dna were obtained from the PDB file "1bd1". This file also contained coordinates for X-ray waters of crystallization and two triethylammonium ions. Neither of these species are included in this representation. Only one nucleotide strand is found in the PDB file since the symmetry operations necessary to generate the second strand are given. These operations are used with the LEAP transform command in order to create the double helix molecule.

This example illustrates:

• ·creating a dna molecule within LEAP
• ·loading a PDB file into a LEAP UNI
• ·using PDB symmetry operations to build a molecule
• ·using the addIons command to place counterions around a molecule
• ·solvate a molecule using the LEAP solvateBox command
• ·creating LEAP OFF libraries
• ·creating topology and coordinate files for AMBER or LEAP.

  #
  dnaStrandA = loadPdb 1bd1.pdb
  #
  # A second strand is generated by making a
  # copy of the first strand and then applying
  # the symmetry operations found in the PDB
  # file to the latter strand. The DNA double helix
  # is generated by combining the two strands.
  # The combine command does not link the UNITs
  # in the LIST argument with bonds:
  #
  dnaStrandB = copy dnaStrandA
  transform dnaStrandB {
  { -1 0 0 }
  { 0 1 0 }
  { 0 0 -1 }
  }
  dnaDoubleHelix = combine { dnaStrandA dnaStrandB }
  #
  # Two complexes are now created. The first is
  # the "in vacuo" association of a DNA double
  # strand with ions. Since this complex will not be
  # solvated, the ions have large VDW radii in order to
  # simulate a shell of water surrounding them. The
  # second complex associates DNA, ions, and a box of
  # TIP3P waters. The ions (Na+) in the later case
  # have VDW radii of "standard" size. In both
  # complexes, enough counterions are added to
  # neutralize the DNA:
  #
  loadOff ions.lib
  dnaDHBigIons = copy dnaDoubleHelix
  addIons dnaDHBigIons IB 0
  #
  dnaDHIons = copy dnaDoubleHelix
  solvateBox dnaDHIons WATBOX216 9
  addIons dnaDHIons Na+ 0
  #
  # The DNA double helix UNITs are saved in
  # an OFF library file. Topology and coordinate
  # files for the double helix molecules with associated
  # counterions/water are created for use in AMBER:
  #
  saveOff dnaDoubleHelix dna.lib
  saveOff dnaDHBigIons dna.lib
  saveOff dnaDHIons dna.lib
  saveAmberParm dnaDHBigIons dnaDHBigIons.top dnaDHBigIons.crd
  saveAmberParm dnaDHIons dnaDHIons.top dnaDHIons.crd
  #
  # Quit the program
  #
  quit
  

A Protein/Ligand Complex: trp Repressor

This example illustrates one method of building the trp repressor protein and solvating it within LEAP. The trp repressor protein is dimeric. Each protein monomer contains one zwitterionic trp amino acid ligand and one polypeptide unit. The two monomeric units are related by symmetry.

In the example, a PDB structure (1wrp) is used to generate the trp repressor Cartesian coordinates. The PDB file contains coordinates for only one monomer; the second monomer must be generated from the symmetry operations given in the PDB file. The original PDB file also contained Cartesian coordinates for water molecules found in the crystal structure. We have removed the waters from the PDB file as we did not want to use them in the example. In addition, we have split the original PDB file into two files: the file "wrp.pdb" contains one polypeptide chain and the file "trp.pdb" contains the trp ligand coordinates.

The example illustrates:

• ·creating a protein within LEAP
• ·loading a PDB file into a LEAP UNIT
• ·using PDB symmetry operations to build a molecule
• ·the use of the copy command to simplify the creation of a molecule
• ·using each of the solvation options of LEAP
• ·creating LEAP OFF libraries
• ·creating topology and coordinate files for AMBER.

  #
  # First, execute the commands in the distribution
  # "leaprc" file. Because of space considerations,
  # these commands have been deleted from this example.
  # Load the Cartesian coordinates of one protein
  # monomer into a LEaP UNIT:
  #
  wrpA = loadPdb wrp.pdb
  #
  # The trp repressor ligand is a zwitterionic amino acid.
  # Since the LEaP libraries do not include
  # zwitterionic amino acids, the ligand will be
  # modeled by making a copy of "NTRP" (N-terminal TRP)
  # and modifying the C-terminal to be an anion. The
  # Unit Editor will be used to make the modification
  # and to add the Atom Properties of the modified ATOMs.
  # First, a copy will be made of "NTRP". The UNIT and
  # residue names of the copy will then be changed from
  # "NTRP" to "WRP". Finally, the "WRP" UNIT will be
  # modified within the Unit Editor. The following changes
  # are made - 1) a second oxygen is drawn and connected to
  # the C-terminal carbon; 2) the C-terminal carboxylate
  # ion is selected and the "Build" command is chosen
  # from the "Build menu"; 3) the "Edit selected atoms"
  # command is chosen from the "Selection menu". An Atom
  # Properties table will appear for the carboxylate ion group.
  # In the # table, the "Name" of the new oxygen is changed to
  # "OXT" and its "Type" is changed to "O". The "Charge"s
  # of both carboxylate ion oxygens are changed to -0.7525.
  # The Atom Properties table and Unit Editor are then exited
  # and the "trpA" Cartesian coordinates are loaded from a PDB
  # file:
  #
  WRP = copy NTRP
  set WRP name "WRP"
  set WRP.1 name "WRP"
  edit WRP
  trpA = loadPdb trp.pdb
  #
  # Create one monomer by associating the protein
  # chain with the zwitterion. Then, to create the
  # second monomer, copy the first monomer and apply
  # the symmetry operations found in the PDB file.
  # Finally, combine the two monomers to form the
  # protein:
  #
  monomerA = combine { wrpA trpA }
  monomerB = copy monomerA
  translate monomerB { 0 0 -24.53333 }
  transform monomerB {
  { 1 0 0 }
  { 0 -1 0 }
  { 0 0 -1 }
  }
  trpRepressor = combine { monomerA monomerB }
  #
  # Solvate the protein using the various LEaP
  # solvation options. Prior to each solvation,
  # a copy of the protein UNIT is made so that
  # the original UNIT is not altered by the water
  # molecules. The solvent used in each of the
  # examples is TIP3P water taken from the
  # water OFF distribution library. In the
  # case of the "solvateCap" command, the cap is
  # centered at the Calpha carbon of the "trpA" ligand
  # in "monomerA". Two of the solvation commands are
  # shown below but not executed in the example:
  #
  # trpRepressorDontClip = copy trpRepressor
  # solvateDontClip trpRepressorDontClip WATBOX216 8.0
  # trpRepressorBox = copy trpRepressor
  # solvateBox trpRepressorBox WATBOX216 8.0
  #
  trpRepressorCap = copy trpRepressor
  solvateCap trpRepressorCap WATBOX216
  trpRepressorCap.106.CA 8.0
  trpRepressorShell = copy trpRepressor
  solvateShell trpRepressorShell WATBOX216 8.0
  #
  # The monomer, ligand, and trp repressor UNITs are
  # saved in an OFF library file. Topology and coordinate
  # files are created for use in AMBER or SPASMS:
  #
  saveOff wrpA trpRepressor.lib
  saveOff trpA trpRepressor.lib
  saveOff trpRepressor trpRepressor.lib
  saveAmberParm trpRepressor trpRepressor.top trpRepressor.crd
  #
  # Quit the program
  #
  quit