How-to¶

This notebook serves as a practical guide to common questions users might have. Some functionalities of MolDrug may not be available on PyPi or Conda yet. This could be because the new version is not yet release, but in brief it is going to be. You could install directly from the repo if some problem pops up. Just paste the following in a code cell:

! pip install git+https://github.com/ale94mleon/MolDrug.git@main

Table of content

1. Execute moldrug from the command line.
2. Create your own cost function.

Execute moldrug from the command line¶

MolDrug it is meant to be for using as a python module. However it has capabilities to work from the command line. Due to the complexity of the input for a normal simulation, yaml structured file are used for the input specification. Let first see the help of moldrug.

As it shows, only one positional arguments is needed yaml_file. Then you could give:

• fitness: A customized fitness python module where your cost function must be implemented.

• outdir: The out directory.

Lets go steep by steep

yaml_file anatomy¶

Now we have all that we need. The SMILES, the definition of the box, the CReM data base and the receptor. But first let see what are the arguments that the function utils.Local needs:

We will pass all of this parameters using config.yaml file.

Note: Here we will not work neither modify with the desirability parameter of the cost function. We will talk about it in a separate tutorial inside of Advance Topics.

Warning!: Check your correct temporal file and change it accordantly. Of course you could also set a normal directory and safe the results of the simulation.

Here is an example yaml file to work with the class utils.Local.

main:
  type: Local
  njobs: 1
  pick: 1
  AddHs: False
  seed_mol: COC(=O)C=1C=CC(=CC1)S(=O)(=O)N
  costfunc: Cost
  costfunc_kwargs:
    vina_executable: vina
    receptor_path: /tmp/your_tmp_dir/x0161.pdbqt
    boxcenter:
      - 12.11
      - 1.84
      - 23.56
    boxsize:
      - 22.5
      - 22.5
      - 22.5
    exhaustiveness: 4
    ncores: 12
    num_modes: 1
  crem_db_path: /tmp/your_tmp_dir/crem.db
  grow_crem_kwargs:
    radius: 3
    min_atoms: 1
    max_atoms: 2
    ncores: 12

The structure of a yaml file is like a python dictionary but more human readable (see this youtube video if you are not familiar).

First we have the directrix:

main:

This is just the name of your main project and could be any word (we will see that moldrug.utils.GA accept also follow projects). Then we have:

type: Local
  njobs: 1
  pick: 1
  AddHs: False

Look how this section is inside of main (one indentation level). type is the name of the class inside moldrug.utils; for now could be GA or Local (case insensitive). In this case we want to use the class Local. njobs, pick and AddHs are the parameters when the class Local (or GA) is call. The rest of the parameters are just the parameters that needs Local for the initialization. All of them are at the same level of the previous ones. Because costfunc_kwargs is a dictionary; we add for its value a new indentation level:

costfunc_kwargs:
    vina_executable: vina
    receptor_path: /tmp/your_tmp_dir/x0161.pdbqt
    boxcenter:
      - 12.11
      - 1.84
      - 23.56
    boxsize:
      - 22.5
      - 22.5
      - 22.5
    exhaustiveness: 4
    ncores: 12
    num_modes: 1

As you could imagine the keyword boxcenter is a list and this is one of the way to represent list in yaml. Cool right?

The next dictionary is the same that the yaml file and we will save it in the temporal file that we created.

Move to the directory and run moldrug¶

As you see the simulation run successfully and the following files were generated:

1. local_result.pbz2; the compresses version of the Local class with all the information of the simulation (binary file).

2. local_pop.sdf; the mol structures. This file is useful to use inside Pymol.

Specify a custom cost function¶

The following is a simple cost function based on the vina score and very similar to moldrug.fitness.Cost used by MolDrug. The details on the parameters and how to correctly implement it will be discussed letter on the following section.

