# Portable rDock Pipeline

## List known-good PDB complexes
```bash
python scripts/prepare_pdb_ligand_dataset.py --list-known-good
```

## Environment and rDock installation
The portable bundle no longer ships `rdock/` binaries.
The runtime is rDock-only. Removed secondary and legacy backends are not part of this portable pack.

`setup_server.sh` creates or updates the conda environment with:
- `conda-forge` + `bioconda`
- `--strict-channel-priority`
- platform-independent Python/runtime packages:
  - `python=3.11`
  - `numpy pandas scipy scikit-learn matplotlib pyyaml tqdm joblib biopython pytest requests`
  - `rdkit`
  - `openbabel`
  - `plip`

PLIP is part of the portable environment. rDock still produces the poses and native
SCORE terms; PLIP is used after docking to count biological interaction classes
such as hydrophobic contacts, hydrogen bonds, salt bridges, pi-stacking and related
contacts. If PLIP fails for a specific pose, the row records the failure and falls
back to the rDock energy/pose interaction proxy instead of silently dropping the
ligand.

For rDock itself, the setup script first tries the conda/bioconda `rdock=24.04.204_legacy`
package. That package is available on common Linux server platforms but may not exist for
every local development platform, including some Apple Silicon setups. If conda cannot
install rDock, the script falls back to an existing `rbdock`/`rbcavity` installation on
`PATH` and validates the detected `RBT_ROOT`.

On Apple Silicon, a typical local setup is:

```bash
brew install rdock
```

On Linux x86_64 servers, use `setup_server.sh`; it will install the conda package when
available through the configured channels.

Setup exports:
- `RDOCK_ROOT=$CONDA_PREFIX`
- `RBT_ROOT=$CONDA_PREFIX/share/rdock` or `$CONDA_PREFIX/share` when that layout is used by the installed package
- `RBT_HOME=$RBT_ROOT`
- `PATH=$CONDA_PREFIX/bin:$PATH`
- `LD_LIBRARY_PATH=$CONDA_PREFIX/lib:$LD_LIBRARY_PATH`

`check_environment.sh` hard-fails if any required dependency is missing, including:
- `rdkit`
- `rdkit.Chem`
- `rdkit.Chem.AllChem`
- `rdkit.Chem.Descriptors`
- `rdkit.Chem.rdMolDescriptors`
- `plip`
- `plip.structure.preparation`
- `obabel`
- `rbdock`
- `rbcavity`
- `RBT_ROOT`
- `dock.prm`

so the pipeline code does not need to be patched on the server.

## Setup on server
```bash
bash setup_server.sh portable-rdock-pipeline
conda activate portable-rdock-pipeline
bash check_environment.sh
```

## Inspect hetero ligands for a PDB
```bash
python scripts/prepare_pdb_ligand_dataset.py \
  --list-hetero \
  --pdb-id 1A52
```

For `1A52`, the parser should confirm that the correct reference ligand is `EST`, not `MK1`. `AU` is treated as an ion/metal and ignored as a docking reference ligand.

## Prepare a dataset for one PDB complex
```bash
python scripts/prepare_pdb_ligand_dataset.py \
  --pdb-id 4HG7 \
  --receptor-chain A \
  --reference-ligand-resname NUT \
  --reference-ligand-chain A \
  --n-ligands 1000 \
  --ligand-source smiles_file \
  --smiles-file data/examples/example_smiles_1000.smi \
  --out /data/datasets/rdock_4hg7_1000 \
  --force
```

## Prepare a dataset by fetching real similar SMILES from PubChem
```bash
python scripts/prepare_pdb_ligand_dataset.py \
  --pdb-id 4HG7 \
  --receptor-chain A \
  --reference-ligand-resname NUT \
  --reference-ligand-chain A \
  --n-ligands 1000 \
  --ligand-source pubchem_similarity \
  --identity-threshold-start 99 \
  --identity-threshold-stop 85 \
  --allow-partial-ligand-set \
  --out /data/datasets/rdock_4hg7_pubchem \
  --force
```

This fetches real PubChem compounds via the official PUG REST `fastsimilarity_2d` endpoint. If PubChem returns fewer than the requested count within the threshold range, the script keeps the available ligands, writes them into the dataset, and records a warning in `dataset_manifest.json`, `qc/preparation_report.md`, and `logs/pubchem_diagnostics.json`.

You can also request larger libraries, for example `25000`, when `--ligand-source pubchem_similarity` is set explicitly. In that case the script still keeps a partial real result set if PubChem does not return enough compounds within the threshold range.

## Auto-detect the reference ligand
```bash
python scripts/prepare_pdb_ligand_dataset.py \
  --pdb-id 1A52 \
  --receptor-chain A \
  --auto-reference-ligand \
  --n-ligands 1000 \
  --ligand-source smiles_file \
  --smiles-file data/examples/example_smiles_1000.smi \
  --out /data/datasets/rdock_1a52_1000 \
  --force
```

