Metadata-Version: 2.4
Name: rustina
Version: 0.4.0
Classifier: Development Status :: 5 - Production/Stable
Classifier: License :: OSI Approved :: MIT License
Classifier: Operating System :: POSIX :: Linux
Classifier: Programming Language :: Python :: 3
Classifier: Programming Language :: Rust
Classifier: Topic :: Scientific/Engineering :: Bio-Informatics
Requires-Dist: rdkit>=2024.3 ; extra == 'prep'
Requires-Dist: molscrub>=0.2.2 ; extra == 'prep'
Requires-Dist: meeko>=0.7.1 ; extra == 'prep'
Requires-Dist: scipy>=1.12.0 ; extra == 'prep'
Requires-Dist: gemmi>=0.6.0 ; extra == 'prep'
Provides-Extra: prep
License-File: LICENSE
Summary: CPU molecular docking with native GNINA-compatible CNN rescoring
Author: Adrien H. Cerdan
Requires-Python: >=3.9
Description-Content-Type: text/markdown; charset=UTF-8; variant=GFM
Project-URL: Homepage, https://github.com/AdrienCerdan/rustina
Project-URL: Issues, https://github.com/AdrienCerdan/rustina/issues
Project-URL: Repository, https://github.com/AdrienCerdan/rustina

# Rustina

Rustina is a CPU molecular-docking engine written in Rust. It provides Vina and
Vinardo empirical scoring, Rustina/Vina/QuickVina2 search, RILC and BFGS local
optimization, and native inference for the GNINA 1.3 default three-model CNN
ensemble and its optional distilled `fast` model. The supported v1 interfaces
are Python and the command-line program.

Rustina v1 targets Linux x86-64 and Python 3.9 or newer.

Search policy is explicit:

- `qvina2` is the default search policy. It uses the Vina budget with the
  official QuickVina 2 BFGS-history significance test.
- `rustina` is the opt-in optimized CPU policy with capped steps, LBFGS,
  early termination, restart/polish-reuse heuristics, and torsion-adaptive
  run scaling (8 runs for rigid ligands, 12 for mid, 16 for flexible).
- `atlas` is the opt-in interaction-field and torsion-beam search. It races a
  diverse candidate population through score-only, shallow, and full local
  refinement stages; existing modes and defaults are unchanged. It remains an
  experimental mode until the documented release-qualification gates pass.
- `hybrid` merges the native and ATLAS basins, deduplicates them by heavy-atom
  RMSD, applies one shared direct-refinement pass, and then performs one shared
  empirical/CNN ranking. Per-pose diagnostics identify `native`, `atlas`, or
  `native+atlas` provenance. It shares ATLAS's experimental qualification
  status.
- `vina` follows the AutoDock Vina 1.2.3 Monte Carlo/BFGS policy and defaults
  to 0.375 A grid spacing.

The compatibility modes reject custom steps, BFGS iteration counts, and RILC.
Use `search_mode="rustina"` for experimental optimizer or budget combinations.

Rustina v1 also includes opt-in template docking for congeneric ligands through
`dock_reference()` and `screen_reference()`.

## Installation

Install the Python package:

```bash
python -m pip install rustina
```

Ligand preparation is an optional extra because it installs RDKit, molscrub,
Meeko, SciPy, and Gemmi:

```bash
python -m pip install "rustina[prep]"
```

For a source checkout, build a production extension with:

```bash
maturin develop --release
```

## Python API

Inputs to docking and scoring are prepared PDBQT paths or raw PDBQT strings.
CNN rescoring defaults to GNINA's distilled `fast` model. Set `cnn="default"`
for Rustina's bundled GNINA 1.3 ensemble, `skip_cnn=True` to disable
rescoring, or supply a native weight path through `cnn`.

```python
import rustina

poses = rustina.dock(
    receptor="data/receptor.pdbqt",
    ligand="data/ligand.pdbqt",
    cx=-14.0,
    cy=18.0,
    cz=-15.0,
    sx=14.0,
    sy=18.0,
    sz=15.0,
    runs=8,
    threads=8,
    seed=42,
    direct_refine_candidates=8,
)
print(poses[0]["affinity"], poses[0]["cnn_score"])
```

`cnn_rotations` (CLI `--cnn-rotations`, 1-24, default 1) averages the CNN pose
score over that many orientations of the voxel grid, the equivalent of GNINA's
`--cnn_rotations`. Only the rescore is affected: the search, the poses, and every
empirical number are identical at any setting. The default of 1 is deliberate.
Measured over 361 redocking targets by rescoring stored poses at all 24
orientations, averaging changed the chosen pose on 4 of 361 targets and never
improved one, while costing about 7 s per target on six cores. See
`CLAUDE.md` for the measurement.

