7.2. Example 02: FreeLigand

This example demonstrates how to use the FreeLigand class to perform a multi-state RESP calculation.

The key difference between this and Example 01: LazyLigand is that this example uses the ligandparam.stages.resp.StageMultiRespFit stage to perform a multi-state RESP calculation after performing gaussian rotations using ligandparam.stages.gaussian.StageGaussianRotation. These types of rotations ared useful for ensuring that the choice of grid points in the RESP calculation is not biased by the initial orientation of the molecule.

7.2.1. Learning Outcomes:

  1. Learn how to use the FreeLigand class to perform a gaussian rotational RESP calculation to parametrize a ligand using multi-state resp fitting.

  2. Learn how to list and execute stages in the pipeline.

7.2.2. Files

The files for this example can be found in the LigandParameterization/examples/02_FreeLigand directory of the source code.

7.2.3. Tutorial

To start with any of the premade recipes, you need to import the recipe from the module. Using the wild-card character will import all the recipes available in the module, whereas, importing a specific recipe will only import that recipe. Just like before, there is also a default stage list, which can be used to disable stages from the recipe, if needed. This only needs to be defined if you plan to disable stages from the recipe.

from ligandparam.recipes import *

default_stage_list = {
    "Initialize": True,
    "Normalize1": True,
    "Minimize": True,
    "Rotate": True,
    "GrabGaussianCharge": True,
    "MultiRespFit": True,
    "UpdateCharge": True,
    "Normalize2": True,
    "UpdateNames": True,
    "UpdateTypes": True,
    "ParmChk": True,
    "Leap": True,
}

The next step is to load the pdb file into the recipe of your choosing, and set various machine parameters. Full parameters that can be selected are available in the documentation for the class (for instance, for FreeLigand, see ligandparam.module.FreeLigand).

inputoptions = {
    'base_name': 'thiophenol',
    'net_charge': 0,
    'mem': '60GB',
    'nproc': 12
}

# Load the pdb as a instance of the FreeLigand class
test = FreeLigand(inputoptions = inputoptions)

Note that here, we have chosen the base_name thiophenol and set the net charge to 0 using the netcharge parameter. If you have a charged ligand, you need to select the proper net charge, as later when Gaussian runs it will need to know the net charge of the molecule. The nproc and mem parameters are used to set the number of processors and the memory to be used by Gaussian, respectively. Note that this class is called nearly identically to the LazyLigand class, but with the FreeLigand class used instead.

The next step is to select the pre-initialized stages for the FreeLigand class. This can be done using the setup method.

test.setup()

The FreeLigand class has a number of stages that are pre-initialized.

In brief, these are:

  1. ligandparam.stages.initialize.StageInitialize - Generates the initial mol2 file with bcc charges using antechamber.

  2. ligandparam.stages.charge.StageNormalizeCharge - Ensures that the bcc charges sum to the net charge, and tries to correct for any errors.

  3. ligandparam.stages.gaussian.StageGaussian - Runs Gaussian to minimize the molecule.

  4. ligandparam.stages.gaussian.StageGaussianRotation - Runs Gaussian to rotate the molecule and generate multiple conformations.

  5. ligandparam.stages.resp.StageMultiResp - Uses ligandparam.multiresp to generate a mol2 with the RESP charges from the Gaussian calculation.

  6. ligandparam.stages.charge.StageNormalizeCharge - Ensures that the bcc charges sum to the net charge, and tries to correct for any errors.

  7. ligandparam.stages.typematching.StageUpdate - Updates the atom names to match the original antechamber atom names in the calculation.

  8. ligandparam.stages.typematching.StageUpdate - Updates the atom types to match the original antechamber atom types in the calculation.

  9. ligandparam.stages.parmchk.StageParmChk - Generates the frcmod file for the ligand using parmchk2.

  10. ligandparam.stages.leap.StageLeap - Runs tleap to generate the final .off parameter files for the ligand.

To list the stages out to the user, you can use the list_stages method.

test.list_stages()
Finally, to execute the stages in order, you can use the execute method. The dry_run parameter is used to test the pipeline

without actually creating any files. This is useful to check if the pipeline is working as expected; however, it has limited functionality as many stages depend on files generated by previous stages.

test.execute(dry_run=False)

This will run the pipeline in order, generating the necessary files for the ligand parameterization.

The output files will be generated in the same directory as the input pdb file, and will have the same name as the pdb file, but with different extensions.

These files are:

  • thiophenol.resp.mol2 - The final mol2 file with the RESP charges.

  • thiophenol.frcmod - The frcmod file for the ligand.

  • thiophenol.off - The off(lib) parameter file for the ligand.

Note that the charges ion the mol2 file should be similar, but not exactly the same as the charges that you obtained from the LazyLigand class.

7.2.4. Full code

# Import the module
from ligandparam.recipes import FreeLigand

# Example default stage list, which could be passed to the disable_stages method to mass remove stages from
# the recipe. To do that, you would uncomment the line marked by a commment.
default_stage_list = {
    "Initialize": True,
    "Normalize1": True,
    "Minimize": True,
    "Rotate": True,
    "GrabGaussianCharge": True,
    "MultiRespFit": True,
    "UpdateCharge": True,
    "Normalize2": True,
    "UpdateNames": True,
    "UpdateTypes": True,
    "ParmChk": True,
    "Leap": True,
}

inputoptions = {
    'base_name': 'thiophenol',
    'net_charge': 0,
    'mem': '60GB',
    'nproc': 12
}

# Load the pdb as a instance of the FreeLigand class
test = FreeLigand(inputoptions = inputoptions)

# Select the pre-initialized stages for Lazy Ligand
test.setup()

# Disable stages from the default list
#test.disable_stages(default_stage_list)

# List the stages out to the user
test.list_stages()

# Execute the stages in order.
test.execute(dry_run=False)