## Validate an existing dataset directory
```bash
python scripts/prepare_pdb_ligand_dataset.py \
  --validate-only \
  --out /data/datasets/rdock_4hg7_1000
```

or

```bash
python -m docking_pipeline validate-dataset \
  --dataset-dir /data/datasets/rdock_4hg7_1000
```

## Resume ligand 3D preparation from an interrupted dataset
If `ligands/all_ligands.smi` already exists and the build stopped during OpenBabel 3D conversion, resume only the ligand preparation step:

```bash
python scripts/prepare_pdb_ligand_dataset.py \
  --prepare-ligands-only \
  --out /data/datasets/rdock_4hg7_pubchem_50k \
  --ligand-batch-size 250 \
  --ligand-jobs auto \
  --ligand-cpu-fraction 0.85
```

This mode:
- uses the existing `ligands/all_ligands.smi`
- writes batch outputs under `ligands/obabel_batches/`
- updates `ligands/ligand_preparation_progress.json`
- skips completed batches on rerun
- if only a partial `ligands/all_ligands.sdf` exists and there are no batch files yet, it resumes from the prepared prefix and builds the missing suffix

Check progress with:

```bash
cat /data/datasets/rdock_4hg7_pubchem_50k/ligands/ligand_preparation_progress.json
```

## Prepare a real 1000-SMILES example file
```bash
python scripts/download_example_smiles.py \
  --n 1000 \
  --out data/examples/example_smiles_1000.smi
```

This file contains real compounds derived from the benchmark ligand library already present in this repository. It is not a synthetic or random SMILES set.

## Large libraries
For large libraries such as `50000` ligands, use a real SMILES file:

```bash
python scripts/prepare_pdb_ligand_dataset.py \
  --pdb-id 1A52 \
  --receptor-chain A \
  --auto-reference-ligand \
  --n-ligands 50000 \
  --ligand-source smiles_file \
  --smiles-file /data/libraries/real_50k_library.smi \
  --out /data/datasets/rdock_1a52_50k \
  --force
```

`PubChem similarity` may fail to collect 1000 analogs for some reference ligands. When that happens, the generator writes `logs/pubchem_diagnostics.json` and you should switch to `--ligand-source smiles_file`.

## Run multi-fidelity docking
```bash
python -m docking_pipeline benchmark-adaptive \
  --dataset-dir /data/datasets/rdock_1a52_50k \
  --strategy multifidelity_adaptive_rdock \
  --fidelity-levels 5,10,15,30,50 \
  --cost-budget-runs 250000 \
  --jobs auto \
  --cpu-fraction 0.85 \
  --out /data/results/rdock_1a52_50k_mf \
  --resume
```

## Production screening without a reference ligand
```bash
python -m docking_pipeline screen-reference-free \
  --dataset-dir /data/datasets/my_prepared_pocket_dataset \
  --strategy reference_free_active_learning_v2 \
  --triage-model classifier \
  --triage-controller auto_recall \
  --triage-retain-fraction 0.02 \
  --triage-target-recall 0.98 \
  --classifier-top-percentile 0.05 \
  --classifier-threshold-mode recall_target \
  --model-fallback-if-worse union_with_cluster_only \
  --survivor-combination-policy union \
  --max-retain-fraction-before-not-useful 0.5 \
  --calibration-size 1000 \
  --fidelity-validation-size 200 \
  --fidelity-levels 5,10,15,30,50 \
  --jobs 64 \
  --chunk-size 50 \
  --cpu-fraction 0.85 \
  --out /data/results/my_screen \
  --resume
```

This mode:
- does not require `reference_ligand.sdf`
- treats any reference ligand present in the dataset as optional benchmark-only metadata
- fails early with `RDKit_REQUIRED_FOR_REFERENCE_FREE_MODEL` if RDKit is unavailable
- writes triage artifacts:
  - `tables/initial_prefilter_decisions.csv`
  - `tables/triage_scores.csv`
  - `tables/triage_survivors.csv`
  - `tables/triage_rejected.csv`
  - `metrics/triage_metrics.json`
  - `metrics/fidelity_reliability.json`

## Compare strategies on one shared sampled reference
```bash
python -m docking_pipeline benchmark-triage-strategies \
  --dataset-dir /data/datasets/test_1500 \
  --reference-mode sampled \
  --reference-sample-size 300 \
  --reference-sample-seed 42 \
  --strategies cluster_only_triage,cheap_descriptor_filter_only,diverse_random_cost_balanced,single_fidelity_cost_balanced,reference_free_triage_bandit_v1,reference_free_active_learning_v2 \
  --triage-target-recall 0.98 \
  --triage-retain-fraction 0.02 \
  --calibration-size 300 \
  --fidelity-validation-size 100 \
  --cost-budget-runs 3000 \
  --seeds 1,2,3 \
  --jobs 10 \
  --chunk-size 20 \
  --rdock-timeout-seconds 14400 \
  --out /data/results/test_1500_strategy_comparison \
  --resume
```

