Skip to content

DETAILS.md

Latest commit

 

History

History
159 lines (142 loc) · 11.9 KB

DETAILS.md

File metadata and controls

159 lines (142 loc) · 11.9 KB
Raw

chemtriage — standalone v22 density-sensitivity predictor · cluster + periodic (PBC)

Minimal deployable set to predict the DC-DFT density sensitivity S = |Σ_i c_i s_i| (kcal/mol) of a reaction from cheap GFN2-xTB descriptors, and to explain why. Model = the v22 XGBoost model (v21 pipeline retrained on the S^SC target; 104 spin-aware tblite GFN2-xTB features). v21 added a periodic (PBC) input path so the same descriptors can be computed for a crystal/surface cell — a research extrapolation (the model is trained on molecules), see Periodic mode below. The predicted S is mapped to a graded triage tier (not a hard yes/no), since a value of 2.01 is not the same call as 5.4 and the model's typical error is ~2 kcal/mol:

The tiers are calibrated on GMTKN55 (10×5-fold cross-validation, out-of-fold; 1505 reactions) — the last two columns are the actual fraction of reactions the model predicted in each band that turned out truly sensitive (true S>2) / strongly sensitive (true S>5):

predicted S (kcal/mol) tier truly S>2 truly S>5 meaning
< 1 SAFE 1% 0% run SC-DFT only
1–2 LIKELY SAFE 17% 2% usually fine, but ~1 in 6 is sensitive — spot-check
2–5 BE CAUTIOUS 64% 14% sensitive more often than not — SC-DFT often unreliable
≥ 5 BE VERY CAUTIOUS 95% 77% almost always sensitive — SC-DFT likely unreliable

So a prediction of 1.80 (17% real risk) is not the same call as 5.4 (95%). Each run prints its band's calibrated risk. Cut points + percentages live in BANDS at the top of predict.py (regenerate the % by binning ../y.npy by ../oof_pred.npy over the same edges) — edit for the paper.

Each run also reports an applicability-domain check — the mean standardized 5-NN distance from this reaction to the 1505 GMTKN55 training reactions, flagged inside / near-edge / outside the training domain, plus the 3 nearest known reactions. This operationalizes the LOSO caveat: outside the domain, trust the flag (AUC ~0.91) but not the exact value (r² ~0.42). Uses train_ref.npz; if that file is absent the check is silently skipped.

Layout

Fully self-contained / portable — no dependency on anything outside this folder (every path resolves relative to the folder; the only external requirement is the conda env in environment.yml). Copy it anywhere and run. The training pipeline that builds the data/ artifacts lives separately in ../train/.

