Improve host/ligand_name/selection

[ale94mleon]

This may be confused. So it would be nice to have only one definition

• *_name: is used for:
  • The creation of the index files and use the syntax of gmx make_ndx. At the moment, it uses the name of the group for the receptor but it can be easily extended.
  • By the function bindflow.utils.tools.center_xtc. This function is called on the equi.boresch and mmbpsa.ana to fix the trajectory before the analysis
• *_selection is used at snakemake level when the Boresch restraint are selected. Here MDAnalysis is used, so the syntax is different to gmx make_ndx.

The best solution: only use one syntax either gmx make_ndx or MDAnalysis, define better naming and work with selection instead of group names. e.g.:

host_ndx_selection and host_boresch_selection

This allows more flexibility for the selection of the host for gmx_MMPBSA. One quick solution is to include the current host_selection keyword for mmpbsa calculations and create the index file at snakemake level (this solution has the problem of different syntax as well). But I would go for standardise to something like host_ndx_selection and host_boresch_selection.

For example on gmx_MMPBSA if some cofactors would like to be included in the calculation, we can use host_name=Protein or resname yCA or resname COF, and this will be internally converted to group Protein or resname COF, but the addition of the group keyword is hidden and prompt to confusion on the user end. Better that users are in total charge and aware to use a full-valid gmx make_ndx. This fix will break the workflow because of the call of center_xtc on the mmpbsa workflow that uses host_name

[ale94mleon]

Temporary Solution

To address this issue, the pipeline can be executed in two stages as follows:

Step 1: Initial Run

First, run the pipeline up to the build_ligand_system rule:

snakemake --until build_ligand_system --jobs 100000 --latency-wait 360 --cluster-cancel scancel --rerun-incomplete --keep-incomplete --keep-going --cluster '<cluster definition>'

This can be done by modifying the command in the job.sh script. Once this step is complete, use the script make_ndx.py to modify the RECEPTOR group to include ions and water cofactors.

Step 2: Modifying the Index File

The make_ndx.py script processes the index files in the pipeline and updates them appropriately. Here’s the script:

from bindflow.preparation.solvent import index_for_soluble_system
from pathlib import Path
from glob import glob

index_files = glob('mmpbsa/*/*/input/complex/index.ndx')

for index_file in index_files:
    index_file = Path(index_file)
    parent_dir = index_file.parent
    print(index_file)
    index_for_soluble_system(
        configuration_file=parent_dir/'complex.gro',
        ndxout=index_file,
        ligand_name='LIG',
        host_name='Protein or resname COF',
        load_dependencies=['export GMX_MAXBACKUP=-1']
    )

Step 3: Resuming the Pipeline

After modifying the index files, remove the --until build_ligand_system option from the command and resume execution of the pipeline:

snakemake --jobs 100000 --latency-wait 360 --cluster-cancel scancel --rerun-incomplete --keep-incomplete --keep-going --cluster '<cluster definition>'

This two-stage approach ensures that the RECEPTOR group is updated correctly before proceeding with the remaining workflow.