Dev benchmark

README.md

@@ -62,8 +62,8 @@ First, prepare a samplesheet with your input data that looks as follows:
 `samplesheet.csv`:
 
 ```csv
-id,binder_name,starting_pdb,chains,target_hotspot_residues,min_length,max_length,number_of_final_designs,settings_advanced,settings_filters
-demo,PDL1,PDL1.pdb,A,"56",65,150,10,default_4stage_multimer.json,default_filters.json
+id,binder_name,starting_pdb,chains,target_hotspot_residues,min_length,max_length,number_of_final_designs,max_trajectories,settings_advanced,settings_filters
+demo,PDL1,PDL1.pdb,A,"56",65,150,10,100,default_4stage_multimer.json,default_filters.json
 ```
 
 Each row represents a single design instance. Detailed documentation describing job parameters can be found in the BindCraft [documentation](https://github.com/martinpacesa/BindCraft).
@@ -76,6 +76,7 @@ Briefly:
 - **min_length** defines the minimum length of the designed binder.
 - **max_length** defines the maximum length of the designed binder.
 - **number_of_final_designs** defines the number of binders required to pass QC criteria before the job is complete.
+- **max_trajectories** defines the maximum number of design trajectories (if provided in **settings_advanced** as well, it will be overridden).
 - **settings_advanced** defines advanced BindCraft settings (JSON format).
 - **settings_filters** defines advanced BindCraft filter settings (JSON format).
 

assets/samplesheet.csv

@@ -1,2 +1,2 @@
-id,binder_name,starting_pdb,chains,target_hotspot_residues,min_length,max_length,number_of_final_designs,settings_filters,settings_advanced
-s1,PDL1,https://raw.githubusercontent.com/Australian-Structural-Biology-Computing/bindflow/refs/heads/main/assets/bindcraft/PDL1.pdb,A,56,65,150,4,https://raw.githubusercontent.com/Australian-Structural-Biology-Computing/bindflow/refs/heads/main/assets/bindcraft/default_filters.json,https://raw.githubusercontent.com/Australian-Structural-Biology-Computing/bindflow/refs/heads/main/assets/bindcraft/default_4stage_multimer.json
+id,binder_name,starting_pdb,chains,target_hotspot_residues,min_length,max_length,number_of_final_designs,max_trajectories,settings_filters,settings_advanced
+s1,PDL1,https://raw.githubusercontent.com/Australian-Structural-Biology-Computing/bindflow/refs/heads/main/assets/bindcraft/PDL1.pdb,A,56,65,150,2,4,https://raw.githubusercontent.com/Australian-Structural-Biology-Computing/bindflow/refs/heads/main/assets/bindcraft/default_filters.json,https://raw.githubusercontent.com/Australian-Structural-Biology-Computing/bindflow/refs/heads/main/assets/bindcraft/default_4stage_multimer.json

assets/schema_input.json

@@ -56,6 +56,13 @@
                 "errorMessage": "",
                 "meta": ["number_of_final_designs"]
             },
+            "max_trajectories": {
+                "type": "integer",
+                "minimum": 1,
+                "default": 1000,
+                "errorMessage": "",
+                "meta": ["max_trajectories"]
+            },
             "settings_filters": {
                 "type": "string",
                 "format": "file-path",
@@ -70,6 +77,6 @@
             }
         },
         "required": ["id", "binder_name", "starting_pdb", "chains", "target_hotspot_residues",
-                    "min_length", "max_length", "number_of_final_designs"]
+                    "min_length", "max_length", "number_of_final_designs", "max_trajectories"]
     }
 }

docs/usage.md

@@ -45,8 +45,8 @@ You will need to create a samplesheet with information about the binder design j
 A final samplesheet file will look something like the one below. This is for a single design job aiming to produce 10 binders that pass QC between lengths 65-150, targeting the PDL1 protein at hotspot residue 56.
 
 ```csv title="samplesheet.csv"
-id,binder_name,starting_pdb,chains,target_hotspot_residues,min_length,max_length,number_of_final_designs,settings_advanced,settings_filters
-demo,PDL1,PDL1.pdb,A,"56",65,150,10,default_4stage_multimer.json,default_filters.json
+id,binder_name,starting_pdb,chains,target_hotspot_residues,min_length,max_length,number_of_final_designs,max_trajectories,settings_advanced,settings_filters
+demo,PDL1,PDL1.pdb,A,"56",65,150,10,100,default_4stage_multimer.json,default_filters.json
 ```
 
 | Column    | Description                                                                                                                                                                            |