## Small benchmark with full reference and the new reference-free strategy
```bash
python -m docking_pipeline benchmark-adaptive \
  --dataset-dir /data/datasets/test_1500 \
  --strategy reference_free_active_learning_v2 \
  --reference-mode full \
  --evaluation-pool-mode same_pool \
  --triage-model classifier \
  --triage-controller auto_recall \
  --triage-retain-fraction 0.05 \
  --triage-target-recall 0.98 \
  --calibration-size 300 \
  --fidelity-validation-size 100 \
  --balanced-baselines true \
  --jobs 16 \
  --chunk-size 25 \
  --out /data/results/test_1500_ref_free_strategy \
  --resume
```

## Large 50k sampled benchmark with the new reference-free strategy
```bash
python -m docking_pipeline benchmark-adaptive \
  --dataset-dir /data/datasets/rdock_50k \
  --strategy reference_free_active_learning_v2 \
  --reference-mode sampled \
  --reference-sample-size 5000 \
  --reference-sample-seed 42 \
  --evaluation-pool-mode same_pool \
  --triage-model classifier \
  --triage-controller auto_recall \
  --triage-retain-fraction 0.02 \
  --triage-target-recall 0.98 \
  --calibration-size 1500 \
  --fidelity-validation-size 300 \
  --balanced-baselines true \
  --jobs 64 \
  --chunk-size 50 \
  --cpu-fraction 0.85 \
  --out /data/results/rdock_50k_ref_free_strategy \
  --resume
```

## Repeated triage safety benchmark
```bash
python -m docking_pipeline benchmark-triage-repeated \
  --dataset-dir /data/datasets/test_1500 \
  --reference-mode full \
  --strategy reference_free_triage_bandit_v1 \
  --baselines cluster_only_triage,cheap_descriptor_filter_only,diverse_random_cost_balanced \
  --triage-model classifier \
  --classifier-top-percentile 0.05 \
  --classifier-threshold-mode recall_target \
  --triage-controller auto_recall \
  --model-fallback-if-worse union_with_cluster_only \
  --survivor-combination-policy union \
  --triage-target-recall 0.98 \
  --triage-retain-fraction 0.02 \
  --calibration-size 300 \
  --fidelity-validation-size 100 \
  --seeds 1,2,3 \
  --jobs 8 \
  --chunk-size 20 \
  --cpu-fraction 0.85 \
  --out /data/results/test_1500_triage_repeated \
  --force
```

## Results layout
- `target/`, `ligands/`, `rdock/`, `poses/`, `tables/`, `metrics/`, `plots/`, `checkpoints/`
- `report.md`, `manifest.json`, `config.yaml`, `commands.log`

## Candidate search for a 500k run
```bash
python scripts/find_pdb_ligand_candidates.py \
  --out configs/benchmarks/pdb_ligand_candidates_500k.tsv \
  --max-candidates 50 \
  --min-heavy-atoms 15 \
  --max-heavy-atoms 60 \
  --min-estimated-hits 500000
```

## Prepare a 500k dataset
```bash
export PDB_ID=<BEST_PDB>
export RECEPTOR_CHAIN=<CHAIN>
export REF_LIGAND=<RESNAME>
export REF_CHAIN=<CHAIN>
export DATASET_DIR=/data/datasets/rdock_${PDB_ID}_pubchem_500k

python scripts/prepare_pdb_ligand_dataset.py \
  --pdb-id $PDB_ID \
  --receptor-chain $RECEPTOR_CHAIN \
  --reference-ligand-resname $REF_LIGAND \
  --reference-ligand-chain $REF_CHAIN \
  --n-ligands 500000 \
  --ligand-source pubchem_similarity \
  --pubchem-threshold-ladder 99,95,90,85,80,75,70,65,60,55,50,45,40,35,30 \
  --allow-partial-ligand-set \
  --ligand-batch-size 100 \
  --ligand-jobs 60 \
  --ligand-cpu-fraction 0.85 \
  --out $DATASET_DIR \
  --force
```

## Validate the prepared dataset
```bash
python -m docking_pipeline validate-dataset \
  --dataset-dir $DATASET_DIR
```

