Skip to content

Commit 2605783

Browse files
committed
Document assignment artifact generator
1 parent 9ab00cd commit 2605783

4 files changed

Lines changed: 84 additions & 0 deletions

File tree

vimms/AssignmentChemicalArtifact.py

Lines changed: 25 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -1,3 +1,5 @@
1+
"""Materialize assignment-noise peak tables as VIMMS chemicals and mzML sidecars."""
2+
13
from __future__ import annotations
24

35
import json
@@ -54,6 +56,8 @@ class AssignmentChemicalArtifactConfig:
5456
def _chromatogram_for_peak(
5557
apex_seconds: float, width_seconds: float
5658
) -> tuple[EmpiricalChromatogram, float]:
59+
"""Build a simple empirical chromatogram centered on one picked feature."""
60+
5761
half_width = max(float(width_seconds), 1.0)
5862
start = max(0.0, float(apex_seconds) - half_width)
5963
apex = max(float(apex_seconds), start + 0.25)
@@ -69,6 +73,8 @@ def _make_ms2_children(
6973
row: pd.Series,
7074
rng: np.random.Generator,
7175
) -> list[MSN]:
76+
"""Create fragment children for chemicals with available/plausible MS2 support."""
77+
7278
score = float(row.get("best_ms2_score", 0.0))
7379
available = float(row.get("ms2_available", 0.0)) > 0.0
7480
if not available and score < 0.25:
@@ -97,6 +103,8 @@ def _materialize_assignment_chemicals(
97103
config: AssignmentChemicalArtifactConfig,
98104
rng: np.random.Generator,
99105
) -> tuple[list[UnknownChemical], pd.DataFrame]:
106+
"""Convert every truth/artifact peak row into a VIMMS UnknownChemical."""
107+
100108
chemicals: list[UnknownChemical] = []
101109
truth_rows: list[dict[str, Any]] = []
102110