@@ -59,6 +59,7 @@ demo,PDL1,PDL1.pdb,A,"56",65,150,10,default_4stage_multimer.json,default_filters
 | `min_length` | The minimum length to sample from potential binders.                                                             |
 | `max_length` | The maximum length to sample for potential binders.                                                             |
 | `number_of_final_designs` | The number of designs that pass all QC criteria before the pipeline will complete.                 |
+| `max_trajectories` | The maximum number of design trajectories (if provided in `settings_advanced` as well, it will be overridden).                          |
 | `settings_advanced` | Full path to json file defining bindcraft advanced settings (Optional: Default is `default_4stage_multimer.json`)                      |
 | `settings_filters` | Full path to json file defining bindcraft QC filter settings (Optional: Default is `default_filters.json`)                  |
 

modules/local/bindcraft/main.nf

@@ -10,7 +10,7 @@ process BINDCRAFT {
     output:
         tuple val(meta), path("*_final_design_stats.csv"), emit: stats
         tuple val(meta), path("*_output/Accepted/Ranked"), emit: accepted_ranked
-        tuple val(meta), path("*_output/Accepted/*pdb"), emit: accepted
+        tuple val(meta), path("*_output/Accepted/*pdb"), emit: accepted,  optional: true
         tuple val(meta), path("*_output"), emit: output_dir
         tuple val(meta), path("*_output/failure_csv.csv") , emit: failure_csv
         tuple val(meta), path("*_output/final_design_stats.csv"), emit: final_design_stats
@@ -26,7 +26,7 @@ process BINDCRAFT {
     def args = task.ext.args ?: ''
 
     """
-    python /work/FreeBindCraft/bindcraft.py \\
+    python -u /work/FreeBindCraft/bindcraft.py \\
         --settings ${target_file} \\
         --filters ${filters} \\
         --advanced ${advanced_settings} \\

modules/local/jsonmanager/main.nf

@@ -11,21 +11,35 @@ process JSONMANAGER {
     val batches
 
     output:
-    path "*.json"      , emit: json
+    path "target_json/*.json"  , emit: json
+    path "advanced_json/*.json", emit: advanced_json
     path "versions.yml", emit: versions
 
     when:
     task.ext.when == null || task.ext.when
 
     script:
     def quote_char = task.ext.args2 ?: "\""
+    def default_advanced = params.settings_advanced ?: "${projectDir}/assets/bindcraft/default_4stage_multimer.json"
     """
     #!/usr/bin/env python3
     import csv
     import json
     import os, sys
     import math
 
+    sample_dir = os.path.dirname(os.path.abspath("${samplesheet}"))
+    os.makedirs("target_json", exist_ok=True)
+    os.makedirs("advanced_json", exist_ok=True)
+
+    def resolve_path(raw_path):
+        if not raw_path:
+            return raw_path
+        if os.path.isabs(raw_path):
+            return raw_path
+        candidate = os.path.join(sample_dir, raw_path)
+        return candidate if os.path.exists(candidate) else raw_path
+
     with open("${samplesheet}", 'r') as csvfile:
         reader = csv.DictReader(csvfile, quotechar="\\${quote_char}")
         for row in reader:
@@ -38,26 +52,46 @@ process JSONMANAGER {
                 row['starting_pdb'] = os.path.basename(row['starting_pdb'])
             row['lengths'] = [int(row['min_length']), int(row['max_length'])]
             number_of_final_designs = int(row['number_of_final_designs'])
-            if batches > number_of_final_designs:
-                batches = number_of_final_designs
+            max_trajectories = int(row['max_trajectories'])
+            if batches > max_trajectories:
+                batches = max_trajectories
+
+            settings_advanced = (row.get('settings_advanced') or '').strip()
+            if not settings_advanced:
+                settings_advanced = "${default_advanced}"
+            settings_advanced = resolve_path(settings_advanced)
+            with open(settings_advanced, 'r') as advanced_file:
+                advanced_template = json.load(advanced_file)
 
-            batch_size = math.ceil(number_of_final_designs / batches)
-            row['number_of_final_designs'] = batch_size
+            max_trajectories_per_batch = math.ceil(max_trajectories / batches)
+            number_of_final_designs_per_batch = math.ceil(number_of_final_designs / batches)
+            row['number_of_final_designs'] = number_of_final_designs_per_batch
+            row.pop('max_trajectories', None)
             for batch_id in range(0, batches):
                 if batch_id == batches - 1:
-                    row['number_of_final_designs'] = number_of_final_designs - batch_id * batch_size
+                    row['number_of_final_designs'] = number_of_final_designs - batch_id * number_of_final_designs_per_batch
                 row['design_path'] = f"{sample_id}_{batch_id}_output"
-                with open(f"{sample_id}-{batch_id}.json", 'w') as jsonfile:
+                with open(f"target_json/{sample_id}-{batch_id}.json", 'w') as jsonfile:
                     json.dump(row, jsonfile, indent=2)
+
+                advanced_settings = dict(advanced_template)
+                advanced_settings['max_trajectories'] = max_trajectories_per_batch
+                if batch_id == batches - 1:
+                    advanced_settings['max_trajectories'] = max_trajectories - batch_id * max_trajectories_per_batch
+                with open(f"advanced_json/{sample_id}-{batch_id}-advanced.json", 'w') as advanced_json_file:
+                    json.dump(advanced_settings, advanced_json_file, indent=2)
 
     with open ("versions.yml", "w") as version_file:
 	    version_file.write("\\"${task.process}\\":\\n    python: {}\\n".format(sys.version.split()[0].strip()))
     """
 
     stub:
     """
-    touch s1-0.json
-    touch s1-1.json
+    mkdir -p target_json advanced_json
+    touch target_json/s1-0.json
+    touch target_json/s1-1.json
+    touch advanced_json/s1-0-advanced.json
+    touch advanced_json/s1-1-advanced.json
 
     cat <<-END_VERSIONS > versions.yml
     "${task.process}":

subworkflows/local/run_bindcraft.nf

@@ -36,8 +36,7 @@ workflow RUN_BINDCRAFT {
        [
             ["id": row.id], 
             file(row.starting_pdb, checkIfExists: true),
-            get_file(row.settings_filters, params.settings_filters, "${projectDir}/assets/bindcraft/default_filters.json"),
-            get_file(row.settings_advanced, params.settings_advanced, "${projectDir}/assets/bindcraft/default_4stage_multimer.json"),    
+            get_file(row.settings_filters, params.settings_filters, "${projectDir}/assets/bindcraft/default_filters.json")
        ]
     }
     .set { ch_settings }
@@ -57,10 +56,19 @@ workflow RUN_BINDCRAFT {
 
     JSONMANAGER.out.json
         .flatten()
-        .map{[["id": it.baseName.split('-')[0..-2].join('-')], it.baseName.split('-')[-1].replace(".json", ""), it]}
+        .map{[["id": it.baseName.split('-')[0..-2].join('-'), "batch": it.baseName.split('-')[-1].replace(".json", "")], it]}
+        .join(
+            JSONMANAGER.out.advanced_json
+                .flatten()
+                .map{
+                    def stripped_base = it.baseName.replaceFirst(/-advanced$/, '')
+                    [["id": stripped_base.split('-')[0..-2].join('-'), "batch": stripped_base.split('-')[-1]], it]
+                }
+        )
+        .map{[["id": it[0].id], it[0].batch, it[1], it[2]]}
         .combine(ch_settings)
-        .filter{it[0].id == it[3].id}
-        .map{[["id": it[0].id, "batch": it[1]], it[2], it[4], it[5], it[6]]}
+        .filter{it[0].id == it[4].id}
+        .map{[["id": it[0].id, "batch": it[1]], it[2], it[5], it[6], it[3]]}
         .set {ch_bindcraft_input}
 
     BINDCRAFT(
@@ -73,7 +81,8 @@ workflow RUN_BINDCRAFT {
     BINDCRAFT.out.accepted
         .map{[["id": it[0].id], it[1]]}
         .groupTuple()
-        .join(ch_final_designs)
+        .join(ch_final_designs, remainder: true)
+        .map { [it[0], it[1] ?: [], it[2]] }
         .subscribe{
             if (it[1].size() < it[2]){
                 log.warn "Sample: ${it[0].id}: The pipeline was unable to generate the target number of successful designs (${it[1].size()} of ${it[2]}) in the allocated time. Please consider changing hotspot residues or design configuration to increase design success rates"
@@ -137,4 +146,4 @@ def get_file(String sheet_file, String general_file, String assets_file) {
     }
 
     return file(sheet_file)
-}
+}