docs(aroma): record the aroma evaluation + defer decision

Honest negative result (keller_2016 CV-R² ≤ 0 across all 20 descriptors); defer aroma, keep OpenPOM scaffold for licensed/customer data, lead with taste. + keller_2016 attribution.

Author: rvnminers-A-and-N

docs/AROMA.md

@@ -0,0 +1,66 @@
+# Aroma — evaluated and deferred
+
+Flavor = taste **+** aroma. The taste side ships and is clean (sweet/bitter/umami
+~0.95 AUROC + intensity + sour/salty). The aroma side is **deliberately deferred** —
+not for lack of effort, but because the public, commercially-clean data is
+insufficient. This documents the evaluation, so the decision is legible.
+
+## What an aroma model needs
+
+A model mapping a molecule's **structure → odor descriptors** (sweet, floral,
+woody…). The state of the art is the **Principal Odor Map** (Lee et al., *Science*
+2023; Google/Osmo), reimplemented open-source as **OpenPOM** (`BioMachineLearning/
+openpom`, MIT) — a message-passing GNN trained on ~5,000 **expert-labeled** molecules
+(the GS-LF dataset).
+
+## The licensing wall (audited)
+
+OpenPOM's *code* is MIT, but its *training data* is restricted. A full audit of the
+Pyrfume catalogue (see [`DATA-SOURCES.md`](DATA-SOURCES.md)) found that **every rich
+odor-descriptor dataset is proprietary or NonCommercial** — Leffingwell, GoodScents,
+Arctander, Flavornet (©Datu), FlavorDB / FooDB (NC), OlfactionBase ("all rights
+reserved"), AromaDB (CSIR), Dravnieks (ASTM ©), sharma_2021 (ACS ©), snitz_2019
+(CC-BY-NC). The **only** commercially-clean odor-descriptor set is **`keller_2016`**
+(Keller & Vosshall 2016, *BMC Neuroscience*, CC-BY-4.0; ~480 molecules, 20 descriptors).
+
+## The empirical result
+
+We aggregated `keller_2016` ([`training/build_aroma_dataset.py`](../training/build_aroma_dataset.py))
+and trained one RandomForest regressor per descriptor on 2048-bit Morgan fingerprints
+([`training/train_aroma.py`](../training/train_aroma.py)), scored by **honest 5-fold
+cross-validation** (400 trees):
+
+| descriptor | CV-R² | descriptor | CV-R² | descriptor | CV-R² | descriptor | CV-R² |
+|---|---|---|---|---|---|---|---|
+| acid | −0.24 | cold | −0.19 | fruit | −0.05 | sour | −0.05 |
+| ammonia | −0.15 | decayed | −0.20 | garlic | −0.12 | spices | −0.12 |
+| bakery | −0.22 | edible | −0.19 | grass | −0.45 | sweaty | +0.03 |
+| burnt | −0.13 | fish | −0.04 | musky | −0.23 | sweet | −0.05 |
+| chemical | −0.10 | flower | −0.12 | warm | −0.28 | wood | −0.08 |
+
+**All 20 descriptors scored CV-R² ≤ 0** (range −0.45 to +0.03) — every model is
+*worse than predicting the mean*. **0/20 usable heads.**
+
+## Root cause
+
+`keller_2016` is **naive-subject** data: random volunteers rating *unfamiliar*
+molecules on a 0–100 scale. People can't reliably name what a molecule smells like,
+so the labels are noise (the per-descriptor means compress to ~22–30 for nearly every
+molecule). The *learnable* odor data uses **expert** labels — which is exactly the
+data that's restricted. **That is the structural reason the entire field trains on
+GS-LF**, and why a clean public aroma model isn't currently possible.
+
+## Decision
+
+We do **not** ship a negative-R² model — it would output confident, wrong smells, and
+a flavor chemist would catch it instantly (worse than nothing). `predict_aroma()`
+returns an honest "not available" marker, and **we lead with the taste engine**.
+
+The **OpenPOM scaffold (`training/train_odor.py`) is kept** as the aroma engine for
+when clean fuel exists:
+1. **License PMP 2001** (~$2,775, Leffingwell & Associates) → re-curate → train.
+2. **A customer's own odor data** (the paid pilot) → train on-prem.
+3. A future open *expert-labeled* dataset, if one emerges.
+
+Until then: **aroma comes with your data** — which is on-thesis (public data proves
+the method; the customer's data unlocks the rest).

docs/README.md

@@ -8,5 +8,6 @@ Architecture, capability catalogue, and design decisions.
 - `SOURCES.md` — data sources, libraries, research, and license attribution
 - `API-CONTRACT.md` — the fixed JSON contract the .NET API exposes and the React UI consumes
 - `DATA-SOURCES.md` — data acquisition tracker: what each source unlocks and how to get it
+- `AROMA.md` — why the aroma model is deferred (clean-data audit + the empirical evaluation)
 
 These document *what we're building and why*, so decisions don't get re-argued.

docs/SOURCES.md

@@ -42,6 +42,7 @@ record of provenance, not legal advice. Get an IP/OSS-license review before ship
 
 **Aroma**
 - **Pyrfume** + **Leffingwell / GoodScents (GS-LF)** odor datasets — the usual training data behind OpenPOM, but **RESTRICTED and NOT USED**: Leffingwell's manifest cites use restrictions (*John Leffingwell & Google*); GoodScents/Arctander/Flavornet (© Datu Inc.) are likewise proprietary. We **exclude all of them** (the demo may go to a customer / commercial use). The aroma model will use only commercial-clean **open** odor data (CC-BY sets like `keller_2016`; smaller — see `DATA-SOURCES.md`). The OpenPOM *code* is MIT.
+- **keller_2016** — Keller & Vosshall (2016), *BMC Neuroscience*, **CC-BY-4.0**; ~480 molecules with naive-subject odor-descriptor ratings. The only commercially-clean odor-descriptor set — evaluated for the aroma model and found too noisy to learn from (CV-R² ≤ 0 across all 20 descriptors; see `docs/AROMA.md`).
 
 **Safety / regulatory (lookups — data-gated)**
 - **FEMA GRAS list** — usual/maximum use levels for the dosing analyzer. (FEMA.)

training/build_aroma_dataset.py

@@ -0,0 +1,60 @@
+"""
+build_aroma_dataset.py — assemble the clean aroma (odor-descriptor) training table.
+
+Source: keller_2016 (Keller & Vosshall 2016, BMC Neuroscience, CC-BY-4.0) — the only
+commercially-clean odor-descriptor dataset available (~480 molecules, 20 descriptors).
+Each molecule gets a mean 0-100 panel rating per descriptor.
+
+Output: aroma_master.parquet  (smiles + 20 descriptor columns, 0-100)
+"""
+from pathlib import Path
+
+import pandas as pd
+from rdkit import Chem
+
+KELLER = Path("aroma/keller_2016")
+DESCRIPTORS = ["ACID", "AMMONIA/URINOUS", "BAKERY", "BURNT", "CHEMICAL", "COLD",
+               "DECAYED", "EDIBLE", "FISH", "FLOWER", "FRUIT", "GARLIC", "GRASS",
+               "MUSKY", "SOUR", "SPICES", "SWEATY", "SWEET", "WARM", "WOOD"]
+
+
+def canon(smiles):
+    if not isinstance(smiles, str):
+        return None
+    m = Chem.MolFromSmiles(smiles)
+    return Chem.MolToSmiles(m) if m else None
+
+
+def build():
+    mol = pd.read_csv(KELLER / "molecules.csv")
+    sti = pd.read_csv(KELLER / "stimuli.csv")
+    beh = pd.read_csv(KELLER / "behavior.csv")
+
+    beh = beh[beh["MeasurementValue"].isin(DESCRIPTORS)].copy()
+    beh["Value"] = pd.to_numeric(beh["Value"], errors="coerce")
+    beh = beh.dropna(subset=["Value"])
+
+    # mean rating per (Stimulus, descriptor) across the panel
+    sm = beh.groupby(["Stimulus", "MeasurementValue"])["Value"].mean().unstack()
+
+    # Stimulus -> CID (single-molecule) -> canonical SMILES
+    sti = sti[["Stimulus", "CIDs"]].copy()
+    sti["CID"] = pd.to_numeric(sti["CIDs"], errors="coerce")
+    sm = sm.join(sti.set_index("Stimulus")[["CID"]]).dropna(subset=["CID"])
+    sm["CID"] = sm["CID"].astype(int)
+    mol_map = mol.dropna(subset=["CID"]).drop_duplicates("CID").set_index("CID")["CanonicalSMILES"]
+    sm["smiles"] = sm["CID"].map(mol_map).map(canon)
+    sm = sm.dropna(subset=["smiles"])
+
+    # molecule level: mean across the molecule's stimuli (concentrations)
+    agg = sm.groupby("smiles")[DESCRIPTORS].mean().reset_index()
+    agg.to_parquet("aroma_master.parquet")
+
+    print(f"aroma_master.parquet: {len(agg)} molecules x {len(DESCRIPTORS)} descriptors")
+    for d in DESCRIPTORS:
+        print(f"  {d:16s} mean={agg[d].mean():5.1f}  max={agg[d].max():5.1f}  "
+              f">=15: {int((agg[d] >= 15).sum())}")
+
+
+if __name__ == "__main__":
+    build()

training/train_aroma.py

@@ -0,0 +1,66 @@
+"""
+train_aroma.py — odor-descriptor regressors (one per descriptor) on Morgan fingerprints.
+
+Reads aroma_master.parquet (keller_2016, CC-BY); trains a RandomForestRegressor per
+descriptor predicting its 0-100 panel rating; reports HONEST 5-fold CV R2 (small,
+noisy panel-mean data — expect modest, descriptor-dependent scores). Saves the ones
+that clear a minimum CV-R2 to aroma_models/.
+"""
+from pathlib import Path
+
+import joblib
+import numpy as np
+import pandas as pd
+from rdkit import Chem
+from rdkit.Chem import DataStructs, rdFingerprintGenerator
+from sklearn.ensemble import RandomForestRegressor
+from sklearn.model_selection import cross_val_score
+
+FP_BITS, FP_RADIUS = 2048, 2
+_MORGAN = rdFingerprintGenerator.GetMorganGenerator(radius=FP_RADIUS, fpSize=FP_BITS)
+MIN_R2 = 0.10  # below this a descriptor is too noisy to ship an honest head
+OUT = Path("aroma_models")
+OUT.mkdir(exist_ok=True)
+for s in OUT.glob("*_reg.joblib"):
+    s.unlink()
+
+DESCRIPTORS = ["ACID", "AMMONIA/URINOUS", "BAKERY", "BURNT", "CHEMICAL", "COLD",
+               "DECAYED", "EDIBLE", "FISH", "FLOWER", "FRUIT", "GARLIC", "GRASS",
+               "MUSKY", "SOUR", "SPICES", "SWEATY", "SWEET", "WARM", "WOOD"]
+
+
+def fp(smiles):
+    m = Chem.MolFromSmiles(smiles)
+    if m is None:
+        return None
+    bv = _MORGAN.GetFingerprint(m)
+    arr = np.zeros((FP_BITS,), dtype=np.int8)
+    DataStructs.ConvertToNumpyArray(bv, arr)
+    return arr
+
+
+df = pd.read_parquet("aroma_master.parquet")
+feats, keep = [], []
+for i, s in enumerate(df["smiles"]):
+    f = fp(s)
+    if f is not None:
+        feats.append(f)
+        keep.append(i)
+X = np.array(feats)
+df = df.iloc[keep].reset_index(drop=True)
+
+print(f"training odor-descriptor regressors on {len(df)} molecules:")
+kept = 0
+for d in DESCRIPTORS:
+    y = df[d].values
+    mask = ~np.isnan(y)
+    Xd, yd = X[mask], y[mask]
+    r2 = cross_val_score(RandomForestRegressor(n_estimators=400, n_jobs=-1, random_state=42),
+                         Xd, yd, cv=5, scoring="r2").mean()
+    flag = "kept" if r2 >= MIN_R2 else "drop (too noisy)"
+    if r2 >= MIN_R2:
+        reg = RandomForestRegressor(n_estimators=400, n_jobs=-1, random_state=42).fit(Xd, yd)
+        joblib.dump(reg, OUT / f"{d.replace('/', '_')}_reg.joblib")
+        kept += 1
+    print(f"  {d:16s} n={int(mask.sum()):3d}  CV-R2={r2:+.2f}  -> {flag}")
+print(f"\nkept {kept}/{len(DESCRIPTORS)} descriptor heads (CV-R2 >= {MIN_R2}) -> aroma_models/")