chemtriage/
├── predict.py            CLI predictor + core (SHAP decomposition, triage tiers, applicability domain)
├── compute_gfn2_xtb.py   runtime GFN2-xTB engine (tblite; cluster + periodic descriptors; SMILES geom-opt)
├── webui.py              Flask web UI (examples tree, cluster/periodic tabs, hover previews)
├── molviz.py             ASE ball-and-stick renderer (incl. periodic unit-cell drawing)
├── environment.yml       conda env spec  → conda env create -f environment.yml
├── README.md
├── data/                 model + bundled reference data (all read-only)
│   ├── xgb_model.joblib      trained XGBoost model + feature list
│   ├── train_ref.npz         standardized training descriptors + kNN percentiles (applicability domain)
│   ├── refdb.pkl.gz          GMTKN55 geometries (ref/info/xyz for 1505 rxns) — cluster examples + hover images
│   ├── reference_S.dat       reference DC-DFT sensitivities S (shown vs the model's prediction)
│   ├── gmtkn_pred.json        precomputed model S for every GMTKN55 reaction
│   └── pbc_examples.json      the periodic systems from ../periodic (X23/ICE7/CO·NaCl/graphite/covalent)
└── web/                  browser assets
    ├── model_card.html       self-contained model card (metrics, ROC, calibration, importances) — /model-card
    ├── logo_tccl.png
    ├── brandserif.woff2 / brandserif-italic.woff2   title serif (Source Serif 4, OFL), served locally
    └── html2canvas.min.js    bundled DOM-screenshot lib (MIT) — powers the web "⤓ Save PNG" button

Regenerate the data/ artifacts with the ../train/ scripts (train.py → build_refdb.py / build_domain_ref.py / model_card.py). Runtime logs are written to a logs/ folder created on demand.

  • logs/ — one folder per run (created automatically), see below.

Setup (one command)

conda env create -f environment.yml     # creates env "dsxtb" with tblite + all python deps
conda activate dsxtb

This pins the exact stack it was built with (tblite 0.6.0, xgboost, flask, ase, …). Then jump to Usage / Web UI below. If you'd rather install into an existing env, see Requirements.

Requirements

  • the tblite python API — pin tblite-python=0.6.0 (conda install -c conda-forge tblite-python=0.6.0, or pip install tblite==0.6.0). It's a standalone GFN2-xTB library (the engine CP2K uses); it does not depend on Grimme's xtb program, so you do not need xtb/xtb-python installed. Portability note: the descriptor engine reads tblite result attributes (bond-orders, orbital-energies, dipole, …) whose names differ across tblite versions — a different tblite on another machine is the usual cause of "only N/104 descriptors could be computed". The engine is version-tolerant (tries key aliases, isolates each descriptor group so one missing key doesn't nuke the rest) and the error now names the missing keys + the tblite version + the fix, but the clean solution is to match the pinned 0.6.0.
  • Both predict.py and webui.py run a smoke test first (a water molecule through the full pipeline): if tblite is absent, the wrong version, or produces fewer descriptors than expected, you get a clear error naming the problem (missing descriptors, tblite version, pin tblite-python=0.6.0) instead of a silent/meaningless S = 0.
  • python: numpy, joblib, xgboost (xgboost is needed to load the model; no shap needed — the per-reaction decomposition uses xgboost's native pred_contribs).
  • web UI only (all optional, feature degrades gracefully if absent): flask, ase+matplotlib (structure pictures), rdkit (SMILES input — type CC=O, click build 3D; RDKit embeds the molecule and it's GFN2-optimized to match the descriptors, with charge & spin auto-detected).

Usage

conda activate xtb
python predict.py  MOL1:charge:uhf:coeff  MOL2:charge:uhf:coeff ...
# reactants get a negative coeff, products positive. uhf = # unpaired electrons (2S).
# example — F2 -> 2 F•  :
python predict.py  f2.xyz:0:0:-1   f.xyz:0:1:2

Output: predicted S + its graded triage tier, then a WHY decomposition — S = base + Σ pushes split into the three physical levers (bond rearrangement / near degeneracy / charge anisotropy) + size, and the top individual descriptor pushes (exact TreeSHAP).

Web UI (optional)

A Flask browser front-end (styled after claude.ai — warm palette, coral accent) over the same predictor (pip install flask):

conda activate xtb
python webui.py [port]        # default 8000
# open http://127.0.0.1:8000  (binds to localhost; if on a remote box, SSH-tunnel:
#   ssh -L 8000:127.0.0.1:8000 user@host )

The page is branded "Chemical Triage for Density Sensitivity" (TCCL / Yonsei logo + contacts). Paste each species' XYZ + set its charge / unpaired-electrons / coefficient, or pick from the left-hand Examples tree — a foldable hierarchy of every GMTKN55 reaction (cluster → subset → reaction number; hover a number to preview the reaction's molecules, click to load) and every periodic system from ../periodic (periodic → X23 / ICE7 / CO·NaCl / graphite / covalent → system). Species load on demand (/reaction; cluster geometries from the bundled refdb.pkl.gz, periodic from pbc_examples.json), then click Predict. It renders the reaction as ball-and-stick molecules (same ASE drawing as the waterfalls figures), the predicted S_predicted, the colour-coded tier, the GMTKN55-calibrated risk, the lever decomposition (bars), and the top-10 drivers table. A ⤓ Save PNG button (top-right of the result card) captures the result exactly as styled — reaction · gauge · tier · WHY decomposition, minus the raw tblite output — as a single image, rendered client-side from the live DOM by the bundled html2canvas (no server round-trip, no external network). Each prediction is logged under logs/web_<timestamp>/. Needs flask + numpy/joblib/xgboost + tblite; the geometry opt for SMILES inputs and the molecule pictures use ase (if ase is absent the UI keeps the MMFF geometry and omits pictures).

What gets logged (so you can look inside)

Every run creates logs/run_<timestamp>/:

  • run.log — the full console output (prediction + WHY decomposition).
  • features.txt — the aggregated model feature vector, column = value (* imputed = was NaN→0).
  • descriptors.json — every per-compound tblite descriptor for each species, plus the decomposition.
  • species_<i>_<name>/ — the tblite GFN2 results of the 2 single-points per species (1_property_300K.out, fod_5000K.out) plus the raw arrays (*__arrays.json: orbital energies/occupations, charges, dipole, quadrupole) — inspect these to see exactly what tblite produced.

Plus a running one-line-per-run index at logs/history.log (timestamp S flag command).

Logs accumulate; delete old logs/run_* folders whenever you like — they're never read back.

Periodic (PBC) mode — research extrapolation

Each species can be a cluster (molecule) or a periodic cell. In the web UI, switch a species block to periodic (a warning is shown first); a lattice field (3 cell vectors in Å) and a supercell field (n1 n2 n3, Γ k-fold) appear. In the CLI, put an extended-XYZ Lattice="ax ay az bx by bz cx cy cz" comment on line 2 of the crystal .xyz (auto-detected). Descriptors are computed per provided cell via tblite GFN2 with periodic boundary conditions (Makov–Payne alignment for charged cells; supercell → band-folded k-sampling). Set coefficients so atoms balance across the reaction's cells & molecules.

Intensive S/Vmon×100 for molecular crystals. When the reaction has exactly one periodic species (the crystal; in the CLI: crystal coeff −1 + its Z molecules with Σ positive coeffs = Z), the run also reports S/Vmon×100 with Vmon = unit-cell volume / Z — the per-volume metric the DC-DFT paper's Fig. 4 actually uses for molecular-crystal cohesion (computed from the cell alone, no external tables; reproduces the paper's S/Vmon*100 columns exactly, see ../periodic/real_raw_data_x23/_ice7). Compare crystals on that scale, not the molecular S bands: dense small-molecule crystals (ice, Vmon ≈ 20–31 ų) sit structurally higher in S/Vmon than large-molecule crystals at equal per-molecule S.

This is out of domain. The model is trained only on GMTKN55 molecules — for solids trust the flag, not the S value; adsorption (CO/NaCl) is a known failure; metals/semimetals need a supercell; periodic dipole/quad are ill-defined. The periodic path is validated separately in ../periodic/. NOTE the v22 model predicts S^SC, not the r2SCAN S the paper reference uses, so this is a cross-target check: per-molecule S tracks the paper's r2SCAN sensitivities at Spearman +0.78 ex-CO/NaCl (v19, trained on r2SCAN, reached +0.88); in the paper's own S/Vmon metric the molecular crystals pool to +0.79 (X23 +0.54, ice +0.36). Magnitudes underpredicted (~2x low on X23). Trust rank direction, not the value.

Provenance

Byte-identical model/pipeline to ../ (the full v22 folder; the model is v22 = v21 retrained on S^SC); this is just the prediction-time subset with logging added. See ../README.md and ../train/FEATURES.md for the model, features, and validation (10×5-fold cross-validation r² = 0.729, triage AUC = 0.958; 104 features).