The low-level CLI and Python API require an explicit pocket center and size.
The product-style benchmark driver defaults to ligand-sized dynamic boxes,
using the crystallographic heavy-atom extent plus adaptive padding
`max(4 Å, 8 Å - 0.2 Å × torsions)`. The fixed-25 release qualification remains
explicitly pinned to its published protocol.

Screen multiple ligands while reusing the receptor grid:

```python
results = rustina.screen(
    receptor="data/receptor.pdbqt",
    ligands=["ligand-1.pdbqt", "ligand-2.pdbqt"],
    cx=-14.0,
    cy=18.0,
    cz=-15.0,
    sx=14.0,
    sy=18.0,
    sz=15.0,
    threads=8,
)
```

Experimental negative-image screening compiles typed pocket hotspots once,
generates rigid pharmacophore alignments, and selects/refines them with the
atom-specific Vina grids. It is a retrieval prototype, not a replacement for
final docking:

```python
pocket = rustina.prepare_negative_image(
    "data/receptor.pdbqt",
    cx=-14.0, cy=18.0, cz=-15.0,
    sx=14.0, sy=18.0, sz=15.0,
)
matches = pocket.screen(["ligand-1.pdbqt", "ligand-2.pdbqt"], threads=6)

# Pose-generation diagnostics: refine and serialize up to 32 distinct basins.
candidates = pocket.screen_candidates(
    "ligand-1.pdbqt", limit=16, refine_candidates=32
)
```

Each successful record reports physical grid energy, pharmacophore coverage,
clash diagnostics, pre/post-refinement transforms, selection provenance, and a
ready-to-score PDBQT pose. Malformed ligands return an isolated
`{"error": ...}` record. Use `scripts/benchmark_negative_image.py` for
versioned teacher-recall experiments; do not interpret the prototype score as
experimental binding affinity.

Before promoting RIFT as a docking prescreen, run
`scripts/run_rift_pose_qualification.py`. It evaluates randomized crystal,
single-ETKDG, and eight-ETKDG inputs on the fixed-box Astex diagnostic panel,
reporting raw proposal RMSD separately from grid-refined RMSD and physical
PoseBusters validity. Qualification and retrieval benchmark drivers print
configuration/target progress as they run and keep resumable raw results under
`scratch/`.

Optional ligand preparation accepts SMILES, structure files/blocks, or an
RDKit molecule:

```python
states = rustina.prepare_ligand("CC(=O)Nc1ccc(O)cc1", random_seed=42)
print(states[0]["pdbqt"])
```

### Template docking

Template docking is a supported, opt-in v1 feature for congeneric series when a
ligand with known coordinates is already in the receptor coordinate frame. It requires the
`prep` extra because maximum-common-substructure matching is performed with
RDKit:

```python
poses = rustina.dock_reference(
    receptor="receptor.pdbqt",
    ligand="query.sdf",
    reference="co-crystal-ligand.sdf",
    cx=-14.0,
    cy=18.0,
    cz=-15.0,
    sx=14.0,
    sy=18.0,
    sz=15.0,
    runs=8,
    seed=42,
)
print(
    poses[0]["relaxed_core_rmsd"],
    poses[0]["physical_score"],
    poses[0]["reference_satisfied"],
)
```

The mapped heavy-atom core is restrained by a soft flat-bottom potential during
search, followed by a short unrestrained relaxation. Results report constrained
and relaxed core RMSD, physical and guided scores, the restraint penalty,
mapping identity, both MCS coverage fractions, and whether the relaxed core
remains within the default 1.0 Angstrom satisfaction threshold.

Automatic guidance requires at least six mapped heavy atoms and 50% query
coverage. An explicit `atom_map` can define a smaller anchor. Rustina fails
rather than silently switching to free docking when automatic guidance does not
meet these gates. `screen_reference()` applies the same workflow to a series
while reusing one receptor grid.

Up to eight symmetry-distinct MCS mappings are evaluated deterministically.
The exact requested run budget is distributed across them globally; mappings do
not multiply the run count. Candidates are merged and deduplicated before one
CNN rescore using the requested `cnn_pool_size`.

#### Timing smoke test

A release-mode `5SAK_ZRY` methyl-analog smoke test used four total runs, four
mapping hypotheses, empirical scoring, matched pose-pool budgets, and three
seeds. Median wall times were:

| CPU threads | Normal docking | Template docking | Difference |
|---:|---:|---:|---:|
| 2 | 0.711 s | 0.732 s | +3% |
| 8 | 0.620 s | 0.682 s | +10% |