def MyAwesomeCost(Individual:utils.Individual, wd:str = '.vina_jobs', vina_executable:str = 'vina', receptor_path:str = None, boxcenter:List[float] = None, boxsize:List[float] =None, exhaustiveness:int = 8, ncores:int = 1,  num_modes:int = 1):

    # Getting Vina score
    cmd = f"{vina_executable} --receptor {receptor_path} --ligand {os.path.join(wd, f'{Individual.idx}.pdbqt')} "\
        f"--center_x {boxcenter[0]} --center_y {boxcenter[1]} --center_z {boxcenter[2]} "\
        f"--size_x {boxsize[0]} --size_y {boxsize[1]} --size_z {boxsize[2]} "\
        f"--out {os.path.join(wd, f'{Individual.idx}_out.pdbqt')} --cpu {ncores} --exhaustiveness {exhaustiveness} --num_modes {num_modes}"

    # Creating the ligand pdbqt
    with open(os.path.join(wd, f'{Individual.idx}.pdbqt'), 'w') as l:
        l.write(Individual.pdbqt)
    utils.run(cmd)

    # Getting the information
    best_energy = utils.VINA_OUT(os.path.join(wd, f'{Individual.idx}_out.pdbqt')).BestEnergy()
    # Changing the 3D conformation by the conformation of the binding pose
    Individual.pdbqt = ''.join(best_energy.chunk)

    # Getting the Scoring function of Vina
    Individual.vina_score = best_energy.freeEnergy

    # Getting the Scoring function of Vina and assign it to the cost attribute of the Individual
    Individual.cost = best_energy.freeEnergy
    return Individual

Because this cost function accept the same parameters we can use the same configuration yaml file; but the name of the function is different. So we have to modify that parameter. Let's save the function in a .py file with the corresponded import statements

Now we only need to call MolDrug and specify the new fitness module to use. It is recommendable to output in on folder the results; therefore also use the option outdir of the command line.

As you can see the simulation run successfully! But be aware of a possible issue! Because the cost function used is not actually inside of MolDrug. If we would like to use on the future the binary file local_result.pbz2. We have first say to Python where to find the used fitness function. For reproducibility and also peace of mind, MolDrug will generate a new file: CustomMolDrugFitness.py this is nothing more than a copy of your implemented python fitness module. This module is the one used for the internal calculation of MolDrug. Therefore, we have to say to python where is CustomMolDrugFitness and append to sys.path in order to open the binary file. If we do not do that we will get an error. Let see the example.

You could see the difference between the two cost functions and the error printed in the try-except block. Finally let's take a look to the sys.path; check that at the end it is our path with our new fitness module.

Create your own cost function¶

Create your own cost function should be straightforward as soon as we understand how MolDrug works with it. So, let's take a look again to the function defined previously:

MyAwesomeCost function it will only:

1. Perform Docking
2. Update:
   • pdbqt attribute of Individual
   • cost attribute of Individual (the must important in order to MolDrug works properly)
3. Assign vina_score attribute to Individual

Basically the idea is give a number to the attribute cost of the passed Individual instance. MolDrug works based on the class moldrug.utils.Individual; this is just a representation of each chemical structure. Therefore all the custom cost function must work based on this class; in oder words, "the chemical structure input" must be coded as an Individual and it will be the first positional argument to pass. The rest is completely optional. If the function perform I/O operations you can also provided as keyword wd. By default moldrug.utils.GA and moldrug.utils.Local create a temporal directory (in /tmp) even when you specify some values for it. This could be a problem if your function create big files (bigger than the capacity of your /tmp directory) or if you would like to preserve some of this files. In this case just create an alternative keyword (like wd_dont_touch) and work based on it. The other important idea to have in mind is that moldrug.utils.GA perform minimization (the best element is the one with lower cost attribute).

WARNING: MolDrug uses multiprocessing for the parallelization. This module uses itself pickle. In our case MolDrug parallelize the cost function; therefore the COST FUNCTION MUST BE PICKLEABLE. To see what can be pickled and unpickled go here.

Common error

multiprocessing.pool.MaybeEncodingError: Error sending result: '<multiprocessing.pool.ExceptionWithTraceback object at 0x7efd82c72b80>'. Reason: 'PicklingError("Can't pickle <class 'Boost.Python.ArgumentError'>: import of module 'Boost.Python' failed")'

The last error could be due to some user defined functions (or classes) that uses your custom cost function that are not a the top level of the module. In other words, you are using user defined functions (or classes) that are outside of your custom_fitness_module.py (in case command line is used) or your main script where the cost function was defined. Solution: put all your code for the cost function in only one .py ((in case command line is used)) file or in your main script.

The following is probably the easiest and dummy cost function possible. But I bring it here to show that there is not limitation on how to create your own cost function.

moldrug.utils.GA exports an sdf file that contains the individual of the current generations every save_pop_every_gen. For that, it uses the information contained in the pdbqt or pdbqts attributes (depending if single or multi-receptor is used respectively). Therefore if we use docking we must update the pdbqt attribute with the pose obtained by vina. If we don't do that, we will just get a conformation generated by RDKit that is automatically created when a instance of utils.Individual is made it.