## Smoke ablation benchmark
```bash
export SMOKE_OUT=/data/results/${PDB_ID}_500k_smoke_ablation

python -m docking_pipeline benchmark-regressor-ablation \
  --dataset-dir $DATASET_DIR \
  --reference-mode sampled \
  --reference-sample-size 1000 \
  --reference-sample-seed 42 \
  --fidelity-levels 5,10,15,30,50 \
  --cost-budget-runs 3000 \
  --calibration-size 500 \
  --fidelity-validation-size 150 \
  --seeds 1 \
  --jobs 32 \
  --chunk-size 50 \
  --cpu-fraction 0.85 \
  --out $SMOKE_OUT \
  --resume
```

## Full regressor ablation benchmark
```bash
export ABLATION_OUT=/data/results/${PDB_ID}_500k_regressor_ablation

python -m docking_pipeline benchmark-regressor-ablation \
  --dataset-dir $DATASET_DIR \
  --reference-mode sampled \
  --reference-sample-size 5000 \
  --reference-sample-seed 42 \
  --fidelity-levels 5,10,15,30,50 \
  --cost-budget-runs 200000 \
  --calibration-size 5000 \
  --fidelity-validation-size 1000 \
  --seeds 1,2,3 \
  --jobs 64 \
  --chunk-size 100 \
  --cpu-fraction 0.85 \
  --rdock-timeout-seconds 14400 \
  --out $ABLATION_OUT \
  --resume
```

## Full strategy benchmark
```bash
export STRATEGY_OUT=/data/results/${PDB_ID}_500k_strategy_benchmark

python -m docking_pipeline benchmark-triage-strategies \
  --dataset-dir $DATASET_DIR \
  --reference-mode sampled \
  --reference-sample-size 5000 \
  --reference-sample-seed 42 \
  --strategies cluster_only_triage,cheap_descriptor_filter_only,random_cost_balanced,diverse_random_cost_balanced,single_fidelity_cost_balanced,reference_free_active_learning_v2 \
  --adaptive-policies cluster_bandit \
  --triage-target-recall 0.98 \
  --triage-retain-fraction 0.02 \
  --calibration-size 5000 \
  --fidelity-validation-size 1000 \
  --cost-budget-runs 200000 \
  --fidelity-levels 5,10,15,30,50 \
  --seeds 1,2,3 \
  --jobs 64 \
  --chunk-size 100 \
  --cpu-fraction 0.85 \
  --rdock-timeout-seconds 14400 \
  --out $STRATEGY_OUT \
  --resume
```

## Comparable adaptive benchmark
Use this when you need a frozen evaluation universe and cost-comparable reporting.
```bash
export COMPARABLE_OUT=/data/results/${PDB_ID}_500k_comparable_benchmark

python -m docking_pipeline benchmark-comparable-adaptive \
  --dataset-dir $DATASET_DIR \
  --evaluation-universe-size 25000 \
  --evaluation-universe-seed 42 \
  --reference-mode full \
  --strategies adaptive_classifier_regressor,classifier_only,diverse_random_cost_balanced,cluster_only,single_fidelity_cost_balanced \
  --cost-budget-runs 100000 \
  --fidelity-levels 5,10,15,30,50 \
  --promotion-policy quota_ladder \
  --min-final-ligands 200 \
  --model-validation-split cluster \
  --jobs 120 \
  --chunk-size 10 \
  --cpu-fraction 0.85 \
  --rdock-timeout-seconds 1800 \
  --out $COMPARABLE_OUT \
  --resume
```

Prosty wrapper shell:
```bash
export DATASET_DIR=$SCRATCH/plgproma/data/rdock_${PDB_ID}_pubchem_500k
export STRATEGY_OUT=$SCRATCH/plgproma/results/${PDB_ID}_500k_strategy_benchmark
export JOBS=120
export CHUNK_SIZE=30
export CALIBRATION_SIZE=50000
export FIDELITY_VALIDATION_SIZE=50
export COST_BUDGET_RUNS=100000
export REFERENCE_SAMPLE_SIZE=500
export RDOCK_TIMEOUT_SECONDS=4000

bash scripts/run_strategy_benchmark_server.sh
```

## Audits
```bash
python -m docking_pipeline audit-strategy-benchmark \
  --run-dir $STRATEGY_OUT

python -m docking_pipeline audit-regressor-ablation \
  --run-dir $ABLATION_OUT
```

## Resume
- Re-run the same `benchmark-adaptive` command with `--resume`.
- Checkpoints are written under `checkpoints/`.

## Output interpretation
- `full_docking_scores.csv`: exhaustive final-fidelity docking reference.
- `random_baseline_scores.csv`: random baseline at matched cost budget.
- `single_fidelity_adaptive_scores.csv`: direct final-fidelity model-priority comparator.
- `multifidelity_trace.csv`: per-ligand trace across fidelity levels.
- `promotion_decisions.csv`: promotion and rejection reasons.
- `final_hits.csv`: hit table after multi-fidelity promotion logic.
- `adaptive_benchmark_metrics.json`: aggregate benchmark metrics and cost accounting.