|
| 1 | +from __future__ import annotations |
| 2 | + |
| 3 | +import argparse |
| 4 | +import gzip |
| 5 | +import hashlib |
| 6 | +import json |
| 7 | +import os |
| 8 | +from pathlib import Path |
| 9 | +import subprocess |
| 10 | +import sys |
| 11 | +import tempfile |
| 12 | +from typing import Any |
| 13 | + |
| 14 | +import pandas as pd |
| 15 | +import torch |
| 16 | + |
| 17 | +_THIS_FILE = Path(__file__).resolve() |
| 18 | +_REPO_ROOT = _THIS_FILE.parents[2] |
| 19 | +if str(_REPO_ROOT) not in sys.path: |
| 20 | + sys.path.insert(0, str(_REPO_ROOT)) |
| 21 | + |
| 22 | +from catpred.inference import PredictionRequest, run_prediction_pipeline |
| 23 | +from catpred.security import load_torch_artifact |
| 24 | + |
| 25 | +_VALID_AAS = set("ACDEFGHIKLMNPQRSTVWY") |
| 26 | +_TARGET_TO_PARAMETER = { |
| 27 | + "kcat": "kcat", |
| 28 | + "Km": "km", |
| 29 | + "km": "km", |
| 30 | + "ki": "ki", |
| 31 | + "Ki": "ki", |
| 32 | +} |
| 33 | +_PREDICTION_COLUMN = { |
| 34 | + "kcat": "Prediction_(s^(-1))", |
| 35 | + "km": "Prediction_(mM)", |
| 36 | + "ki": "Prediction_(mM)", |
| 37 | +} |
| 38 | + |
| 39 | + |
| 40 | +def _repo_root() -> Path: |
| 41 | + return _REPO_ROOT |
| 42 | + |
| 43 | + |
| 44 | +def _contains_model_checkpoints(path: Path) -> bool: |
| 45 | + return path.exists() and path.is_dir() and any(path.rglob("model.pt")) |
| 46 | + |
| 47 | + |
| 48 | +def _discover_checkpoint_root(repo_root: Path) -> Path: |
| 49 | + production_root = (repo_root / ".e2e-assets" / "pretrained" / "production").resolve() |
| 50 | + checkpoints_root = (repo_root / "checkpoints").resolve() |
| 51 | + if _contains_model_checkpoints(production_root): |
| 52 | + return production_root |
| 53 | + return checkpoints_root |
| 54 | + |
| 55 | + |
| 56 | +def _env_path(name: str, default: Path) -> Path: |
| 57 | + value = os.environ.get(name) |
| 58 | + return Path(value).resolve() if value else default.resolve() |
| 59 | + |
| 60 | + |
| 61 | +def _resolve_parameter(payload: dict[str, Any]) -> str: |
| 62 | + target = payload.get("target") |
| 63 | + if isinstance(target, str) and target in _TARGET_TO_PARAMETER: |
| 64 | + return _TARGET_TO_PARAMETER[target] |
| 65 | + |
| 66 | + params = payload.get("params") or {} |
| 67 | + kinetics_type = str(params.get("kinetics_type", "")).strip().upper() |
| 68 | + if kinetics_type == "KCAT": |
| 69 | + return "kcat" |
| 70 | + if kinetics_type == "KM": |
| 71 | + return "km" |
| 72 | + if kinetics_type == "KI": |
| 73 | + return "ki" |
| 74 | + raise RuntimeError(f"Unsupported target in payload: {target!r}") |
| 75 | + |
| 76 | + |
| 77 | +def _stable_seq_id(sequence: str) -> str: |
| 78 | + digest = hashlib.sha1(sequence.encode("utf-8")).hexdigest()[:16] |
| 79 | + return f"seq_{digest}" |
| 80 | + |
| 81 | + |
| 82 | +def _resolve_seq_ids(sequences: list[str], tools_path: Path, media_path: Path) -> list[str]: |
| 83 | + seqmap_cli = tools_path / "seqmap" / "main.py" |
| 84 | + seqmap_db = media_path / "sequence_info" / "seqmap.sqlite3" |
| 85 | + if not seqmap_cli.exists() or not seqmap_db.exists(): |
| 86 | + return [_stable_seq_id(sequence) for sequence in sequences] |
| 87 | + |
| 88 | + seqmap_python = os.environ.get("CATPRED_SEQMAP_PYTHON", sys.executable) |
| 89 | + payload = "\n".join(sequences) + "\n" |
| 90 | + cmd = [ |
| 91 | + seqmap_python, |
| 92 | + str(seqmap_cli), |
| 93 | + "--db", |
| 94 | + str(seqmap_db), |
| 95 | + "batch-get-or-create", |
| 96 | + "--stdin", |
| 97 | + ] |
| 98 | + proc = subprocess.run( |
| 99 | + cmd, |
| 100 | + input=payload, |
| 101 | + text=True, |
| 102 | + stdout=subprocess.PIPE, |
| 103 | + stderr=subprocess.PIPE, |
| 104 | + ) |
| 105 | + if proc.returncode != 0: |
| 106 | + raise RuntimeError( |
| 107 | + f"seqmap failed (rc={proc.returncode})\nSTDOUT:\n{proc.stdout}\nSTDERR:\n{proc.stderr}" |
| 108 | + ) |
| 109 | + seq_ids = [line.strip() for line in proc.stdout.splitlines() if line.strip()] |
| 110 | + if len(seq_ids) != len(sequences): |
| 111 | + raise RuntimeError(f"seqmap returned {len(seq_ids)} ids for {len(sequences)} sequences") |
| 112 | + return seq_ids |
| 113 | + |
| 114 | + |
| 115 | +def _load_or_compute_esm2_feature(sequence: str, cache_file: Path) -> Path: |
| 116 | + cache_file.parent.mkdir(parents=True, exist_ok=True) |
| 117 | + if cache_file.exists(): |
| 118 | + load_torch_artifact( |
| 119 | + cache_file, |
| 120 | + purpose="CatPred ESM2 cache entry", |
| 121 | + map_location="cpu", |
| 122 | + roots=[cache_file.parent], |
| 123 | + ) |
| 124 | + return cache_file |
| 125 | + |
| 126 | + os.environ.setdefault("PROTEIN_EMBED_USE_CPU", "1") |
| 127 | + from catpred.data.esm_utils import get_single_esm_repr |
| 128 | + |
| 129 | + embedding = get_single_esm_repr(sequence).cpu() |
| 130 | + torch.save(embedding, cache_file) |
| 131 | + return cache_file |
| 132 | + |
| 133 | + |
| 134 | +def _build_input_dataframe(rows: list[dict[str, Any]], seq_ids: list[str]) -> pd.DataFrame: |
| 135 | + formatted_rows = [] |
| 136 | + for row, seq_id in zip(rows, seq_ids): |
| 137 | + substrate = row.get("substrates", row.get("substrate", row.get("Substrate", ""))) |
| 138 | + if isinstance(substrate, list): |
| 139 | + if len(substrate) != 1: |
| 140 | + raise RuntimeError("CatPred expects exactly one substrate per row.") |
| 141 | + substrate = substrate[0] |
| 142 | + substrate = str(substrate).strip() |
| 143 | + sequence = str(row.get("sequence", "")).strip() |
| 144 | + formatted_rows.append( |
| 145 | + { |
| 146 | + "SMILES": substrate, |
| 147 | + "sequence": sequence, |
| 148 | + "pdbpath": seq_id, |
| 149 | + } |
| 150 | + ) |
| 151 | + return pd.DataFrame(formatted_rows) |
| 152 | + |
| 153 | + |
| 154 | +def _write_protein_records( |
| 155 | + rows: list[dict[str, Any]], |
| 156 | + seq_ids: list[str], |
| 157 | + parameter: str, |
| 158 | + media_path: Path, |
| 159 | + out_path: Path, |
| 160 | +) -> None: |
| 161 | + records: dict[str, dict[str, Any]] = {} |
| 162 | + needs_esm = parameter in {"kcat", "km"} |
| 163 | + esm_cache_dir = media_path / "sequence_info" / "esm2_last" / "per_residue" |
| 164 | + |
| 165 | + for row, seq_id in zip(rows, seq_ids): |
| 166 | + sequence = str(row.get("sequence", "")).strip() |
| 167 | + record: dict[str, Any] = {"name": seq_id, "seq": sequence} |
| 168 | + if needs_esm: |
| 169 | + cache_file = _load_or_compute_esm2_feature(sequence, esm_cache_dir / f"{seq_id}.pt") |
| 170 | + record["esm2_feats_path"] = str(cache_file.resolve()) |
| 171 | + records[seq_id] = record |
| 172 | + |
| 173 | + with gzip.open(out_path, "wt", encoding="utf-8") as handle: |
| 174 | + json.dump(records, handle) |
| 175 | + |
| 176 | + |
| 177 | +def run_from_payload(payload: dict[str, Any]) -> dict[str, Any]: |
| 178 | + rows = payload.get("rows") or [] |
| 179 | + if not isinstance(rows, list): |
| 180 | + raise RuntimeError("'rows' must be a list in input payload.") |
| 181 | + |
| 182 | + valid_rows: list[dict[str, Any]] = [] |
| 183 | + invalid_indices: list[int] = [] |
| 184 | + for idx, row in enumerate(rows): |
| 185 | + sequence = str(row.get("sequence", "")).strip() |
| 186 | + substrate = row.get("substrates", row.get("substrate", row.get("Substrate", ""))) |
| 187 | + if isinstance(substrate, list): |
| 188 | + substrate = substrate[0] if len(substrate) == 1 else "" |
| 189 | + substrate = str(substrate).strip() |
| 190 | + if not sequence or not substrate or not set(sequence).issubset(_VALID_AAS): |
| 191 | + invalid_indices.append(idx) |
| 192 | + continue |
| 193 | + valid_rows.append(row) |
| 194 | + |
| 195 | + predictions: list[float | None] = [None] * len(rows) |
| 196 | + if not valid_rows: |
| 197 | + for idx in range(len(rows)): |
| 198 | + print(f"Progress: {idx + 1}/{len(rows)}", flush=True) |
| 199 | + return {"predictions": predictions, "invalid_indices": invalid_indices} |
| 200 | + |
| 201 | + repo_root = _env_path("CATPRED_REPO_ROOT", _repo_root()) |
| 202 | + media_path = _env_path("CATPRED_MEDIA_PATH", repo_root / "media") |
| 203 | + tools_path = _env_path("CATPRED_TOOLS_PATH", repo_root / "tools") |
| 204 | + checkpoint_root = _env_path("CATPRED_CHECKPOINT_ROOT", _discover_checkpoint_root(repo_root)) |
| 205 | + parameter = _resolve_parameter(payload) |
| 206 | + |
| 207 | + seq_ids = _resolve_seq_ids( |
| 208 | + [str(row.get("sequence", "")).strip() for row in valid_rows], |
| 209 | + tools_path=tools_path, |
| 210 | + media_path=media_path, |
| 211 | + ) |
| 212 | + |
| 213 | + with tempfile.TemporaryDirectory(prefix="catpred_webkinpred_") as tmp_dir_str: |
| 214 | + tmp_dir = Path(tmp_dir_str).resolve() |
| 215 | + input_csv = tmp_dir / "input.csv" |
| 216 | + protein_records = tmp_dir / "protein_records.json.gz" |
| 217 | + results_dir = tmp_dir / "results" |
| 218 | + |
| 219 | + _build_input_dataframe(valid_rows, seq_ids).to_csv(input_csv, index=False) |
| 220 | + _write_protein_records( |
| 221 | + rows=valid_rows, |
| 222 | + seq_ids=seq_ids, |
| 223 | + parameter=parameter, |
| 224 | + media_path=media_path, |
| 225 | + out_path=protein_records, |
| 226 | + ) |
| 227 | + |
| 228 | + request = PredictionRequest( |
| 229 | + parameter=parameter, |
| 230 | + input_file=str(input_csv), |
| 231 | + checkpoint_dir=str((checkpoint_root / parameter).resolve()), |
| 232 | + use_gpu=False, |
| 233 | + repo_root=str(repo_root), |
| 234 | + python_executable=sys.executable, |
| 235 | + protein_records_file=str(protein_records), |
| 236 | + ) |
| 237 | + output_file = run_prediction_pipeline(request=request, results_dir=str(results_dir)) |
| 238 | + output_df = pd.read_csv(output_file) |
| 239 | + value_col = _PREDICTION_COLUMN[parameter] |
| 240 | + if value_col not in output_df.columns: |
| 241 | + raise RuntimeError(f"CatPred output is missing expected column: {value_col}") |
| 242 | + |
| 243 | + valid_predictions = output_df[value_col].tolist() |
| 244 | + if len(valid_predictions) != len(valid_rows): |
| 245 | + raise RuntimeError( |
| 246 | + f"CatPred produced {len(valid_predictions)} predictions for {len(valid_rows)} rows." |
| 247 | + ) |
| 248 | + |
| 249 | + valid_iter = iter(valid_predictions) |
| 250 | + for idx in range(len(rows)): |
| 251 | + if idx in invalid_indices: |
| 252 | + continue |
| 253 | + predictions[idx] = float(next(valid_iter)) |
| 254 | + |
| 255 | + total = len(rows) |
| 256 | + for idx in range(total): |
| 257 | + print(f"Progress: {idx + 1}/{total}", flush=True) |
| 258 | + |
| 259 | + return { |
| 260 | + "predictions": predictions, |
| 261 | + "invalid_indices": sorted(set(invalid_indices)), |
| 262 | + } |
| 263 | + |
| 264 | + |
| 265 | +def main() -> None: |
| 266 | + parser = argparse.ArgumentParser(description="Run CatPred via the webKinPred subprocess contract.") |
| 267 | + parser.add_argument("--input", required=True, help="Input JSON path.") |
| 268 | + parser.add_argument("--output", required=True, help="Output JSON path.") |
| 269 | + args = parser.parse_args() |
| 270 | + |
| 271 | + with open(args.input, "r", encoding="utf-8") as handle: |
| 272 | + payload = json.load(handle) |
| 273 | + |
| 274 | + result = run_from_payload(payload) |
| 275 | + |
| 276 | + with open(args.output, "w", encoding="utf-8") as handle: |
| 277 | + json.dump(result, handle) |
| 278 | + |
| 279 | + |
| 280 | +if __name__ == "__main__": |
| 281 | + try: |
| 282 | + main() |
| 283 | + except Exception as exc: |
| 284 | + print(f"[CatPred] ERROR: {exc}", file=sys.stderr, flush=True) |
| 285 | + raise |
0 commit comments