Global Matching Sphere Optimization: Difference between revisions

From DISI
Jump to navigation Jump to search
No edit summary
mNo edit summary
 
(5 intermediate revisions by the same user not shown)
Line 1: Line 1:
== Goal ==
To optimize your matching sphere (MS) setups getting faster docking and more high-scoring ligands with fewer spheres.


== Goal ==
To optimize your matching sphere (MS) setups getting more enrichment with fewer spheres.
== Description ==
== Description ==
The program performs optimization of matching spheres using genetic algorithm. It selects spheres from two sets:
The program performs optimization of matching spheres by pruning and stochastic optimization. It selects spheres from two sets:
* heavy atoms of xtal-lig  
* heavy atoms of xtal-lig  
* spheres prepared by SPHGEN program
* spheres prepared by SPHGEN program
At each generation, N matching sphere sets are created, containing a maximum of M spheres each. Then retrospective docking is done for each set, and sets are ranked by the


* enrichment (normalized logAUC, see [http://arxiv.org/abs/2210.10905 Ian's paper]),  
Juggler generates an initial MS set consisting of 100 spheres (maximum in DOCK 3.8). This set is used for retrospective docking, and then KDTree algorithm is used to prune the set to the required number of spheres by discarding all spheres that were not used in generation of the poses of the known binders ("actives"). This procedure is repeated to account for any differences in matching produced by reducing the MS set.
* RMSD of the docked pose to the experimental one.  


After that, a quarter of sets "survive" and produce a new generation by direct transfer, mutations, and crossover. This process is repeated until enrichment, RMSD and minimum number of spheres do not change substantially in 10 generations.
After this, the resulting set is transferred to the stepwise optimization procedure which conducts random perturbations of the sphere sets. Retrospective docking is done for each set, and sets are ranked by the
 
* enrichment (normalized logAUC, see [http://arxiv.org/abs/2210.10905 Ian's paper]),
* the average score of the top 1% of ligands.


The program consists of two main modules:
The program consists of two main modules:
* a Python script (<code>juggler.py</code>) that performs MS generation, optimization, and ranking.
* a Python script (<code>juggler.py</code>) that performs MS generation, optimization, and ranking.
* a Bash script, that watches created directory structure, runs docking and processes docking results
* a Bash script (<code>rundockd.sh</code>), that watches created directory structure, runs docking and processes docking results.
 
== Setup & Running ==
== Setup & Running ==
So far, the program is running on Wynton. LMK if you are interested in launching it on Gimel or other clusters.


The scripts and example config file are in <code>/wynton/group/bks/soft/juggler</code>
=== Setup ===
Dependencies:
* python 3.8.1 or higher
* [https://dock.docking.org/DOCK3.8/ DOCK3.8]
* [https://github.com/docking-org/pydock3 pydock3]
* [https://github.com/docking-org/SUBDOCK subdock]
* rdkit
* pandas
* numpy
* yaml


=== Preparation ===
=== Preparation ===
Prepare <code>dockfiles</code> directory with any tools of your liking (blastermaster, dockopt etc). You will also need <code>rec.pdb</code>, <code>rec.crg.pdb</code>, <code>xtal-lig.pdb</code>, <code>ligands.names</code>, <code>decoys.names</code> and a <code>sdi</code> directory with the paths to ligand .tgz files. To get RMSD of xtal-lig docked poses to the experimental pose, your xtal-lig.pdb must have correct bond orders and atom valences. You can edit it in Schrodinger and save as <code>xtal-lig.pdb</code>.
What you need to prepare:


Prepare <code>juggler_config.yml</code> file. Put it into an empty directory.
* <code>dockfiles</code> directory with any tools of your liking (blastermaster, dockopt etc).
* <code>rec.crg.pdb</code>
* <code>xtal-lig.pdb</code>: To get RMSD of xtal-lig docked poses to the experimental pose, your xtal-lig.pdb must have correct bond orders and atom valences. You can edit it in Schrodinger and save as <code>xtal-lig.pdb</code>.
* <code>ligands.names</code>
* <code>decoys.names</code>
* <code>sdi</code> file with the paths to ligand <code>.tgz</code> files.


=== Running ===
Prepare <code>juggler_config.yml</code> file. Put the config into an empty directory.<syntaxhighlight lang="yaml">
Enter a screen environment so your run is not interrupted if you disconnect your SSH session. Then do:
################################################
# Paths for your target
receptor_file_path: "/test/rec.crg.pdb"
xtal_lig_file_path: "/test/xtal-lig.pdb"
dock_files_dir_path: "/test/dockfiles"
lig_names_file_path: "/test/ligands.names"
dec_names_file_path: "/test/decoys.names"
sdi_file_path: "test/ligands_sdi"
 
################################################
# Executables and running
dockbase: "/path/to/DOCK"
dock64_bin: "path/to/dock64"
subdock_bash_file_path: "/path/to/subdock.bash"
queue_type: "sge" # "slurm" or "sge"
 
###############################################
# Max and min number of spheres
min_sph: 4 # min is 4
max_sph: 10 # max is 100
</syntaxhighlight>


<code>source /wynton/group/bks/soft/python_envs/python3.8.5.sh</code>
The <code>dock64_bin</code> parameter is optional; if absent, <code>{dockbase}/docking/DOCK/bin/dock64</code> will be used.


<code>python /wynton/group/bks/soft/juggler/juggler.py</code>
=== Running ===
You can either:
* Enter a screen environment so your run is not interrupted if you disconnect your SSH session, or
* Run Juggler using a queuing system. See example files for the slurm and sge below.


You can detach from the screen (Ctrl-A d).
In both cases you need to launch Juggler and the docking daemon simultaneously.


Open a new screen. In the same directory launch a docking daemon
==== In a screen ====
<syntaxhighlight lang="bash">
source /path/to/python/env
# or
conda activate pydock3
# rundockd should run in the background to manage docking jobs
sh rundockd.sh 2>&1 > rundockd.log &
python juggler.py 2>&1 > juggler.log
</syntaxhighlight>


<code>/wynton/group/bks/soft/juggler/rundockd-wynton-taskid.sh</code>
==== Via a queue ====


You can run other calculations on Wynton in the meantime, as Juggler will track the task IDs that it launched.
===== SGE =====
<syntaxhighlight lang="bash">
#! /bin/bash
#$ -cwd
#$ -q long.q
#$ -o stdout_juggler
#$ -e stdout_juggler
#$ -l s_rt=72:58:00
#$ -l h_rt=73:00:00
#$ -l mem_free=10G
#$ -pe smp 2
source /path/to/pydock3/env.sh
# or conda activate pydock3
sh /path/to/juggler/rundockd.sh 2>&1 > rundockd.log &
python /path/to/juggler/juggler.py 2>&1 > juggler.log
</syntaxhighlight>


The run will take few hours to ~2 days depending on the number of actives and decoys and the load of Wynton.
===== SLURM =====
<syntaxhighlight lang="bash">
#! /bin/bash
#$ -cwd
#$ -q long.q
#$ -o stdout_juggler
#$ -e stdout_juggler
#$ -l s_rt=23:58:00
#$ -l h_rt=24:00:00
#$ -l mem_free=10G
source /path/to/pydock3/env.sh
# or conda activate pydock3
sh /path/to/juggler/rundockd.sh 2>&1 > rundockd.log &
python /path/to/juggler/juggler.py 2>&1 > juggler.log
</syntaxhighlight>


=== Processing results ===
=== Processing results ===
The script will print the paths to where three best matching sphere sets are:
At the end of a run you will get a message that convergence was reached. You will see the directory <code>best_set</code> that contains <code>dockfiles</code> and docking results for the best matching sphere set found. This directory is updated at each step, so if the run fails or convergence is not reached, you can still access the optimal set.
 
Other output files:
* <code>stepwise_opt_best_sets.dat</code> — lists the IDs and the nlogAUC values for the best set in each stepwise optimization round.
* <code>stepwise_opt_metrics.dat</code> — lists IDs, nlogAUC, RMSD and average scores for the top 1% ligands for all sets tested during the stepwise optimization.
* optional: <code>juggler.log</code> or <code>stdout</code> — contains the log of the run.


* best enrichment
=== For BKS lab users ===
* best RMSD
==== Gimel ====
* best balanced metrics (highest enrichment, lowest RMSD and lowest Nsph).
Juggler is in <code>/mnt/nfs/exa/work/ak87/UCSF/JUGGLER/SCRIPTS/JUGGLER</code>


You can use <code>dockfiles</code> from the listed directories.
Subdock is in <code>/mnt/nfs/exa/work/ak87/PROGRAM/SUBDOCK/SUBDOCK</code>


You can also track the optimization progress running the following script in your working directory:
You can use this file to submit to SLURM queue


<code>/wynton/group/bks/soft/juggler/plot_all_metrics.py</code>
<syntaxhighlight lang="bash">
#! /bin/bash
#$ -cwd
#$ -q long.q
#$ -o stdout_juggler
#$ -e stdout_juggler
#$ -l s_rt=23:58:00
#$ -l h_rt=24:00:00
#$ -l mem_free=10G
source /nfs/soft/ian/python3.8.5.sh
sh /mnt/nfs/exa/work/ak87/UCSF/JUGGLER/SCRIPTS/JUGGLER/rundockd.sh 2>&1 > rundockd.log & #/dev/null &
python /mnt/nfs/exa/work/ak87/UCSF/JUGGLER/SCRIPTS/JUGGLER/juggler.py 2>&1 > orbebb.log
</syntaxhighlight>


It produces <code>combined_metrics.png</code>
==== Wynton ====
[[File:Combined metrics plot from a GA run.png|thumb|Combined metrics plot from a GA run]]


'''If not converged'''
Juggler is in <code>/wynton/group/bks/work/ak87/UCSF/JUGGLER/SCRIPTS/JUGGLER</code>


The program will stop after 200 generations if convergence is not reached. In case it takes too long you can stop it any time by pressing Ctrl-C. It doesn't mean that you have no results, though. Juggler generates <code>combined_metrics.dat</code> file in the working directory, which contains metrics for all sets explored. It contains the following columns:
Subdock is in <code>/wynton/group/bks/work/ak87/UCSF/JUGGLER/SCRIPTS/SUBDOCK</code>


Generations Set# NormLogAUC RMSD Nsph Combined_metrics
You can use this file to submit to SLURM queue


You can paste its content into Excel, sort by the highest NormLogAUC and pick a MS set of your liking.
<syntaxhighlight lang="bash">
#! /bin/bash
#$ -cwd
#$ -q long.q
#$ -o stdout_juggler
#$ -e stdout_juggler
#$ -l s_rt=72:58:00
#$ -l h_rt=73:00:00
#$ -l mem_free=10G
#$ -pe smp 2
source /wynton/group/bks/soft/python_envs/env.sh
sh /wynton/group/bks/work/ak87/UCSF/JUGGLER/SCRIPTS/JUGGLER/rundockd.sh 2>&1 > rundockd.log & #/dev/null &
python /wynton/group/bks/work/ak87/UCSF/JUGGLER/SCRIPTS/JUGGLER/juggler.py 2>&1 > orbebb.log
</syntaxhighlight>

Latest revision as of 22:18, 17 April 2026

Goal

To optimize your matching sphere (MS) setups getting faster docking and more high-scoring ligands with fewer spheres.

Description

The program performs optimization of matching spheres by pruning and stochastic optimization. It selects spheres from two sets:

  • heavy atoms of xtal-lig
  • spheres prepared by SPHGEN program

Juggler generates an initial MS set consisting of 100 spheres (maximum in DOCK 3.8). This set is used for retrospective docking, and then KDTree algorithm is used to prune the set to the required number of spheres by discarding all spheres that were not used in generation of the poses of the known binders ("actives"). This procedure is repeated to account for any differences in matching produced by reducing the MS set.

After this, the resulting set is transferred to the stepwise optimization procedure which conducts random perturbations of the sphere sets. Retrospective docking is done for each set, and sets are ranked by the

  • enrichment (normalized logAUC, see Ian's paper),
  • the average score of the top 1% of ligands.

The program consists of two main modules:

  • a Python script (juggler.py) that performs MS generation, optimization, and ranking.
  • a Bash script (rundockd.sh), that watches created directory structure, runs docking and processes docking results.

Setup & Running

Setup

Dependencies:

Preparation

What you need to prepare:

  • dockfiles directory with any tools of your liking (blastermaster, dockopt etc).
  • rec.crg.pdb
  • xtal-lig.pdb: To get RMSD of xtal-lig docked poses to the experimental pose, your xtal-lig.pdb must have correct bond orders and atom valences. You can edit it in Schrodinger and save as xtal-lig.pdb.
  • ligands.names
  • decoys.names
  • sdi file with the paths to ligand .tgz files.

Prepare juggler_config.yml file. Put the config into an empty directory.

################################################
# Paths for your target
receptor_file_path: "/test/rec.crg.pdb"
xtal_lig_file_path: "/test/xtal-lig.pdb"
dock_files_dir_path: "/test/dockfiles"
lig_names_file_path: "/test/ligands.names"
dec_names_file_path: "/test/decoys.names"
sdi_file_path: "test/ligands_sdi"

################################################
# Executables and running
dockbase: "/path/to/DOCK"
dock64_bin: "path/to/dock64"
subdock_bash_file_path: "/path/to/subdock.bash"
queue_type: "sge" # "slurm" or "sge"

###############################################
# Max and min number of spheres
min_sph: 4 # min is 4
max_sph: 10 # max is 100

The dock64_bin parameter is optional; if absent, {dockbase}/docking/DOCK/bin/dock64 will be used.

Running

You can either:

  • Enter a screen environment so your run is not interrupted if you disconnect your SSH session, or
  • Run Juggler using a queuing system. See example files for the slurm and sge below.

In both cases you need to launch Juggler and the docking daemon simultaneously.

In a screen

source /path/to/python/env
# or
conda activate pydock3
# rundockd should run in the background to manage docking jobs
sh rundockd.sh 2>&1 > rundockd.log &
python juggler.py 2>&1 > juggler.log

Via a queue

SGE
#! /bin/bash
#$ -cwd
#$ -q long.q
#$ -o stdout_juggler
#$ -e stdout_juggler
#$ -l s_rt=72:58:00
#$ -l h_rt=73:00:00
#$ -l mem_free=10G
#$ -pe smp 2
source /path/to/pydock3/env.sh
# or conda activate pydock3
sh /path/to/juggler/rundockd.sh 2>&1 > rundockd.log &
python /path/to/juggler/juggler.py 2>&1 > juggler.log
SLURM
#! /bin/bash
#$ -cwd
#$ -q long.q
#$ -o stdout_juggler
#$ -e stdout_juggler
#$ -l s_rt=23:58:00
#$ -l h_rt=24:00:00
#$ -l mem_free=10G
source /path/to/pydock3/env.sh
# or conda activate pydock3
sh /path/to/juggler/rundockd.sh 2>&1 > rundockd.log &
python /path/to/juggler/juggler.py 2>&1 > juggler.log

Processing results

At the end of a run you will get a message that convergence was reached. You will see the directory best_set that contains dockfiles and docking results for the best matching sphere set found. This directory is updated at each step, so if the run fails or convergence is not reached, you can still access the optimal set.

Other output files:

  • stepwise_opt_best_sets.dat — lists the IDs and the nlogAUC values for the best set in each stepwise optimization round.
  • stepwise_opt_metrics.dat — lists IDs, nlogAUC, RMSD and average scores for the top 1% ligands for all sets tested during the stepwise optimization.
  • optional: juggler.log or stdout — contains the log of the run.

For BKS lab users

Gimel

Juggler is in /mnt/nfs/exa/work/ak87/UCSF/JUGGLER/SCRIPTS/JUGGLER

Subdock is in /mnt/nfs/exa/work/ak87/PROGRAM/SUBDOCK/SUBDOCK

You can use this file to submit to SLURM queue

#! /bin/bash
#$ -cwd
#$ -q long.q
#$ -o stdout_juggler
#$ -e stdout_juggler
#$ -l s_rt=23:58:00
#$ -l h_rt=24:00:00
#$ -l mem_free=10G
source /nfs/soft/ian/python3.8.5.sh
sh /mnt/nfs/exa/work/ak87/UCSF/JUGGLER/SCRIPTS/JUGGLER/rundockd.sh 2>&1 > rundockd.log & #/dev/null &
python /mnt/nfs/exa/work/ak87/UCSF/JUGGLER/SCRIPTS/JUGGLER/juggler.py 2>&1 > orbebb.log

Wynton

Juggler is in /wynton/group/bks/work/ak87/UCSF/JUGGLER/SCRIPTS/JUGGLER

Subdock is in /wynton/group/bks/work/ak87/UCSF/JUGGLER/SCRIPTS/SUBDOCK

You can use this file to submit to SLURM queue

#! /bin/bash
#$ -cwd
#$ -q long.q
#$ -o stdout_juggler
#$ -e stdout_juggler
#$ -l s_rt=72:58:00
#$ -l h_rt=73:00:00
#$ -l mem_free=10G
#$ -pe smp 2
source /wynton/group/bks/soft/python_envs/env.sh
sh /wynton/group/bks/work/ak87/UCSF/JUGGLER/SCRIPTS/JUGGLER/rundockd.sh 2>&1 > rundockd.log & #/dev/null &
python /wynton/group/bks/work/ak87/UCSF/JUGGLER/SCRIPTS/JUGGLER/juggler.py 2>&1 > orbebb.log