The two-thread end-to-end template call, including MCS generation and final
relaxation, took 0.758 s median. These numbers characterize one small smoke
case, not expected performance across ligand series. The solved-congeneric
qualification requirements are documented in
[`docs/benchmarks/REFERENCE_DOCKING_PROTOCOL.md`](docs/benchmarks/REFERENCE_DOCKING_PROTOCOL.md).

Template docking assumes the reference and target receptor use the same
coordinate frame. It does not align receptor structures or provide shape-only
or pharmacophore guidance. `5SAK_ZRY` demonstrates restraint behavior and
physical validity only; its methyl analog has no experimental pose and is not
an accuracy benchmark.

`rustina.build_profile()` returns `debug` or `release`. Performance results are
valid only when the actually imported extension reports `release`.

## Command line

```bash
rustina dock \
  --receptor data/receptor.pdbqt \
  --ligand data/ligand.pdbqt \
  --output docked.pdbqt \
  --cx -14 --cy 18 --cz -15 \
  --sx 14 --sy 18 --sz 15 \
  --runs 8 --threads 8 --seed 42

rustina score \
  --receptor data/receptor.pdbqt \
  --ligand data/ligand.pdbqt \
  --cnn
```

Use `rustina dock --help` and `rustina score --help` for the complete supported
options.

## Reproducibility and benchmarks

Docking is deterministic when `seed` is provided. Release benchmarks must
record the Rustina version, build profile, input dataset revision, complete
arguments, CPU model, and raw per-target results. Historical pre-v1 research is
preserved by the `pre-v1-research` tag; v1 benchmark qualification lives under
`docs/benchmarks/`.

### PoseBusters + Astex Benchmark

All numbers below come from a single binary (`97d0dd73`), seed 42, budget 8,
393 targets (308 PoseBusters + 85 Astex), pinned to six physical cores on a
Ryzen 5 5600X. Wilson 95% intervals are on the all-target denominator; at ~390
targets the minimum detectable difference is roughly ±4 points, so do not read
a smaller gap as real without a paired test.

Superseded numbers: this section previously reported 79.1% over "485 evaluated
targets" and 68.6% over "433" — both larger than the 393 targets that exist.
Those rows pooled records from several different binaries into one figure. See
`docs/benchmarks/MERGE_QUALIFICATION_RESULTS.md`.

#### Controlled protocol — the headline

Generated ETKDGv3+UFF start conformers and the pinned fixed 25 Å box;
crystallographic coordinates are evaluation-only. This is the protocol
comparable to published PoseBusters results.

| Engine          | Mode | RMSD<=2A (95% CI)  | PB-valid | Mean s |
|-----------------|------|--------------------|----------|--------|
| rustina         | fast | 53.7% (48.7-58.6)  | 95.7%    | 2.18   |
| rustina_qvina2  | fast | 63.9% (59.0-68.5)  | 96.9%    | 4.83   |

On this protocol the `qvina2` compatibility mode outperforms native `rustina`
by 10.2 points — but it also takes 2.2x as long at the same `--runs 8`, so an
unknown part of that gap is simply more search. A budget-matched comparison has
not been run. Earlier releases claimed native mode *surpassed* qvina2; that
comparison was made on a ligand-sized box, where the two modes differ by 2.6
points, which is inside the noise floor.

#### Experimental protocol — a labelled upper bound, not a headline

Crystal ligand conformers fed as docking input and a box sized from the crystal
ligand extent. Both leak the answer, so this runs about 28 points high. It is
published because it is the corpus used for offline re-ranking research, not
because it is a fair number.

| Engine  | Mode | RMSD<=2A | Top-5 <=2A | PB-valid | Any pose <=2A | Mean s |
|---------|------|----------|------------|----------|---------------|--------|
| rustina | fast | 81.4%    | 92.6%      | 98.0%    | 97.5%         | 1.84   |

The gap between 81.4% top-1 and the 97.5% sampling ceiling is the engine's real
limitation: near-native poses are usually generated and then ranked below #1.

CNN `fast` is the recommended mode: it costs about 1.1x empirical. CNN
`default` (the GNINA 1.3 ensemble) is marginally better and 1.3-1.9x slower.

Rustina is research software. Docking scores and predicted poses are not a
substitute for experimental evidence or clinical decision-making.

## Attribution and license

Rustina is MIT licensed. Its algorithms and bundled CNN parameters build on
AutoDock Vina, Smina, QuickVina2, GNINA, Vina-GPU, and related published work.
See [THIRD_PARTY_NOTICES.md](THIRD_PARTY_NOTICES.md) and
[CITATION.cff](CITATION.cff) for provenance and citations. Model conversion is
documented in [docs/models.md](docs/models.md).