@@ -109,6 +117,9 @@ def _materialize_assignment_chemicals(
109117
float(row.get("peak_width", 0.08)) * 60.0 * 1.5,
110118
)
111119
chromatogram, start_rt = _chromatogram_for_peak(target_rt_seconds, width_seconds)
120+
# UnknownChemical stores a neutral mass, but the picked-table contract is
121+
# already m/z based. Subtracting a proton gives a chemical that emits at
122+
# the target m/z without changing model-facing values.
112123
neutral_for_unknown = max(1.0, target_mz - PROTON_MASS)
113124
max_intensity = max(
114125
float(row.get("intensity", np.exp(row.get("log_intensity", 12.0)))), 1.0
@@ -149,6 +160,8 @@ def _materialize_assignment_chemicals(
149160

150161

151162
def _output_path(config: AssignmentChemicalArtifactConfig) -> Path | None:
163+
"""Return an output directory only when mzML writing is requested."""
164+
152165
if not config.write_mzml:
153166
return None
154167
if config.output_dir is None:
@@ -162,6 +175,8 @@ def _run_scan_simulation(
162175
chemicals: list[UnknownChemical],
163176
config: AssignmentChemicalArtifactConfig,
164177
) -> tuple[Environment, Path | None]:
178+
"""Run VIMMS TopN scan simulation for optional audit/export mzML output."""
179+
165180
out_dir = _output_path(config)
166181
out_file = f"{config.prefix}.mzML" if out_dir is not None else None
167182
mass_spec = IndependentMassSpectrometer(
@@ -195,6 +210,8 @@ def _run_scan_simulation(
195210

196211

197212
def _scan_summary(scans: Mapping[int, list[Any]]) -> pd.DataFrame:
213+
"""Flatten VIMMS controller scans into a lightweight audit table."""
214+
198215
rows: list[dict[str, Any]] = []
199216
for ms_level, level_scans in sorted(scans.items()):
200217
for scan in level_scans:
@@ -211,12 +228,16 @@ def _scan_summary(scans: Mapping[int, list[Any]]) -> pd.DataFrame:
211228

212229

213230
def _empty_scan_summary() -> pd.DataFrame:
231+
"""Return the scan-summary schema used when mzML writing is disabled."""
232+
214233
return pd.DataFrame(
215234
columns=["scan_id", "ms_level", "rt_seconds", "n_peaks", "tic"]
216235
)
217236

218237

219238
def _truth_table_from_peaks(peak_table: pd.DataFrame) -> pd.DataFrame:
239+
"""Extract peak-level truth columns from the model-facing peak table."""
240+
220241
columns = [
221242
"peak_id",
222243
"true_label",
@@ -238,6 +259,8 @@ def _truth_table_from_peaks(peak_table: pd.DataFrame) -> pd.DataFrame:
238259

239260

240261
def _jsonable(value: Any) -> Any:
262+
"""Recursively convert dataclass metadata values to JSON-safe objects."""
263+
241264
if isinstance(value, Path):
242265
return str(value)
243266
if isinstance(value, tuple):
@@ -276,6 +299,8 @@ def generate_assignment_chemical_artifact(
276299
rng=rng,
277300
)
278301
if config.write_mzml:
302+
# Scans are sidecars only. The returned peak table stays bit-identical
303+
# to the direct picked-table generator used by IonClassifier.
279304
env, mzml_path = _run_scan_simulation(chemicals, config)
280305
scan_summary = _scan_summary(env.controller.scans)
281306
else:

vimms/AssignmentNoise.py

Lines changed: 18 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -1,3 +1,5 @@
1+
"""Generate labelled picked-peak assignment scenarios for IonClassifier."""
2+
13
from __future__ import annotations
24

35
import json
@@ -161,6 +163,8 @@ def role_properties(role: str, carbon_count: float) -> dict[str, float]:
161163

162164

163165
def _validate_config(config: AssignmentScenarioConfig) -> None:
166+
"""Validate dimensions and explicit candidate-presence overrides."""
167+
164168
profile = config.profile
165169
if profile.n_candidates < 1:
166170
raise ValueError("n_candidates must be at least 1")
@@ -180,12 +184,16 @@ def _presence_pattern(
180184
profile: AssignmentNoiseProfile,
181185
present_pattern: tuple[int, ...] | None,
182186
) -> np.ndarray:
187+
"""Sample which candidates are truly present in a local assignment problem."""
188+
183189
if present_pattern is not None:
184190
return np.asarray(present_pattern, dtype=np.int8)
185191
present = (rng.random(profile.n_candidates) < profile.p_present).astype(np.int8)
186192
if not present.any():
187193
present[0] = 1
188194
if profile.n_candidates > 1 and present.all() and profile.p_matched_decoy > 0.0:
195+
# Keep one plausible absent candidate available so matched-decoy
196+
# clusters can test false positive rejection.
189197
present[-1] = 0
190198
return present
191199

@@ -195,6 +203,8 @@ def _make_candidate_table(
195203
profile: AssignmentNoiseProfile,
196204
present: np.ndarray,
197205
) -> pd.DataFrame:
206+
"""Create candidate metadata around one target-like local mass/RT region."""
207+
198208
base_mass = float(rng.uniform(profile.mass_min, profile.mass_max))
199209
base_rt = float(rng.uniform(profile.rt_min, profile.rt_max))
200210
base_c = int(rng.integers(5, 45))
@@ -236,6 +246,8 @@ def _make_candidate_table(
236246

237247

238248
def _make_ion_table(candidates: pd.DataFrame, roles: tuple[str, ...]) -> pd.DataFrame:
249+
"""Expand candidates into theoretical role/adduct/isotope ions."""
250+
239251
rows: list[dict[str, Any]] = []
240252
for _, candidate in candidates.iterrows():
241253
candidate_index = int(candidate["candidate_index"])
@@ -264,6 +276,8 @@ def _make_ion_truth_table(
264276
candidate_table: pd.DataFrame,
265277
ion_table: pd.DataFrame,
266278
) -> pd.DataFrame:
279+
"""Record observed and missing true companion ions for every candidate role."""
280+
267281
observed = (
268282
peak_table.loc[peak_table["true_label"] > 0, ["true_label", "peak_id", "source_type"]]
269283
.drop_duplicates("true_label")
@@ -321,13 +335,17 @@ def generate_assignment_peak_table(
321335
latent_profiles: dict[int, np.ndarray] = {}
322336
for _, candidate in candidate_table.iterrows():
323337
candidate_index = int(candidate["candidate_index"])
338+
# Candidate ions share a latent abundance profile so isotope/adduct
339+
# companions look correlated, while artifacts can partially mimic it.
324340
latent_profiles[candidate_index] = rng.normal(
325341
float(candidate["expected_log_intensity"]),
326342
rng.uniform(0.35, 0.95),
327343
size=profile.n_bio_samples,
328344
)
329345

330346
peaks: list[dict[str, Any]] = []
347+
# Add true evidence first, then layer increasingly adversarial background
348+
# profiles around the same theoretical ions.
331349
_add_candidate_ions(peaks, rng, candidate_table, ion_table, profile, roles, latent_profiles)
332350
_add_structured_interferents(peaks, rng, candidate_table, ion_table, profile)
333351
_add_coeluting_isobars(peaks, rng, candidate_table, ion_table, profile)

vimms/AssignmentNoisePeaks.py

Lines changed: 35 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -1,3 +1,5 @@
1+
"""Peak-row generation helpers for assignment-noise scenarios."""
2+
13
from __future__ import annotations
24

35
from typing import Any, Mapping
@@ -9,6 +11,8 @@
911

1012

1113
def _skew(values: np.ndarray) -> float:
14+
"""Return sample skewness with a stable value for near-constant profiles."""
15+
1216
arr = np.asarray(values, dtype=float)
1317
sd = arr.std()
1418
if sd < 1e-8:
@@ -17,6 +21,8 @@ def _skew(values: np.ndarray) -> float:
1721

1822

1923
def _beta_score(rng: np.random.Generator, true: bool, fp_prob: float) -> float:
24+
"""Sample an MS2-like support score for true and false candidate evidence."""
25+
2026
if true:
2127
return float(rng.beta(9.0, 2.2))
2228
if rng.random() < fp_prob:
@@ -49,6 +55,8 @@ def _peak_row(
4955
ms2_true_score: bool = False,
5056
metadata: Mapping[str, Any] | None = None,
5157
) -> dict[str, Any]:
58+
"""Create one normalized picked-peak row with feature and truth metadata."""
59+
5260
row: dict[str, Any] = {
5361
"mz": float(mz),
5462
"rt": float(rt),
@@ -93,6 +101,8 @@ def _peak_row(
93101
available_prob = 0.70 if ms2_true_score else min(0.85, 0.25 + 0.55 * fp_prob)
94102
available = bool(rng.random() < available_prob)
95103
if not available:
104+
# Keep weak residual scores so the feature exists, but make it clear
105+
# that no useful spectrum should support this row.
96106
for candidate_index in range(n_candidates):
97107
row[f"ms2_score_cand_{candidate_index}"] *= float(rng.uniform(0.0, 0.25))
98108
ms2_scores = [
@@ -107,6 +117,8 @@ def _peak_row(
107117

108118

109119
def _finalize_peak_table(peaks: list[dict[str, Any]], rng: np.random.Generator) -> pd.DataFrame:
120+
"""Compute local context features and assign stable peak identifiers."""
121+
110122
df = pd.DataFrame(peaks)
111123
if df.empty:
112124
return df
@@ -151,6 +163,8 @@ def _add_diffuse_background(
151163
profile: AssignmentNoiseProfile,
152164
n: int,
153165
) -> None:
166+
"""Append unstructured clutter and blank-contaminant background peaks."""
167+
154168
mz_min = float(candidates["neutral_mass"].min() - 25.0)
155169
mz_max = float(candidates["neutral_mass"].max() + 45.0)
156170
rt_center = float(candidates["pred_rt"].mean())
@@ -193,6 +207,8 @@ def _add_candidate_ions(
193207
roles: tuple[str, ...],
194208
latent_profiles: dict[int, np.ndarray],
195209
) -> None:
210+
"""Append true candidate-role peaks plus split/merged artifacts."""
211+
196212
for _, candidate in candidates.iterrows():
197213
candidate_index = int(candidate["candidate_index"])
198214
if int(candidate["present"]) <= 0:
@@ -272,6 +288,8 @@ def _add_candidate_ions(
272288
emitted_any = True
273289

274290
if rng.random() < profile.p_split:
291+
# Split/shoulder artifacts resemble the parent ion in m/z/RT and
292+
# abundance profile but remain background labels.
275293
split_profile = peak_profile + rng.normal(0.0, 0.15, size=profile.n_bio_samples)
276294
peaks.append(
277295
_peak_row(
@@ -300,6 +318,9 @@ def _add_candidate_ions(
300318
)
301319

302320
if rng.random() < profile.p_merge:
321+
# Merged peaks intentionally have good-looking evidence but are
322+
# labelled background because the picked feature is not the clean
323+
# theoretical ion.
303324
merged_profile = peak_profile + rng.normal(0.0, 0.18, size=profile.n_bio_samples)
304325
peaks.append(
305326
_peak_row(
@@ -335,6 +356,8 @@ def _add_structured_interferents(
335356
ions: pd.DataFrame,
336357
profile: AssignmentNoiseProfile,
337358
) -> None:
359+
"""Append wrong-source peaks near theoretical ions and candidate RTs."""
360+
338361
for _, ion in ions.iterrows():
339362
if rng.random() >= profile.p_interferent:
340363
continue
@@ -408,6 +431,8 @@ def _add_coeluting_isobars(
408431
ions: pd.DataFrame,
409432
profile: AssignmentNoiseProfile,
410433
) -> None:
434+
"""Append high-quality background peaks that coelute with theoretical ions."""
435+
411436
for _, ion in ions.iterrows():
412437
if rng.random() >= profile.p_coeluting_isobar:
413438
continue
@@ -460,6 +485,8 @@ def _add_matched_decoy_clusters(
460485
profile: AssignmentNoiseProfile,
461486
latent_profiles: dict[int, np.ndarray],
462487
) -> None:
488+
"""Append plausible isotope/adduct clusters for absent or wrong candidates."""
489+
463490
candidate_indices = list(candidates["candidate_index"].astype(int))
464491
absent = list(
465492
candidates.loc[candidates["present"].astype(int) == 0, "candidate_index"].astype(int)
@@ -479,6 +506,8 @@ def _add_matched_decoy_clusters(
479506
cand_ions = ions[ions["candidate_index"] == candidate_index].sort_values("role_index")
480507
decoy_log_a = float(candidate["expected_log_intensity"] + rng.normal(-0.05, 0.55))
481508
if latent_profiles and rng.random() < 0.65:
509+
# Reusing a true latent profile makes a decoy cluster look like a
510+
# coherent biological signal instead of independent random clutter.
482511
base_profile = latent_profiles[int(rng.choice(list(latent_profiles.keys())))].copy()
483512
base_profile = base_profile - base_profile.mean() + decoy_log_a
484513
else:
@@ -541,6 +570,8 @@ def _add_matched_decoy_clusters(
541570
},
542571
)
543572
if rng.random() < profile.p_decoy_false_ms2:
573+
# Some decoy clusters get candidate-specific false MS2 support,
574+
# matching the failure mode the classifier must reject.
544575
row[f"ms2_score_cand_{candidate_index}"] = float(
545576
max(row.get(f"ms2_score_cand_{candidate_index}", 0.0), rng.beta(5.0, 2.0))
546577
)
@@ -590,6 +621,8 @@ def _trim_peaks(
590621
rng: np.random.Generator,
591622
max_peaks: int,
592623
) -> list[dict[str, Any]]:
624+
"""Limit peak count while preserving true ions and hard negative artifacts."""
625+
593626
if len(peaks) <= max_peaks:
594627
return peaks
595628
true_rows = [p for p in peaks if int(p["true_label"]) > 0]
@@ -605,6 +638,8 @@ def _trim_peaks(
605638
easy_rows = [
606639
p for p in peaks if int(p["true_label"]) == 0 and p["source_type"] not in hard_types
607640
]
641+
# True ions are always retained. Hard negatives are retained before easy
642+
# diffuse clutter so stress profiles remain adversarial after trimming.
608643
remaining = max(0, max_peaks - len(true_rows))
609644
if len(hard_rows) > remaining:
610645
idx = rng.choice(len(hard_rows), size=remaining, replace=False)

vimms/AssignmentScanSimulation.py

Lines changed: 6 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -1,3 +1,5 @@
1+
"""Backward-compatible wrappers for the old assignment scan-simulation API."""
2+
13
from __future__ import annotations
24

35
import warnings
@@ -13,6 +15,8 @@
1315

1416

1517
def generate_assignment_scan_artifact(*args, **kwargs):
18+
"""Deprecated alias for ``generate_assignment_chemical_artifact``."""
19+
1620
warnings.warn(
1721
"generate_assignment_scan_artifact is deprecated; use "
1822
"generate_assignment_chemical_artifact instead.",
@@ -23,6 +27,8 @@ def generate_assignment_scan_artifact(*args, **kwargs):
2327

2428

2529
def write_assignment_scan_artifact(*args, **kwargs):
30+
"""Deprecated alias for ``write_assignment_chemical_artifact``."""
31+
2632
warnings.warn(
2733
"write_assignment_scan_artifact is deprecated; use "
2834
"write_assignment_chemical_artifact instead.",

0 commit comments

Comments
 (0)