Pharmacological Lattice Quantisation without Number-Theory
Description
This dataset supports a two-layer computational platform for neglected-disease drug discovery, focused on antimicrobial resistance and malaria: (A) a ligand-based drug-triage and de-risking platform, and (B) a structure-based generative design pipeline (SBDD) built on top of it. It provides cleaned, reproducible, honestly-evaluated bioactivity data together with the results needed to judge model reliability, plus the full structure-based pipeline and the candidate molecules it produced. Data were derived from two open resources: ChEMBL (bioactivity, CC BY-SA 3.0) and the Open Targets Platform (human target-disease genetics, CC BY 4.0). Bioactivity records were standardised to a p-scale (p = -log10 of the molar concentration; pIC50 for binding assays, pMIC for whole-cell antibacterial assays, converting ug/mL via molecular weight), deduplicated by InChIKey with replicate measurements aggregated by median, and canonicalised with RDKit. Each compound carries a Bemis-Murcko scaffold split (train/test) so the held-out evaluation can be reproduced exactly, with no analog leakage between splits. Building on the ligand layer, a pathogen-agnostic structure-based pipeline discovers essential-but- undrugged protein targets and designs candidate binders for them. For a chosen pathogen it: (1) finds targets that are essential (real gene-essentiality data - BV-BRC for bacteria; PlasmoDB piggyBac Mutagenesis Index Scores for Plasmodium), undrugged, host-selective (no close human homolog), and druggable (site-anchored pocket hydrophobicity/enclosure), with resistance barriers from CARD; (2) fetches the AlphaFold structure and UniProt catalytic residues; (3) scores small-molecule binding with a positive-control-validated Boltz-2 co-folding oracle returning a binding probability and a 3D pose - validated per target against known drug/target pairs (e.g. DHFR/methotrexate 0.998 vs decoy 0.13; A. baumannii FabI/triclosan lifted 0.577->0.912 in the FabI+NAD+ ternary complex), where rigid docking (AutoDock Vina) failed the same control and the oracle is shown deterministic (sigma ~0.003); (4) checks whether the co-folded pose engages the catalytic site; (5) generates candidate molecules with a Markov-chain / fragment-crossover / substrate-growing generator seeded from each target's native substrate or inhibitor pharmacophore, under a baked-in reactive-group/PAINS/charge/catechol/ thiocarbonyl filter that prevents oracle-gaming (documented and rejected artifacts: boron and thiolate exploits); and (6) diagnoses the binding ceiling as a search, oracle, target, or cofactor limit (e.g. the FabI ceiling was pre-registered as cofactor-limited and confirmed by the NAD+ lift).
Files
Steps to reproduce
Environment. Python 3.13. Install: pip install rdkit scikit-learn scipy numpy matplotlib. All steps are deterministic (random_state=42); network calls to ChEMBL and Open Targets are cached locally, so a second run is offline and identical. 1. Acquire bioactivity data (ChEMBL REST API). For each target/organism, query activity.json filtered to standard_relation="=": malaria = CHEMBL364 (IC50), S. aureus = CHEMBL352 (MIC), E. coli = CHEMBL354 (MIC), acetylcholinesterase = CHEMBL220 (IC50), EGFR = CHEMBL203 (IC50), hERG = CHEMBL240 (IC50). Whole-organism targets are queried with the organism name, not filtered to human. Retain canonical_smiles, standard_value, standard_units. 2. Clean and standardise. Parse each SMILES with RDKit; drop unparseable structures. Convert potency to a p-scale: p = -log10(molar), where molar units (nM, µM…) convert directly and mass-concentration units (µg/mL) convert via RDKit molecular weight molar = value·1e-3 / MW. Keep only 0 ≤ p ≤ 14. 3. Deduplicate. Compute each compound's InChIKey (RDKit). Group records by InChIKey and aggregate replicate p-values by median. This yields one row per unique compound with n_measurements. 4. Featurize. 2048-bit ECFP4 count fingerprints (RDKit GetMorganGenerator, radius 2), kept as full vectors. For similarity/applicability domain, binarize (>0) and use max Tanimoto to the training set. 5. Scaffold split (no leakage). Assign each compound its Bemis–Murcko scaffold (RDKit MurckoScaffold); order scaffolds largest-first and fill 80% into train, remainder into test, so whole scaffolds never straddle the split. This scaffold_split column reproduces the reported metrics exactly. 6. Train and evaluate. RandomForestRegressor(n_estimators=400, min_samples_leaf=2, random_state=42) on the train split; report on the held-out test split: Spearman, Pearson, RMSE, versus a mean-prediction baseline. 7. Abstention / risk–coverage. Per-prediction uncertainty = standard deviation across the 400 trees. Rank test compounds by ascending uncertainty; recompute Spearman/RMSE on the most-confident top-k for k = 100%…10% coverage. 8. hERG liability (classification). Label blocker = pIC50 > 6 (IC50 < 1 µM). Reclaim censored > records with value ≥ 1000 nM as confident non-blockers. RandomForestClassifier(400, class_weight="balanced", random_state=42); report ROC-AUC and PR-AUC on the scaffold-held-out set. 9. Biology axis (Open Targets GraphQL, api.platform.opentargets.org/api/v4/graphql). For each molecular target (Ensembl ID), query associatedDiseases for overall and genetic_association datatype scores; query evidences on gene_burden/gwas_credible_sets/eva for rows with directionOnTrait="protect" (protective loss-of-function). 10. Regenerate the whole package. Running python -m affinity.export_dataset executes steps 1–9 from cache and writes all data/ tables (with scaffold-split columns), results/model_metrics.csv, the figures, checksums, and this protocol, deterministically.