diff --git a/docs/installation.md b/docs/installation.md index 1be39ce..584198e 100644 --- a/docs/installation.md +++ b/docs/installation.md @@ -96,40 +96,49 @@ Installation instructions for running the database on bare metal, rather than Do The ProCogGraph pipeline is built using Nextflow for workflow management. To run the pipeline, follow these steps: -1. The pipeline utilises data from a number of different sources to build the ProCogGraph database. To begin, prepare a data files directory with the following: - - | File | Description | Download | - | ---- | ---- | ---- | - | pdb_chain_enzyme.tsv.gz | Protein chain EC ID annotation from SIFTS for PDe structures. | [SIFTS](https://ftp.ebi.ac.uk/pub/databases/msd/sifts/csv/pdb_chain_enzyme.csv) | - | assemblies_data.csv.gz | Assembly data for PDBe structures from PDBe-KB | [PDBe-KB](https://ftp.ebi.ac.uk/pub/databases/pdbe-kb/complexes/assemblies_data.csv) | - | enzclass.txt | Enzyme classification hierarchy | [ExPASy](https://ftp.expasy.org/databases/enzyme/) | - | enzyme.dat | ENZYME database records | [ExPASy](https://ftp.expasy.org/databases/enzyme/) | - | cath-names.txt | CATH domain names | [CATH](http://download.cathdb.info/cath/releases/latest-release/cath-classification-data/) | - | cath-domain-description-file.txt | CATH domain descriptions | [CATH](http://download.cathdb.info/cath/releases/latest-release/cath-classification-data/) | - | dir.des.scop.1_75.txt | SCOP domain descriptions | [SCOP](https://scop.berkeley.edu/downloads/) | - | dir.cla.scop.1_75.txt | SCOP domain classifications | [SCOP](https://scop.berkeley.edu/downloads/) | - | clan_membership.txt.gz | Pfam clan membership | [InterPro](https://ftp.ebi.ac.uk/pub/databases/Pfam/current_release/database_files/) | - | clan.txt.gz | Pfam clan descriptions | [InterPro](https://ftp.ebi.ac.uk/pub/databases/Pfam/current_release/database_files/) | - | interpro.xml.gz | InterPro domain annotations | [InterPro](https://ftp.ebi.ac.uk/pub/databases/interpro/current_release/) | - | rhea-reaction-smiles.tsv | RHEA reaction smiles strings | [RHEA](https://www.rhea-db.org/help/download) | - | rhea2ec.tsv | RHEA to EC number mappings | [RHEA](https://www.rhea-db.org/help/download) | - | rhea-directions.tsv | RHEA reaction directions | [RHEA](https://www.rhea-db.org/help/download) | - | chebi_names.tsv.gz | ChEBI names | [ChEBI](https://ftp.ebi.ac.uk/pub/databases/chebi/Flat_file_tab_delimited/) | - | relation.tsv | ChEBI relations | [ChEBI](https://ftp.ebi.ac.uk/pub/databases/chebi/Flat_file_tab_delimited/) | - | ChEBI_Results.tsv | ChEBI records with database cross references to KEGG GLYCAN and KEGG COMPOUND, where a structure exists for the record, generated with advanced search function. | [ChEBI](https://www.ebi.ac.uk/chebi/advancedSearchForward.do) | - | scop2-cla-latest.txt | SCOP2 domain classifications | [EBI](https://www.ebi.ac.uk/pdbe/scop/download) | - | scop2-des-latest.txt | SCOP2 domain descriptions | [EBI](https://www.ebi.ac.uk/pdbe/scop/download) | - | ccd.cif | Chemical Component Dictionary Structures | [CCD](https://www.wwpdb.org/data/ccd) | - | pubchem_substance_id_mapping.txt | PubChem substance ID mappings from PubChem search for KEGG data source. | [PubChem](https://www.ncbi.nlm.nih.gov/pcsubstance?term=%22KEGG%22%5BSourceName%5D%20AND%20hasnohold%5Bfilt%5D) | - -2. Clone this repository and install dependencies: +1. Clone this repository and install dependencies: ``` bash git clone m-crown/ProCogGraph cd ProCogGraph - conda env create -f nextflow/envs/environment.yml + conda env create -f nextflow/envs/procoggraph.yaml + conda activate procoggraph ``` +2. The pipeline utilises data from a number of different sources to build the ProCogGraph database. Rather than downloading these manually, run the reference data download script, which fetches, validates, and (where needed) reformats every file the pipeline requires into a single data directory: + + ``` bash + python3 nextflow/bin/download_reference_data.py --data_dir /PATH/TO/DATA_DIR + ``` + + Useful flags: + - `--dry-run` — list what would be fetched without downloading anything. + - `--only NAME[,NAME...]` — fetch a subset of entries (see `nextflow/bin/reference_data_manifest.yaml` for names). + - `--force` — re-download and refresh files that already exist. + + Every file below is fetched automatically (the script will print instructions and skip, rather than fail the run, for any future entry that turns out to have no stable bulk-download source): + + | File | Description | + | ---- | ---- | + | pdb_chain_enzyme.tsv.gz | Protein chain EC ID annotation from SIFTS for PDBe structures. | + | assemblies_data.csv | Assembly data for PDBe structures from PDBe-KB. | + | enzclass.txt | Enzyme classification hierarchy (ExPASy). | + | enzyme.dat | ENZYME database records (ExPASy). | + | cath-names.txt | CATH domain names. | + | cath-domain-description-file.txt | CATH domain descriptions. | + | dir.des.scop.1_75.txt / dir.cla.scop.1_75.txt | SCOP domain descriptions/classifications. | + | scop2-cla-latest.txt / scop2-des-latest.txt | SCOP2 domain classifications/descriptions. | + | Pfam-A.clans.tsv.gz | Pfam-A domain descriptions and clan membership (consolidated by Pfam into a single file in 2026). | + | interpro.xml.gz | InterPro domain annotations. | + | rhea2ec.tsv / rhea-directions.tsv / rhea-reaction-smiles.tsv | RHEA reaction data. | + | chebi_names.tsv.gz / relation.tsv | ChEBI names and relations. | + | ChEBI_Results.tsv | ChEBI records with a KEGG COMPOUND cross-reference and a known structure — derived from bulk ChEBI flat files rather than ChEBI's advanced search UI (see `docs/reference_data_download_plan.md`). | + | ccd.cif | Chemical Component Dictionary structures. | + | pubchem_substance_id_mapping.txt | PubChem CID ↔ KEGG code cross-references — derived from NCBI's bulk `CID-Identifiers.tsv.gz` rather than an interactive PubChem substance search (see `docs/reference_data_download_plan.md`). | + | rd/ | Directory of Rhea per-reaction (`.rd`) files, extracted from Rhea's `rhea-rd.tar.gz` bulk archive. | + + `pfamA.txt.gz`/`clan_membership.txt.gz`/`clan.txt.gz` in older versions of this doc have been consolidated upstream into the single `Pfam-A.clans.tsv.gz` above; the download script fetches it and the pipeline's parsing code (`utils.get_pfam_annotations`) reads it directly, so no separate files are needed. `rd/` — a directory of Rhea reaction (`.rd`) files, needed by `preprocess_rhea.py --rd_dir` — is also fetched automatically (extracted from Rhea's `rhea-rd.tar.gz` bulk archive). + 3. Preprocess RHEA reaction files: ``` bash @@ -140,7 +149,7 @@ The ProCogGraph pipeline is built using Nextflow for workflow management. To run 4. Produce final manifest file of structures to be processed: ``` bash - python3 download_mmcif.py --sifts_file /PATH/TO/DATA_DIR/pdb_chain_enzyme.tsv.gz --assemblies_file /PATH/TO/DATA_DIR/assemblies_data.csv.gz --chunk_size 100 --output_dir /PATH/TO/STRUCTURES_DIR + python3 download_mmcif.py --sifts_file /PATH/TO/DATA_DIR/pdb_chain_enzyme.tsv.gz --assemblies_file /PATH/TO/DATA_DIR/assemblies_data.csv --chunk_size 100 --output_dir /PATH/TO/STRUCTURES_DIR ``` 5. Run the nextflow pipeline: @@ -155,7 +164,7 @@ The ProCogGraph pipeline is built using Nextflow for workflow management. To run - params.output_dir - the desired output directory. - - params.manifest - the path to the manifest file created in step 3. + - params.manifest - the path to the manifest file created in step 4. ``` bash cd /PATH/TO/PROCOGGRAPH_REPOSITORY/nextflow diff --git a/docs/reference_data_download_plan.md b/docs/reference_data_download_plan.md new file mode 100644 index 0000000..a17ca4d --- /dev/null +++ b/docs/reference_data_download_plan.md @@ -0,0 +1,303 @@ +# Plan: Automate External Reference Data Download + +**Status (2026-08-02, branch `feature/reference-data-download`, off +`fix/pipeline-wiring-bugs`):** Implemented as +`nextflow/bin/reference_data_manifest.yaml` + +`nextflow/bin/download_reference_data.py`, with `direct_url`/`derived`/ +`manual`/`needs_code_update` source types, tested end-to-end against live +sources (not just written and assumed correct) during implementation. +Notable deviations/fixes found only by actually running it against real +data, not by re-reading this plan: + +- The `ChEBI_Results.tsv` derivation logic in this plan's §2 was wrong in + three ways only caught by inspecting the real files: ChEBI flat-file + columns are lowercase, not uppercase as assumed; + `database_accession.tsv.gz`'s `type` column has no KEGG-specific value + at all (KEGG cross-references are `type == "MANUAL_X_REF"` rows filtered + by `source_id` against `source.tsv.gz`, a file this plan never + mentioned); and `structures.tsv.gz` has a direct `smiles` column per + compound, not a melted `TYPE`/`STRUCTURE` format. Also: the derived `ID` + column must be `CHEBI:`-prefixed and KEGG-COMPOUND-only (not also + KEGG GLYCAN, which this plan's §2 originally included) to match what + `get_ec_information.py` actually merges on and what its ChEBI + `has_role`-based cofactor-tagging regex (`:682`) expects — a bare + integer ID would have silently zeroed out that regex match rather than + erroring. Verified by diffing derived output against the repo's existing + manually-exported `data_files/ChEBI_Results.tsv`. +- SIFTS's only live source is a comma-separated, uncompressed `.csv`, not + the tab-separated gzipped `.tsv.gz` `download_mmcif.py` reads — added a + `csv_to_tsv_gz` post-process step (not anticipated by this plan) rather + than changing the pipeline script. +- SCOP (`dir.cla.scop.1_75.txt`/`dir.des.scop.1_75.txt`), SCOP2 + (`scop2-cla-latest.txt`/`scop2-des-latest.txt`), the CCD + (`components.cif.gz`), and the three RHEA TSVs could not be confirmed + live during implementation (inconclusive search results / a directory + listing rather than the file itself) — the manifest uses best-known + conventional URLs for these, each marked `confidence: unverified` with a + note. Run `--dry-run --only ` and then a real fetch for these + specifically before depending on them; update the manifest if any 404s. +- `preprocess_rhea.py --rd_dir` needs a directory of Rhea `.rd` files with + no bulk source identified — left as a `manual`-style gap (`rhea_rd_dir` + entry) until the user pointed at Rhea's `rhea-rd.tar.gz` bulk archive + (sibling to the 37,100-file `rd/` directory Rhea's own README says to + avoid browsing directly); resolved via a new `extract_tar_gz` + post-process and directory-aware size checking added to the script. +- **Update:** the Pfam clan-file consolidation (`Pfam-A.clans.tsv.gz`) is + no longer a gap — `utils.get_pfam_annotations` was rewritten to consume + the single consolidated file directly (`pfam_a_file`/`pfam_clan_rels`/ + `pfam_clans` collapsed into one `pfam_clans_file` param throughout + `nextflow.config`/`main.nf`/`process_all_pdb_contacts.py`). One real + data loss found along the way: the consolidated file has no equivalent + of the old `clan.txt.gz`'s free-text `clan_comment` field (confirmed — + no supplementary clan-detail file exists on Pfam's FTP anymore); + `clan_description` now holds the short clan id (e.g. "GPCR_A") and + `clan_comment` is left null rather than faking it. Verified end-to-end + against the real file. +- All `scop_domains_info`/`scop_descriptions`/`scop2_domains_info`/ + `scop2_descriptions`/`ccd`/`rhea2ec`/`rhea_directions`/ + `rhea_reaction_smiles`/`rhea_rd_dir`/`pfam_clans` entries are now + `confidence: verified` (user supplied several of the real source pages; + each was then re-verified end-to-end via the actual script and, where + applicable, the exact downstream pipeline parsing logic). +- **Update:** `pubchem_substance_mapping` is also resolved — user found + NCBI's bulk `CID-Identifiers.tsv.gz` (12.7M CID/identifier/type rows); + filtering to `identifier_type == "KEGG ID"` gives exactly the CID/KEGG + pairs needed (26,933 rows), turned into a `derived` entry + (`derive_pubchem_kegg_mapping`) that reproduces + `get_ec_information.py`'s expected `"CID: KEGG:"` line format. + Verified by feeding the derived output through the real parser and its + downstream CID-matching logic. **Every entry in the manifest is now + either `direct_url`/`derived` with `confidence: verified`, or a + self-contained code-only fix (Pfam) — zero remaining `manual` entries.** + +`docs/installation.md` updated to point at the script in place of the old +manual download table. + +Currently, running the Nextflow pipeline from scratch requires manually +downloading ~20 external reference files (SIFTS, PDBe-KB, ExPASy, CATH, +SCOP, SCOP2, Pfam/InterPro, RHEA, ChEBI, CCD, PubChem) into a single +`data_dir`, per the table in [`docs/installation.md`](installation.md#L99). +This is tedious, error-prone (wrong filename/location silently breaks a +pipeline stage later), and undocumented in a machine-checkable way. This +plan proposes a single script that fetches and validates all of them. + +## 1. Build a authoritative manifest of files needed + +Cross-referencing the installation doc's table against actual +`nextflow.config` params (`nextflow/nextflow.config:86-102`) and script +argparse definitions turned up two gaps to fix as part of this work, not +after: + +- **`pfamA.txt.gz`** is required by `params.pfam_a_file` + (`nextflow.config:95`, consumed by `utils.get_pfam_annotations`), but is + **not listed** in the installation doc's download table at all. +- **`scop2-des-latest.txt`** is listed in the doc's table but has **no + corresponding `params.scop2_descriptions_file` entry in + `nextflow.config`**, even though `process_all_pdb_contacts.py` and + `produce_neo4j_files.py` both take `--scop2_descriptions_file` as a + required argument (this is fix_plan.md item 2 — will need a + `nextflow.config` param added regardless of this plan, but the download + script must fetch the file either way). + +So step 1 is: build one manifest (e.g. +`nextflow/bin/reference_data_manifest.yaml`) that is the single source of +truth for "what files does the pipeline need, and where do they come +from," covering all params referenced in `nextflow.config:86-102` plus +`scop2_descriptions_file`. The installation doc's table and +`nextflow.config` should both be generated from / checked against this +manifest, not maintained separately, so they can't drift again. + +Each manifest entry needs: `name`, `target_filename`, `param_name` (the +matching `nextflow.config` param, for cross-checking), `source_type` +(`direct_url` / `latest_release_dir` / `manual`), `url`, and optionally +`post_process` (e.g. `gunzip`, `none`). + +## 2. Handle three distinct source patterns + +The current file list isn't uniformly "download this URL" — group by +pattern so the script can treat each correctly: + +- **Direct, stable URL** (majority of files): SIFTS, PDBe-KB assemblies, + RHEA files, ChEBI names/relations, CCD — just `curl`/`requests.get`. +- **"Latest release" directory listings**: CATH + (`cath-classification-data/`), SCOP2 (`pdbe/scop/download`), Pfam/InterPro + (`current_release/`) don't have one fixed filename/version — the script + needs to resolve the actual file to fetch (e.g. parse a directory index + or use a known-stable "latest" alias URL if the source provides one). + Where the source has no stable "latest" alias, log the resolved + version/date so re-runs are reproducible and it's visible when upstream + data has moved on. + + **Path correction (checked live, 2026-08):** Pfam has also restructured. + `clan_membership.txt.gz` and `clan.txt.gz` — both named explicitly in + `docs/installation.md`'s download table — **no longer exist** at + `https://ftp.ebi.ac.uk/pub/databases/Pfam/current_release/`. Clan data is + now consolidated into a single file, **`Pfam-A.clans.tsv.gz`** (541K). + `pfamA.txt.gz` (13M, already flagged as missing from the doc's table in + step 1) is still present under that name. The manifest entries for both + `pfam_clan_rels`/`pfam_clans` params (`nextflow.config:93-94`) need to + point at the one consolidated file instead of two, and + `utils.get_pfam_annotations` (`nextflow/bin/utils.py:232-241`), which + currently expects two separate files (`clan_membership_file`, + `clan_info_file`), will need updating to parse the merged format — this + is a pipeline-code change, not just a download-path fix, so it should be + scoped alongside `fix_plan.md` rather than silently patched inside the + download script. +- **Derived from bulk files** (not manual): `ChEBI_Results.tsv` was + originally produced via ChEBI's advanced search UI + (`docs/installation.md:119`), which looked non-scriptable at first + glance. Tracing its actual consumer + (`nextflow/bin/get_ec_information.py:531-536`) shows it only needs four + columns — ChEBI ID, ChEBI name, SMILES, and KEGG COMPOUND + ACCESSION — restricted to ChEBI entries that (a) have a KEGG + COMPOUND/GLYCAN cross-reference and (b) have a structure on file. Both + of those facts are published by ChEBI as plain bulk flat files. + + **Path correction (checked live, 2026-08):** ChEBI has since + restructured its FTP layout entirely. The + `Flat_file_tab_delimited/` directory referenced by + `docs/installation.md` no longer exists — it's now + `https://ftp.ebi.ac.uk/pub/databases/chebi/flat_files/`, and the + filenames changed too: + - `chebi_names.tsv.gz` → **`names.tsv.gz`** (9.0M compressed) + - `relation.tsv` → **`relation.tsv.gz`** (now gzipped, 2.6M) + - the cross-reference table needed here is **`database_accession.tsv.gz`** + (now gzipped, 3.8M) — filter `TYPE` to `KEGG COMPOUND accession` / + `KEGG GLYCAN accession`. + - the structure table is **`structures.tsv.gz`** (not `.csv.gz` — + renamed, and larger than assumed at 88M compressed) — filter `TYPE` + to `SMILES`. + + Any manifest entry or doc text written before this check needs the + corrected path/filename/extension, not just the derivation logic. This + is exactly the kind of drift the manifest in step 1 is meant to catch + early via a periodic `--dry-run` / size-sanity check, rather than + silently downloading a 404 page or failing a `pd.read_csv`. + + So this isn't a `source_type: manual` entry — it's a `derived` entry: + download both bulk files (now `.gz`), inner-join them on ChEBI compound + ID, merge in the name from `names.tsv.gz`, and write out + `ChEBI_Results.tsv` in the same shape the pipeline already expects. This + is strictly better than the original manual export too: it's + reproducible and pinned to a specific ChEBI release rather than a + point-in-time UI query result. +- **Manual-only, no stable download URL**: none currently identified. Keep + `source_type: manual` as a category in the manifest schema for future + entries, but the manifest should end up with zero of them once the + `derived` handling above is implemented. If a genuinely manual entry + ever shows up, the script should skip it, print clear instructions, and + continue rather than fail the whole run. + +## 3. Script design + +New script: `nextflow/bin/download_reference_data.py` + +- `--data_dir` (required) — output directory, matches `params.data_dir`. +- `--manifest` (default: bundled `reference_data_manifest.yaml`). +- `--only NAME[,NAME...]` — fetch a subset (useful when only one upstream + source has updated). +- `--skip-existing` (default true) / `--force` — re-download even if the + target file already exists, for refreshing stale reference data. +- `--dry-run` — print what would be fetched/skipped without downloading. +- Behaviour per file: + 1. Skip if present and `--skip-existing` (default), unless `--force`. + 2. Download with retry/backoff (a handful of these hosts — EBI FTP, + Zenodo-adjacent mirrors — are occasionally flaky) and a sane timeout. + 3. Verify the download is non-empty / matches expected minimum size as a + sanity check (no fixed checksums available from most of these + upstream sources, so this is a best-effort integrity check, not a + hash). + 4. Apply `post_process` (e.g. leave `.gz` files as-is where the pipeline + expects gzip, since several params like `sifts_file`/`interpro_xml` + are consumed compressed). For `derived` entries, `post_process` is the + join/filter logic described in step 2 (e.g. a small pandas function + registered per entry) rather than a simple decompress. + 5. Log source URL + resolved version/date actually fetched, into a + `download_manifest.log` written to `data_dir`, so a user can tell + later exactly what version of each upstream file they're running + against — useful for reproducing / debugging a specific pipeline run. +- At the end, print a summary table: fetched / derived / skipped / + manual-required, and exit non-zero only if a required (non-manual) file + failed. + +## 3a. New source: PDBe-KB / PDBeChem v2 ligand-function data (for cofactor coverage) + +This is *not* one of the ~20 files currently required by the pipeline — +it's a new source needed only once +[`docs/cofactor_coverage_plan.md`](cofactor_coverage_plan.md) is +implemented, added here so the manifest stays the single place tracking +every external file the pipeline (present or planned) depends on. + +`docs/installation.md` already lists one PDBe-KB file +(`assemblies_data.csv`, from `pdbe-kb/complexes/`) but that's a different +resource — assembly preference data, unrelated to ligand function. The +relevant PDBe-KB/PDBeChem v2 resource for cofactor coverage is separate: +`https://ftp.ebi.ac.uk/pub/databases/msd/pdbechem_v2/additional_data/pdb_ligand_interactions/`, +containing `interacting_chains_with_ligand_functions.tsv` (202M) and +`pdb_bound_molecules.tsv` (137M). These are bulk, structure/bound-entity +level outputs of PDBe's own **RelLig** pipeline, which classifies every +observed PDB ligand as reactant-like, cofactor-like, or drug-like using +the same PARITY method ProCogGraph already implements — see +`cofactor_coverage_plan.md` for why this is the primary source to use for +that work, ahead of the UniProt-REST approach originally proposed there. +- `source_type: direct_url`, two files, no per-accession chunking needed + (unlike the UniProt REST approach this replaces as primary source). +- Given the file sizes (well over the small reference tables elsewhere in + this manifest — 200M+ each), the size-sanity check described in the + script design below (§3) should use a floor generous enough not to + false-flag on minor version growth, but tight enough to catch a + truncated/failed download. + +## 4. Integration points + +- `docs/installation.md` — replace the current manual download table + (lines 99-124) with instructions to run the new script, keeping a + smaller reference table of what it fetches (generated from the manifest, + see step 1) for transparency. No manual-step instructions should be + needed once `ChEBI_Results.tsv` is handled as a `derived` entry. +- Do **not** wire this into `main.nf` as a pipeline process — reference + data changes far less often than pipeline runs, and folding it into the + Nextflow DAG would force a re-check/re-download on every run. Keep it as + a separate, explicitly-invoked setup step, analogous to + `setup_docker_linux.sh` for the Docker path. +- Add the new script's dependencies (`requests`, `pyyaml` if not already + present) to `nextflow/envs/procoggraph.yaml` if missing. + +## 5. Verification + +- `--dry-run` against the manifest should list exactly the files currently + documented in `docs/installation.md`'s table, plus `pfamA.txt.gz`, minus + `clan_membership.txt.gz`/`clan.txt.gz` (replaced by the single + `Pfam-A.clans.tsv.gz`), and with `chebi_names.tsv.gz`/`relation.tsv` + pointing at their corrected `flat_files/names.tsv.gz` / + `flat_files/relation.tsv.gz` locations (see §2 above). +- A real run into a clean `data_dir` should produce a directory usable + directly by `nextflow run main.nf` (module fix_plan items 1-2 aside), + i.e. every `params.*_file`/`*_dir` referenced in + `nextflow.config:86-102` (plus `scop2_descriptions_file` once added) + resolves to an existing, non-empty file. +- Re-running with no flags should skip everything (fast, idempotent); + re-running with `--force` should refresh everything and update + `download_manifest.log`. + +## Suggested order of work + +1. Build the manifest (step 1) — this alone surfaces/documents the + `pfamA.txt.gz` and `scop2_descriptions_file` gaps precisely, which is + useful even before the script exists. +2. Implement the script for the direct-URL files first (majority of the + list, lowest risk). +3. Add "latest release" resolution for CATH/SCOP2/Pfam/InterPro. +4. Add the `derived` handling for `ChEBI_Results.tsv` (fetch + `database_accession.tsv` + `structures.csv.gz`, join, write output) and + validate its output against a copy of the original manually-exported + file if one is still available, to catch any subtle format/column + differences before relying on it. +5. Update `docs/installation.md` to point at the script. + +This plan only covers *downloading* the reference data — it doesn't touch +the pipeline wiring bugs already tracked in +[`docs/fix_plan.md`](fix_plan.md). Fixing fix_plan.md item 2 +(`scop2_descriptions_file` missing from `nextflow.config`) can happen +independently, but should land before or alongside this work so the +manifest's `param_name` cross-check has something to check against. diff --git a/nextflow/bin/download_reference_data.py b/nextflow/bin/download_reference_data.py new file mode 100644 index 0000000..654e68a --- /dev/null +++ b/nextflow/bin/download_reference_data.py @@ -0,0 +1,350 @@ +#!/usr/bin/env python + +""" +Downloads and validates the external reference data files needed to run +the ProCogGraph Nextflow pipeline from scratch (see docs/installation.md +and nextflow/nextflow.config:85-103), driven by +reference_data_manifest.yaml. See docs/reference_data_download_plan.md +for the design rationale. + +Example usage: + python3 download_reference_data.py --data_dir /path/to/data_dir + python3 download_reference_data.py --data_dir /path/to/data_dir --dry-run + python3 download_reference_data.py --data_dir /path/to/data_dir --only enzyme_dat,rhea2ec + python3 download_reference_data.py --data_dir /path/to/data_dir --force +""" + +import argparse +import csv +import gzip +import shutil +import sys +import tarfile +import time +from datetime import datetime, timezone +from pathlib import Path + +import pandas as pd +import requests +import yaml + +DEFAULT_MANIFEST = Path(__file__).parent / "reference_data_manifest.yaml" +RETRY_ATTEMPTS = 3 +RETRY_BACKOFF_SECONDS = 5 +REQUEST_TIMEOUT_SECONDS = 120 + + +def load_manifest(manifest_path): + with open(manifest_path) as handle: + manifest = yaml.safe_load(handle) + return manifest["entries"] + + +def download_with_retry(url, dest_path): + last_exception = None + for attempt in range(1, RETRY_ATTEMPTS + 1): + try: + with requests.get(url, stream=True, timeout=REQUEST_TIMEOUT_SECONDS) as response: + response.raise_for_status() + tmp_path = dest_path.with_suffix(dest_path.suffix + ".part") + with open(tmp_path, "wb") as out_file: + for chunk in response.iter_content(chunk_size=1024 * 1024): + out_file.write(chunk) + tmp_path.rename(dest_path) + return + except (requests.RequestException, OSError) as exc: + last_exception = exc + if attempt < RETRY_ATTEMPTS: + print(f" download failed (attempt {attempt}/{RETRY_ATTEMPTS}): {exc} - retrying in {RETRY_BACKOFF_SECONDS}s") + time.sleep(RETRY_BACKOFF_SECONDS) + raise RuntimeError(f"failed to download {url} after {RETRY_ATTEMPTS} attempts") from last_exception + + +def check_size(dest_path, entry): + """For file entries, checks entry['min_size_bytes']. For directory + entries (post_process: extract_tar_gz), checks entry['min_file_count'] + instead, since a directory's own stat().st_size isn't meaningful.""" + if dest_path.is_dir(): + min_file_count = entry.get("min_file_count") + if min_file_count is None: + return + actual_count = sum(1 for _ in dest_path.iterdir()) + if actual_count < min_file_count: + raise RuntimeError( + f"{dest_path.name}/ contains {actual_count} files, below the expected " + f"minimum of {min_file_count} - likely a truncated/partial extraction. " + f"Not treating this as a successful fetch." + ) + return + + min_size_bytes = entry.get("min_size_bytes") + if min_size_bytes is None: + return + actual_size = dest_path.stat().st_size + if actual_size < min_size_bytes: + raise RuntimeError( + f"{dest_path.name} is {actual_size} bytes, below the expected minimum of " + f"{min_size_bytes} bytes - likely a truncated download, an error page, or a " + f"stale URL. Not treating this as a successful fetch." + ) + + +def post_process_gunzip(dest_path): + gz_path = dest_path.with_name(dest_path.name + ".gz.tmp") + dest_path.rename(gz_path) + with gzip.open(gz_path, "rb") as src, open(dest_path, "wb") as dst: + shutil.copyfileobj(src, dst) + gz_path.unlink() + + +def post_process_extract_tar_gz(dest_path): + """For entries whose target_filename is a directory (e.g. Rhea's rd/ + folder of per-reaction files): the tarball downloads to dest_path as a + plain file first (same as every other entry), then this extracts it. + + Assumes the tarball's own top-level entry is a directory with the same + name as target_filename (verified true for rhea-rd.tar.gz, which + contains rd/10001.rd etc., not bare 10001.rd) - so extraction target is + dest_path's parent, letting the tarball's own top-level folder become + dest_path. Extracting *into* dest_path itself would double-nest + (dest_path/rd/10001.rd instead of dest_path/10001.rd).""" + tar_path = dest_path.with_name(dest_path.name + ".tar.gz.tmp") + dest_path.rename(tar_path) + with tarfile.open(tar_path, "r:gz") as tar: + tar.extractall(path=dest_path.parent, filter="data") + tar_path.unlink() + if not dest_path.is_dir(): + raise RuntimeError( + f"expected extracting {tar_path.name} to produce a directory at {dest_path}, " + f"but it didn't - the tarball's internal layout may have changed." + ) + + +def post_process_csv_to_tsv_gz(dest_path): + """SIFTS' only live source is comma-separated and uncompressed, but the + pipeline (download_mmcif.py) reads it tab-separated from a .gz file.""" + csv_path = dest_path.with_name(dest_path.name + ".csv.tmp") + dest_path.rename(csv_path) + with open(csv_path, newline="") as src, gzip.open(dest_path, "wt", newline="") as dst: + reader = csv.reader(src) + writer = csv.writer(dst, delimiter="\t") + for row in reader: + writer.writerow(row) + csv_path.unlink() + + +def derive_chebi_results(data_dir, entry, manifest_by_name, force, skip_existing): + """Reconstruct ChEBI_Results.tsv (ChEBI ID/name/SMILES/KEGG COMPOUND + ACCESSION for ChEBI entries with a KEGG cross-reference and a known + structure) from bulk ChEBI flat files, instead of ChEBI's advanced + search UI. See docs/reference_data_download_plan.md #2. + + Schema notes (verified against live ChEBI flat_files, 2026-08 - these + are NOT documented anywhere obvious, so don't "simplify" this without + re-checking against the real files first): + - all ChEBI flat file columns are lowercase. + - database_accession.tsv.gz's `type` column is a mix of CAS/CITATION/ + MANUAL_X_REF/REGISTRY_NUMBER - it does NOT contain a KEGG-specific + type value. KEGG COMPOUND cross-references are `type == + "MANUAL_X_REF"` rows whose `source_id` is 45 in source.tsv.gz (not + downloaded separately here - confirmed by inspecting source.tsv.gz + directly; a stable ChEBI-internal primary key, not expected to + change). Only KEGG COMPOUND is used here, not KEGG GLYCAN (47) or + KEGG DRUG (46), to match get_ec_information.py's actual merge key + (nextflow/bin/get_ec_information.py:536, `right_on = "KEGG COMPOUND + ACCESSION"` - the original manual export also had separate + KEGG GLYCAN/KEGG DRUG columns, but neither is ever read). + - structures.tsv.gz has a `smiles` column directly (one row per + compound_id, not melted by structure type) - no TYPE=="SMILES" + filtering needed or possible. + - the "ID" column MUST be formatted as "CHEBI:" (not a bare + integer) - get_ec_information.py:543 sets `ligand_db` directly + from this column, and :682 later does + `ligand_db.str.extractall("CHEBI:([0-9]+)")` to feed ChEBI + has_role-based cofactor tagging. A bare integer silently produces + zero regex matches there rather than an error. + """ + KEGG_COMPOUND_SOURCE_ID = 45 + + dependency_names = ["chebi_names", "chebi_database_accession", "chebi_structures"] + dependency_paths = {} + for dep_name in dependency_names: + dep_entry = manifest_by_name[dep_name] + dep_path = data_dir / dep_entry["target_filename"] + if not dep_path.exists(): + print(f" dependency {dep_name} not present, fetching it first") + fetch_entry(data_dir, dep_entry, manifest_by_name, force, skip_existing) + dependency_paths[dep_name] = dep_path + + names_df = pd.read_csv(dependency_paths["chebi_names"], sep="\t", compression="gzip") + names_df = names_df.groupby("compound_id").agg({"name": "first"}).reset_index() + + accession_df = pd.read_csv(dependency_paths["chebi_database_accession"], sep="\t", compression="gzip") + kegg_accession_df = accession_df.loc[ + (accession_df["type"] == "MANUAL_X_REF") & (accession_df["source_id"] == KEGG_COMPOUND_SOURCE_ID) + ][["compound_id", "accession_number"]].rename(columns={"accession_number": "KEGG COMPOUND ACCESSION"}) + + structures_df = pd.read_csv(dependency_paths["chebi_structures"], sep="\t", compression="gzip") + smiles_df = structures_df.dropna(subset=["smiles"])[["compound_id", "smiles"]].drop_duplicates( + subset="compound_id", keep="first" + ).rename(columns={"smiles": "SMILES"}) + + merged = kegg_accession_df.merge(smiles_df, on="compound_id", how="inner") + merged = merged.merge(names_df, on="compound_id", how="left") + merged["ID"] = "CHEBI:" + merged["compound_id"].astype(str) + merged = merged.rename(columns={"name": "NAME"}) + merged = merged[["ID", "NAME", "SMILES", "KEGG COMPOUND ACCESSION"]].drop_duplicates() + + dest_path = data_dir / entry["target_filename"] + merged.to_csv(dest_path, sep="\t", index=False) + return dest_path + + +def derive_pubchem_kegg_mapping(data_dir, entry, manifest_by_name, force, skip_existing): + """Reconstruct pubchem_substance_id_mapping.txt (PubChem CID <-> KEGG + code pairs) from NCBI's bulk CID-Identifiers.tsv.gz cross-reference + table, instead of an interactive PubChem substance search export. + + get_ec_information.py's parser (nextflow/bin/get_ec_information.py: + 556-570) is a hand-rolled "KEY:VALUE KEY:VALUE ..." per-line format, + not a real TSV - each output line here is "CID: KEGG:", + matching that exactly. It needs at minimum CID and KEGG keys per line; + no other fields from the original manual export are read anywhere + downstream, so nothing else needs to be included. + """ + dep_entry = manifest_by_name["pubchem_cid_identifiers"] + dep_path = data_dir / dep_entry["target_filename"] + if not dep_path.exists(): + print(f" dependency pubchem_cid_identifiers not present, fetching it first") + fetch_entry(data_dir, dep_entry, manifest_by_name, force, skip_existing) + + cid_identifiers = pd.read_csv( + dep_path, sep="\t", header=None, compression="gzip", + names=["cid", "identifier", "identifier_type"], + ) + kegg_mapping = cid_identifiers.loc[cid_identifiers["identifier_type"] == "KEGG ID"] + + dest_path = data_dir / entry["target_filename"] + with open(dest_path, "w") as out_file: + for cid, kegg_code in zip(kegg_mapping["cid"], kegg_mapping["identifier"]): + out_file.write(f"CID:{cid} KEGG:{kegg_code}\n") + return dest_path + + +DERIVED_FUNCTIONS = { + "derive_chebi_results": derive_chebi_results, + "derive_pubchem_kegg_mapping": derive_pubchem_kegg_mapping, +} + +POST_PROCESS_FUNCTIONS = { + "gunzip": post_process_gunzip, + "csv_to_tsv_gz": post_process_csv_to_tsv_gz, + "extract_tar_gz": post_process_extract_tar_gz, +} + + +def fetch_entry(data_dir, entry, manifest_by_name, force, skip_existing): + name = entry["name"] + source_type = entry["source_type"] + dest_path = data_dir / entry["target_filename"] if entry.get("target_filename") else None + + if dest_path is not None and dest_path.exists() and skip_existing and not force: + print(f"[skip] {name}: {dest_path.name} already exists") + return "skipped", dest_path + + if source_type == "manual": + print(f"[manual] {name}: no automated source available.") + if entry.get("note"): + print(f" {entry['note'].strip()}") + return "manual", None + + if source_type == "needs_code_update": + print(f"[fetched, needs code update] {name}: downloading {dest_path.name}, but pipeline code isn't ready to consume it yet.") + if entry.get("note"): + print(f" {entry['note'].strip()}") + download_with_retry(entry["url"], dest_path) + check_size(dest_path, entry) + return "needs_code_update", dest_path + + if source_type == "derived": + print(f"[derive] {name}: deriving {dest_path.name}") + derive_fn = DERIVED_FUNCTIONS[entry["post_process"]] + derive_fn(data_dir, entry, manifest_by_name, force, skip_existing) + check_size(dest_path, entry) + return "derived", dest_path + + if source_type == "direct_url": + print(f"[fetch] {name}: downloading {entry['url']}") + download_with_retry(entry["url"], dest_path) + post_process = entry.get("post_process", "none") + if post_process != "none": + POST_PROCESS_FUNCTIONS[post_process](dest_path) + check_size(dest_path, entry) + return "fetched", dest_path + + raise ValueError(f"unknown source_type '{source_type}' for entry '{name}'") + + +def main(): + parser = argparse.ArgumentParser(description=__doc__, formatter_class=argparse.RawDescriptionHelpFormatter) + parser.add_argument("--data_dir", required=True, type=str, help="output directory (matches params.data_dir)") + parser.add_argument("--manifest", type=str, default=str(DEFAULT_MANIFEST), help="path to reference_data_manifest.yaml") + parser.add_argument("--only", type=str, default=None, help="comma-separated entry names to fetch (default: all)") + parser.add_argument("--force", action="store_true", help="re-download even if the target file already exists") + parser.add_argument("--include-optional", action="store_true", help="also fetch entries marked optional (e.g. cofactor-coverage data, not required by the current pipeline)") + parser.add_argument("--dry-run", action="store_true", help="print what would be fetched/skipped without downloading") + args = parser.parse_args() + + data_dir = Path(args.data_dir) + data_dir.mkdir(parents=True, exist_ok=True) + + entries = load_manifest(args.manifest) + manifest_by_name = {entry["name"]: entry for entry in entries} + + if args.only: + requested_names = set(args.only.split(",")) + unknown = requested_names - set(manifest_by_name) + if unknown: + print(f"error: unknown entry name(s): {', '.join(sorted(unknown))}", file=sys.stderr) + sys.exit(1) + entries = [entry for entry in entries if entry["name"] in requested_names] + elif not args.include_optional: + entries = [entry for entry in entries if not entry.get("optional", False)] + + if args.dry_run: + print(f"Would process {len(entries)} entries into {data_dir}:\n") + for entry in entries: + target = entry.get("target_filename") or "(derived, no fixed target)" + exists = (data_dir / entry["target_filename"]).exists() if entry.get("target_filename") else False + action = "skip (exists)" if exists and not args.force else entry["source_type"] + print(f" [{action:20s}] {entry['name']:35s} -> {target}") + return + + results = {} + log_lines = [f"# download_reference_data.py run at {datetime.now(timezone.utc).isoformat()}"] + for entry in entries: + try: + status, dest_path = fetch_entry(data_dir, entry, manifest_by_name, args.force, skip_existing=not args.force) + except Exception as exc: + print(f"[FAILED] {entry['name']}: {exc}") + status, dest_path = "failed", None + results.setdefault(status, []).append(entry["name"]) + source_note = entry.get("url") or "(derived/manual)" + log_lines.append(f"{entry['name']}\t{status}\t{source_note}\t{dest_path or ''}") + + log_path = data_dir / "download_manifest.log" + with open(log_path, "a") as log_file: + log_file.write("\n".join(log_lines) + "\n") + + print("\nSummary:") + for status in ["fetched", "derived", "skipped", "needs_code_update", "manual", "failed"]: + names = results.get(status, []) + if names: + print(f" {status}: {len(names)} ({', '.join(names)})") + + if results.get("failed"): + sys.exit(1) + + +if __name__ == "__main__": + main() diff --git a/nextflow/bin/process_all_pdb_contacts.py b/nextflow/bin/process_all_pdb_contacts.py index b65ac16..b47a052 100644 --- a/nextflow/bin/process_all_pdb_contacts.py +++ b/nextflow/bin/process_all_pdb_contacts.py @@ -106,9 +106,7 @@ def main(): parser = argparse.ArgumentParser(description='') parser.add_argument('--ccd_cif', type=str, help='cif file containing chemical component dictionary in mmcif format') parser.add_argument('--contacts_file', type=str, help='csv file containing concatenated contacts from pdbe-arpeggio for all pdbs') - parser.add_argument('--pfam_a_file', type=str, help='pfam a file') - parser.add_argument('--pfam_clan_rels', type=str, help='pfam clan relationships file') - parser.add_argument('--pfam_clans', type=str, help='pfam clans file') + parser.add_argument('--pfam_clans_file', type=str, help='consolidated Pfam-A.clans.tsv.gz file (pfam accession/name/description + clan membership)') parser.add_argument('--scop_domains_info_file', type=str, help='scop domains info file') parser.add_argument('--scop_descriptions_file', type=str, help='scop descriptions file') parser.add_argument('--scop2_domains_info_file', type=str, help='scop2 domains info file') @@ -299,7 +297,7 @@ def main(): scop_contacts.to_csv(f"scop_pdb_residue_interactions.csv.gz", sep = "\t", index = False, compression = "gzip") if len(pfam_contacts) > 0: - pfam_annotations = get_pfam_annotations(args.pfam_a_file, args.pfam_clan_rels, args.pfam_clans) + pfam_annotations = get_pfam_annotations(args.pfam_clans_file) pfam_contacts = pfam_contacts.merge(pfam_annotations, left_on = "xref_db_acc" , right_on = "pfam_accession", how = "left") pfam_contacts.to_csv(f"pfam_pdb_residue_interactions.csv.gz", sep = "\t", index = False, compression = "gzip") else: diff --git a/nextflow/bin/reference_data_manifest.yaml b/nextflow/bin/reference_data_manifest.yaml new file mode 100644 index 0000000..935f365 --- /dev/null +++ b/nextflow/bin/reference_data_manifest.yaml @@ -0,0 +1,354 @@ +# Manifest of external reference data required by the ProCogGraph Nextflow +# pipeline (see nextflow/nextflow.config:85-103 and docs/installation.md). +# This file is the single source of truth for "what does the pipeline need, +# and where does it come from" — consumed by +# nextflow/bin/download_reference_data.py. See docs/reference_data_download_plan.md +# for the rationale behind each entry, especially the ones marked with a +# `confidence` or `note` field below (URLs checked live during planning are +# marked confidence: verified; a few could not be confirmed and are marked +# confidence: unverified — the download script's dry-run/size-check is the +# safety net for those, not blind trust in this file). +# +# Fields: +# name - unique identifier, used with --only +# target_filename- filename written into --data_dir +# param_name - matching nextflow.config param (for cross-checking), or +# null if the file is only used outside nextflow.config +# (e.g. as a CLI arg to download_mmcif.py) +# source_type - direct_url | derived | manual | needs_code_update +# url - source URL (direct_url entries) +# min_size_bytes - best-effort sanity floor for the downloaded file +# post_process - none | gunzip | csv_to_tsv_gz | +# confidence - verified (checked live while building this manifest) | +# unverified (best-known convention, not confirmed live) +# note - free text, e.g. why a source needs special handling + +entries: + - name: sifts_ec_mapping + target_filename: pdb_chain_enzyme.tsv.gz + param_name: sifts_file + source_type: direct_url + url: https://ftp.ebi.ac.uk/pub/databases/msd/sifts/csv/pdb_chain_enzyme.csv + min_size_bytes: 1000000 + post_process: csv_to_tsv_gz + confidence: verified + note: > + The only live SIFTS source is comma-separated and uncompressed, but + the pipeline (download_mmcif.py:182) reads it with sep="\t" from a + .tsv.gz file. post_process converts comma -> tab and gzips it. + + - name: pdbe_kb_assemblies + target_filename: assemblies_data.csv + param_name: null + source_type: direct_url + url: https://ftp.ebi.ac.uk/pub/databases/pdbe-kb/complexes/assemblies_data.csv + min_size_bytes: 10000000 + post_process: none + confidence: verified + note: > + Not a nextflow.config param - consumed directly as a CLI arg by + download_mmcif.py --assemblies_file when building the structure + manifest (docs/installation.md step 4). + + - name: enzyme_dat + target_filename: enzyme.dat + param_name: enzyme_dat_file + source_type: direct_url + url: https://ftp.expasy.org/databases/enzyme/enzyme.dat + min_size_bytes: 5000000 + post_process: none + confidence: verified + + - name: enzyme_class + target_filename: enzclass.txt + param_name: enzyme_class_file + source_type: direct_url + url: https://ftp.expasy.org/databases/enzyme/enzclass.txt + min_size_bytes: 5000 + post_process: none + confidence: verified + + - name: scop_domains_info + target_filename: dir.cla.scop.1_75.txt + param_name: scop_domains_info_file + source_type: direct_url + url: https://scop.berkeley.edu/downloads/parse/dir.cla.scop.1.75.txt + min_size_bytes: 5000000 + post_process: none + confidence: verified + note: > + Confirmed live via https://scop.berkeley.edu/downloads/ver=1.75 + (2026-08). Real upstream filename uses a period (dir.cla.scop.1.75.txt, + 10.5M) not an underscore - dir.cla.scop.1_75.txt 404s. target_filename + keeps the underscore form to match params.scop_domains_info_file. + + - name: scop_descriptions + target_filename: dir.des.scop.1_75.txt + param_name: scop_descriptions_file + source_type: direct_url + url: https://scop.berkeley.edu/downloads/parse/dir.des.scop.1.75.txt + min_size_bytes: 1000000 + post_process: none + confidence: verified + note: > + Confirmed live alongside scop_domains_info (2026-08, 5.8M). Same + period-vs-underscore filename note applies. + + - name: scop2_domains_info + target_filename: scop2-cla-latest.txt + param_name: scop2_domains_info_file + source_type: direct_url + url: https://www.ebi.ac.uk/pdbe/scop/files/scop-cla-latest.txt + min_size_bytes: 1000000 + post_process: none + confidence: verified + note: > + Confirmed live via https://www.ebi.ac.uk/pdbe/scop/download (2026-08). + Note the real upstream filename is scop-cla-latest.txt (no "2") despite + this being the SCOP2 dataset - target_filename keeps the local + scop2-cla-latest.txt name to match params.scop2_domains_info_file. + + - name: scop2_descriptions + target_filename: scop2-des-latest.txt + param_name: scop2_descriptions_file + source_type: direct_url + url: https://www.ebi.ac.uk/pdbe/scop/files/scop-des-latest.txt + min_size_bytes: 100000 + post_process: none + confidence: verified + note: > + Confirmed live alongside scop2_domains_info (2026-08). Same "no 2 in + the upstream filename" note applies. + + - name: interpro_xml + target_filename: interpro.xml.gz + param_name: interpro_xml + source_type: direct_url + url: https://ftp.ebi.ac.uk/pub/databases/interpro/current_release/interpro.xml.gz + min_size_bytes: 20000000 + post_process: none + confidence: verified + + - name: pfam_clans + target_filename: Pfam-A.clans.tsv.gz + param_name: pfam_clans_file + source_type: direct_url + url: https://ftp.ebi.ac.uk/pub/databases/Pfam/current_release/Pfam-A.clans.tsv.gz + min_size_bytes: 100000 + post_process: none + confidence: verified + note: > + Pfam consolidated the three files this used to be (pfamA.txt.gz, + clan_membership.txt.gz, clan.txt.gz) into this single file (2026 + restructure). utils.get_pfam_annotations and + params.pfam_clans_file/nextflow.config updated to match - see + docs/reference_data_download_plan.md #2. No longer excluded from the + "ready to run pipeline" summary - fully wired up now. + + - name: cath_names + target_filename: cath-names.txt + param_name: cath_names + source_type: direct_url + url: http://download.cathdb.info/cath/releases/latest-release/cath-classification-data/cath-names.txt + min_size_bytes: 100000 + post_process: none + confidence: verified + + - name: cath_domain_description_file + target_filename: cath-domain-description-file.txt + param_name: cddf + source_type: direct_url + url: http://download.cathdb.info/cath/releases/latest-release/cath-classification-data/cath-domain-description-file.txt + min_size_bytes: 50000000 + post_process: none + confidence: verified + + - name: ccd + target_filename: ccd.cif + param_name: ccd_cif + source_type: direct_url + url: https://ftp.ebi.ac.uk/pub/databases/msd/pdbechem_v2/ccd/components.cif + min_size_bytes: 500000000 + post_process: none + confidence: verified + note: > + Confirmed live (2026-08) - already uncompressed (856,048,364 bytes, + ~816M via HEAD request), so no gunzip post_process needed. Large + download; min_size_bytes set well below the real size to avoid + false-flagging on minor version growth. + + - name: pubchem_cid_identifiers + target_filename: CID-Identifiers.tsv.gz + param_name: null + source_type: direct_url + url: https://ftp.ncbi.nlm.nih.gov/pubchem/Compound/Extras/CID-Identifiers.tsv.gz + min_size_bytes: 50000000 + post_process: none + confidence: verified + note: > + Bulk NCBI PubChem CID-to-external-identifier cross-reference table + (CID, identifier, identifier_type per row, no header; 12.7M rows, + 2026-08). Intermediate input for the pubchem_substance_mapping + derivation below - user-supplied source. + + - name: pubchem_substance_mapping + target_filename: pubchem_substance_id_mapping.txt + param_name: pubchem_mapping + source_type: derived + url: null + min_size_bytes: 500000 + post_process: derive_pubchem_kegg_mapping + confidence: verified + note: > + Derived from pubchem_cid_identifiers (filtered to + identifier_type == "KEGG ID", 26,933 of 12.7M rows, mixing KEGG + COMPOUND (C#####) and KEGG DRUG (D#####) codes - both valid, KEGG's + REST API handles both uniformly and get_ec_information.py doesn't + distinguish them) instead of the original manual PubChem substance + search export. Output format matches what get_ec_information.py's + hand-rolled "KEY:VALUE KEY:VALUE" line parser + (nextflow/bin/get_ec_information.py:556-570) expects: one line per + CID/KEGG pair, "CID: KEGG:". + + - name: rhea2ec + target_filename: rhea2ec.tsv + param_name: rhea2ec + source_type: direct_url + url: https://ftp.expasy.org/databases/rhea/tsv/rhea2ec.tsv + min_size_bytes: 50000 + post_process: none + confidence: verified + note: Confirmed live (2026-08). + + - name: rhea_directions + target_filename: rhea-directions.tsv + param_name: rhea_directions + source_type: direct_url + url: https://ftp.expasy.org/databases/rhea/tsv/rhea-directions.tsv + min_size_bytes: 50000 + post_process: none + confidence: verified + note: Confirmed live (2026-08). + + - name: rhea_reaction_smiles + target_filename: rhea-reaction-smiles.tsv + param_name: rhea_reactions_smiles + source_type: direct_url + url: https://ftp.expasy.org/databases/rhea/tsv/rhea-reaction-smiles.tsv + min_size_bytes: 1000000 + post_process: none + confidence: verified + note: Confirmed live (2026-08). + + - name: rhea_rd_dir + target_filename: rd + param_name: null + source_type: direct_url + url: https://ftp.expasy.org/databases/rhea/ctfiles/rhea-rd.tar.gz + min_size_bytes: null + min_file_count: 30000 + post_process: extract_tar_gz + confidence: verified + note: > + preprocess_rhea.py --rd_dir expects a directory of Rhea reaction (.rd) + files, one per unidirectional Rhea reaction. The rd/ directory itself + (https://ftp.expasy.org/databases/rhea/ctfiles/rd/) has ~37,100 + individual files and its own README explicitly says "please avoid + browsing the subdirectories" - use the sibling rhea-rd.tar.gz archive + instead (confirmed live, 2026-08, 24.8M), extracted into a local rd/ + directory via the extract_tar_gz post_process. target_filename here is + a directory name, not a single file. + + - name: chebi_names + target_filename: chebi_names.tsv.gz + param_name: null + source_type: direct_url + url: https://ftp.ebi.ac.uk/pub/databases/chebi/flat_files/names.tsv.gz + min_size_bytes: 5000000 + post_process: none + confidence: verified + note: > + ChEBI restructured their FTP layout in 2026; this replaces the old + Flat_file_tab_delimited/chebi_names.tsv.gz path referenced in + docs/installation.md. Renamed on download from upstream's names.tsv.gz + to chebi_names.tsv.gz, since preprocess_rhea.py --chebi_names + (nextflow/bin/preprocess_rhea.py:93) reads it directly under that + local filename with pd.read_csv(sep="\t", compression="gzip"). + + - name: chebi_relation + target_filename: relation.tsv + param_name: null + source_type: direct_url + url: https://ftp.ebi.ac.uk/pub/databases/chebi/flat_files/relation.tsv.gz + min_size_bytes: 500000 + post_process: gunzip + confidence: verified + note: > + Upstream file is now gzipped (relation.tsv.gz); get_ec_information.py + --chebi_relations expects a plain .tsv, so this is gunzipped on + download. + + - name: chebi_database_accession + target_filename: chebi_database_accession.tsv.gz + param_name: null + source_type: direct_url + url: https://ftp.ebi.ac.uk/pub/databases/chebi/flat_files/database_accession.tsv.gz + min_size_bytes: 500000 + post_process: none + confidence: verified + note: Intermediate input for the ChEBI_Results.tsv derivation below. + + - name: chebi_structures + target_filename: chebi_structures.tsv.gz + param_name: null + source_type: direct_url + url: https://ftp.ebi.ac.uk/pub/databases/chebi/flat_files/structures.tsv.gz + min_size_bytes: 50000000 + post_process: none + confidence: verified + note: > + Intermediate input for the ChEBI_Results.tsv derivation below. Large + (~88M compressed) - the size-sanity floor here is set generously + below the real size to avoid false-flagging on minor version growth. + + - name: chebi_results + target_filename: ChEBI_Results.tsv + param_name: null + source_type: derived + url: null + min_size_bytes: 500000 + post_process: derive_chebi_results + confidence: verified + note: > + Derived from chebi_names, chebi_database_accession and + chebi_structures (all fetched above) instead of ChEBI's advanced + search UI - see docs/reference_data_download_plan.md #2. Requires + those three entries to have been fetched first; the script fetches + them automatically as dependencies if missing. + + # --- Not required by the pipeline today - only needed once + # docs/cofactor_coverage_plan.md is implemented. Kept here so this + # manifest stays the single place tracking every external file the + # pipeline depends on, present or planned. Excluded from --dry-run's + # "ready to run pipeline" summary; fetch explicitly with --only. + - name: pdbe_rellig_chain_ligand_functions + target_filename: interacting_chains_with_ligand_functions.tsv + param_name: null + source_type: direct_url + url: https://ftp.ebi.ac.uk/pub/databases/msd/pdbechem_v2/additional_data/pdb_ligand_interactions/interacting_chains_with_ligand_functions.tsv + min_size_bytes: 100000000 + post_process: none + confidence: verified + optional: true + note: For the cofactor-coverage plan (docs/cofactor_coverage_plan.md), not the current pipeline. + + - name: pdbe_rellig_bound_molecules + target_filename: pdb_bound_molecules.tsv + param_name: null + source_type: direct_url + url: https://ftp.ebi.ac.uk/pub/databases/msd/pdbechem_v2/additional_data/pdb_ligand_interactions/pdb_bound_molecules.tsv + min_size_bytes: 50000000 + post_process: none + confidence: verified + optional: true + note: For the cofactor-coverage plan (docs/cofactor_coverage_plan.md), not the current pipeline. diff --git a/nextflow/bin/utils.py b/nextflow/bin/utils.py index 06e0aca..fbab7ff 100644 --- a/nextflow/bin/utils.py +++ b/nextflow/bin/utils.py @@ -229,16 +229,29 @@ def clean_and_merge_scop_col(df, column_id, description_df): scop_domains_info.drop(columns = ["pdb_id", "scop_description"], inplace = True) return scop_domains_info -def get_pfam_annotations(pfam_a_file, clan_membership_file, clan_info_file): - pfam_a = pd.read_csv(pfam_a_file, sep = "\t", header = None, usecols = [0,1,3], names = ["pfam_accession", "pfam_name", "pfam_description"]) - pfam_clan_rels = pd.read_csv(clan_membership_file, sep = "\t", header = None, names = ["clan", "pfam"]) - pfam_clans = pd.read_csv(clan_info_file, sep = "\t", comment = "#", header = None, names = ["clan_acc", "clan_id", "previous_id", "clan_description", "clan_author", "deposited_by", "clan_comment", "updated", "created", "version", "number_structures", "number_archs", "number_species", "number_sequences", "competed", "uniprot_competed"]) - pfam_clan_df = pfam_clan_rels.merge(pfam_clans[["clan_acc", "clan_description", "clan_comment"]], left_on = "clan", right_on = "clan_acc", how = "left", indicator = True) - assert(len(pfam_clan_df.loc[pfam_clan_df._merge != "both"]) == 0) - pfam_clan_df.drop(columns = "_merge", inplace = True) - - pfam_a_clans_merged = pfam_a.merge(pfam_clan_df, left_on = "pfam_accession", right_on = "pfam", how = "left") - return pfam_a_clans_merged +def get_pfam_annotations(pfam_clans_file): + """Pfam consolidated the three files this used to read separately - + pfamA.txt.gz (pfam accession/name/description), clan_membership.txt.gz + (clan-to-pfam membership) and clan.txt.gz (clan details) - into a + single Pfam-A.clans.tsv.gz (2026 restructure, no header row): + pfam_accession, clan_acc, clan_id, pfam_name, pfam_description. Pfam + families with no clan have empty clan_acc/clan_id fields (parsed as + NaN by pandas). See docs/reference_data_download_plan.md #2. + + Schema note: the consolidated file has no equivalent of the old + clan.txt.gz's free-text `clan_comment` field - Pfam no longer + publishes a bulk clan comment/description beyond the short clan id + (e.g. "GPCR_A"). clan_description below is that short id; + clan_comment is left null rather than duplicating clan_description + into it, since that would misrepresent it as distinct information. + """ + pfam_clans = pd.read_csv(pfam_clans_file, sep = "\t", header = None, + names = ["pfam_accession", "clan_acc", "clan_id", "pfam_name", "pfam_description"]) + pfam_clans["clan"] = pfam_clans["clan_acc"] + pfam_clans["pfam"] = pfam_clans["pfam_accession"].where(pfam_clans["clan_acc"].notna()) + pfam_clans["clan_description"] = pfam_clans["clan_id"] + pfam_clans["clan_comment"] = np.nan + return pfam_clans[["pfam_accession", "pfam_name", "pfam_description", "clan", "pfam", "clan_acc", "clan_description", "clan_comment"]] def return_partial_EC_list(ec, total_ec_list): if not isinstance(ec, str) and np.isnan(ec): diff --git a/nextflow/envs/procoggraph.yaml b/nextflow/envs/procoggraph.yaml index bff7cee..27b07bf 100644 --- a/nextflow/envs/procoggraph.yaml +++ b/nextflow/envs/procoggraph.yaml @@ -16,4 +16,5 @@ dependencies: - rdkit - biopython - beautifulsoup4 - - pubchempy \ No newline at end of file + - pubchempy + - pyyaml \ No newline at end of file diff --git a/nextflow/main.nf b/nextflow/main.nf index ebab36c..9e520e6 100644 --- a/nextflow/main.nf +++ b/nextflow/main.nf @@ -65,9 +65,7 @@ process PROCESS_ALL_CONTACTS { input: path combined_contacts path ccd_cif - path pfam_a - path pfam_clan_rels - path pfam_clans + path pfam_clans_file path scop_domains_info path scop_domains_description path scop2_domains_info @@ -95,7 +93,7 @@ process PROCESS_ALL_CONTACTS { script: """ - python3 ${workflow.projectDir}/bin/process_all_pdb_contacts.py --contacts_file ${combined_contacts} --ccd_cif ${ccd_cif} --pfam_a_file ${pfam_a} --pfam_clan_rels ${pfam_clan_rels} --pfam_clans ${pfam_clans} --scop_domains_info_file ${scop_domains_info} --scop_descriptions_file ${scop_domains_description} --scop2_domains_info_file ${scop2_domains_info} --scop2_descriptions_file ${scop2_domains_description} --interpro_xml ${interpro_xml} --cath_names ${cath_names} --cddf ${cddf} --glycoct_cache ${glycoct_cache} --smiles_cache ${smiles_cache} --csdb_linear_cache ${csdb_linear_cache} --enzyme_dat_file ${enzyme_dat_file} --enzyme_class_file ${enzyme_class_file} --sifts_ec_mapping ${sifts_ec_mapping} + python3 ${workflow.projectDir}/bin/process_all_pdb_contacts.py --contacts_file ${combined_contacts} --ccd_cif ${ccd_cif} --pfam_clans_file ${pfam_clans_file} --scop_domains_info_file ${scop_domains_info} --scop_descriptions_file ${scop_domains_description} --scop2_domains_info_file ${scop2_domains_info} --scop2_descriptions_file ${scop2_domains_description} --interpro_xml ${interpro_xml} --cath_names ${cath_names} --cddf ${cddf} --glycoct_cache ${glycoct_cache} --smiles_cache ${smiles_cache} --csdb_linear_cache ${csdb_linear_cache} --enzyme_dat_file ${enzyme_dat_file} --enzyme_class_file ${enzyme_class_file} --sifts_ec_mapping ${sifts_ec_mapping} """ } @@ -227,7 +225,7 @@ workflow { arpeggio_batches = processed_struct_manifest.arpeggio_batch.flatten() arpeggio = RUN_ARPEGGIO( arpeggio_batches ).collect() contacts = PROCESS_CONTACTS( arpeggio, processed_struct_manifest.updated_manifest, params.domain_contact_cutoff ) - all_contacts = PROCESS_ALL_CONTACTS( contacts.combined_contacts, Channel.fromPath("${params.ccd_cif}"), Channel.fromPath("${params.pfam_a_file}"), Channel.fromPath("${params.pfam_clan_rels}"), Channel.fromPath("${params.pfam_clans}"), Channel.fromPath("${params.scop_domains_info_file}"), Channel.fromPath("${params.scop_descriptions_file}"), Channel.fromPath("${params.scop2_domains_info_file}"), Channel.fromPath("${params.scop2_descriptions_file}"), Channel.fromPath("${params.interpro_xml}"), Channel.fromPath("${params.cath_names}"), Channel.fromPath("${params.cddf}"), Channel.fromPath("${params.glycoct_cache}"), Channel.fromPath("${params.smiles_cache}"), Channel.fromPath("${params.csdb_linear_cache}"), Channel.fromPath("${params.enzyme_dat_file}"), Channel.fromPath("${params.enzyme_class_file}"), Channel.fromPath("${params.sifts_file}") ) + all_contacts = PROCESS_ALL_CONTACTS( contacts.combined_contacts, Channel.fromPath("${params.ccd_cif}"), Channel.fromPath("${params.pfam_clans_file}"), Channel.fromPath("${params.scop_domains_info_file}"), Channel.fromPath("${params.scop_descriptions_file}"), Channel.fromPath("${params.scop2_domains_info_file}"), Channel.fromPath("${params.scop2_descriptions_file}"), Channel.fromPath("${params.interpro_xml}"), Channel.fromPath("${params.cath_names}"), Channel.fromPath("${params.cddf}"), Channel.fromPath("${params.glycoct_cache}"), Channel.fromPath("${params.smiles_cache}"), Channel.fromPath("${params.csdb_linear_cache}"), Channel.fromPath("${params.enzyme_dat_file}"), Channel.fromPath("${params.enzyme_class_file}"), Channel.fromPath("${params.sifts_file}") ) score_ligands = SCORE_LIGANDS( all_contacts.bound_entities, Channel.fromPath(params.cognate_ligands), Channel.fromPath(params.parity_cache), Channel.from(params.parity_threshold) ) produce_neo4j_files = PRODUCE_NEO4J_FILES( score_ligands.all_parity_calcs, Channel.fromPath(params.cognate_ligands) , all_contacts.bound_entities, all_contacts.cath, all_contacts.scop, all_contacts.pfam, all_contacts.superfamily, all_contacts.g3dsa, all_contacts.scop2b_sf, all_contacts.scop2b_fa, Channel.fromPath("${params.enzyme_dat_file}"), Channel.fromPath("${params.enzyme_class_file}"), Channel.from(params.parity_threshold), Channel.fromPath("${params.rhea2ec}"), Channel.fromPath("${params.rhea_directions}"), Channel.fromPath("${params.rhea_reactions_smiles}"), Channel.fromPath("${params.scop2_domains_info_file}"), Channel.fromPath("${params.scop2_descriptions_file}") ) } \ No newline at end of file diff --git a/nextflow/nextflow.config b/nextflow/nextflow.config index a154192..a739880 100644 --- a/nextflow/nextflow.config +++ b/nextflow/nextflow.config @@ -91,9 +91,7 @@ params { scop2_domains_info_file = "${params.data_dir}/scop2-cla-latest.txt" scop2_descriptions_file = "${params.data_dir}/scop2-des-latest.txt" interpro_xml = "${params.data_dir}/interpro.xml.gz" - pfam_clans = "${params.data_dir}/clan.txt.gz" - pfam_clan_rels = "${params.data_dir}/clan_membership.txt.gz" - pfam_a_file = "${params.data_dir}/pfamA.txt.gz" + pfam_clans_file = "${params.data_dir}/Pfam-A.clans.tsv.gz" cath_names = "${params.data_dir}/cath-names.txt" cddf = "${params.data_dir}/cath-domain-description-file.txt" ccd_cif = "${params.data_dir}/ccd.cif"