|
| 1 | +""" |
| 2 | +Differential parity test: Python `CalibrationState` vs. the real MATLAB |
| 3 | +staircase algorithm, executed through Octave (or MATLAB) as a reference engine. |
| 4 | +
|
| 5 | +Unlike `test_calibration.py` (which asserts against hand-derived numbers), this |
| 6 | +feeds the *same* scripted win/loss sequences to both implementations and checks |
| 7 | +the per-trial target windows match. The reference lives in |
| 8 | +`tests/matlab_ref/calibration_ref.m`, extracted verbatim from |
| 9 | +`MID-updated/scripts/PresentTarget.m` lines 14-48. |
| 10 | +
|
| 11 | +If neither `octave` nor `matlab` is on PATH, the parity tests skip — but a loud |
| 12 | +warning is emitted so a green run never silently implies parity was verified. |
| 13 | +""" |
| 14 | +from __future__ import annotations |
| 15 | + |
| 16 | +import shutil |
| 17 | +import subprocess |
| 18 | +import warnings |
| 19 | +from pathlib import Path |
| 20 | + |
| 21 | +import pytest |
| 22 | + |
| 23 | +from mid_det import config |
| 24 | +from mid_det.calibration import CalibrationState |
| 25 | + |
| 26 | +_REF_DIR = Path(__file__).parent / "matlab_ref" |
| 27 | + |
| 28 | +# cue id (1..6, as used by var.cues in MATLAB) -> (polarity, magnitude) |
| 29 | +_CUE_TO_PM = {cue_id: pm for pm, cue_id in config.TRIAL_TYPE_MAP.items()} |
| 30 | + |
| 31 | + |
| 32 | +# ── reference-engine discovery ──────────────────────────────────────────────── |
| 33 | + |
| 34 | +def _find_engine() -> tuple[str, str] | None: |
| 35 | + """Return (kind, executable_path) for octave or matlab, else None.""" |
| 36 | + for kind in ("octave", "matlab"): |
| 37 | + path = shutil.which(kind) |
| 38 | + if path: |
| 39 | + return kind, path |
| 40 | + return None |
| 41 | + |
| 42 | + |
| 43 | +_ENGINE = _find_engine() |
| 44 | + |
| 45 | + |
| 46 | +@pytest.fixture(scope="session", autouse=True) |
| 47 | +def _warn_if_no_engine(): |
| 48 | + """Make a missing reference engine impossible to overlook in the summary.""" |
| 49 | + if _ENGINE is None: |
| 50 | + warnings.warn( |
| 51 | + "calibration MATLAB-parity NOT verified — no Octave/MATLAB engine " |
| 52 | + "found on PATH; install octave (`brew install octave`) to run it.", |
| 53 | + stacklevel=2, |
| 54 | + ) |
| 55 | + yield |
| 56 | + |
| 57 | + |
| 58 | +requires_engine = pytest.mark.skipif( |
| 59 | + _ENGINE is None, reason="no octave/matlab engine on PATH (see emitted warning)" |
| 60 | +) |
| 61 | + |
| 62 | + |
| 63 | +# ── drivers ─────────────────────────────────────────────────────────────────── |
| 64 | + |
| 65 | +def _python_windows( |
| 66 | + cues: list[int], wins: list[int], base_rt: float, rt_change: float |
| 67 | +) -> list[float]: |
| 68 | + """Per-trial target window from the Python implementation under test.""" |
| 69 | + cal = CalibrationState(base_rt_s=base_rt, rt_change_s=rt_change) |
| 70 | + out: list[float] = [] |
| 71 | + for cue_id, win in zip(cues, wins): |
| 72 | + polarity, magnitude = _CUE_TO_PM[cue_id] |
| 73 | + out.append(cal.next_target_dur_s(polarity, magnitude)) |
| 74 | + cal.record_outcome(polarity, magnitude, bool(win)) |
| 75 | + return out |
| 76 | + |
| 77 | + |
| 78 | +def _octave_windows( |
| 79 | + cues: list[int], wins: list[int], base_rt: float, rt_change: float |
| 80 | +) -> list[float]: |
| 81 | + """Per-trial target window from the MATLAB/Octave reference engine.""" |
| 82 | + assert _ENGINE is not None |
| 83 | + kind, exe = _ENGINE |
| 84 | + cues_lit = " ".join(str(c) for c in cues) |
| 85 | + wins_lit = " ".join(str(w) for w in wins) |
| 86 | + # %.17g round-trips a double exactly; printf recycles the format per element. |
| 87 | + code = ( |
| 88 | + f"addpath('{_REF_DIR}');" |
| 89 | + f"cv = calibration_ref([{cues_lit}], [{wins_lit}], {base_rt!r}, {rt_change!r});" |
| 90 | + f"printf('%.17g\\n', cv);" |
| 91 | + ) |
| 92 | + cmd = [exe, "-batch", code] if kind == "matlab" else [exe, "--quiet", "--eval", code] |
| 93 | + proc = subprocess.run(cmd, capture_output=True, text=True, timeout=120) |
| 94 | + if proc.returncode != 0: |
| 95 | + raise RuntimeError( |
| 96 | + f"{kind} reference failed (rc={proc.returncode}):\n{proc.stderr}" |
| 97 | + ) |
| 98 | + return [float(line) for line in proc.stdout.split() if line.strip()] |
| 99 | + |
| 100 | + |
| 101 | +def _assert_parity( |
| 102 | + cues: list[int], wins: list[int], base_rt: float, rt_change: float |
| 103 | +) -> list[float]: |
| 104 | + py = _python_windows(cues, wins, base_rt, rt_change) |
| 105 | + ref = _octave_windows(cues, wins, base_rt, rt_change) |
| 106 | + assert len(py) == len(cues) |
| 107 | + assert ref == pytest.approx(py, abs=1e-12) |
| 108 | + return py |
| 109 | + |
| 110 | + |
| 111 | +# ── shared test data ────────────────────────────────────────────────────────── |
| 112 | + |
| 113 | +def _random_sequence(seed: int, n: int) -> tuple[list[int], list[int]]: |
| 114 | + """Interleaved (cue, win) trials spanning all 6 cues, reproducible by seed.""" |
| 115 | + import random |
| 116 | + |
| 117 | + rng = random.Random(seed) |
| 118 | + cue_ids = sorted(_CUE_TO_PM) |
| 119 | + cues = [rng.choice(cue_ids) for _ in range(n)] |
| 120 | + wins = [rng.randint(0, 1) for _ in range(n)] |
| 121 | + return cues, wins |
| 122 | + |
| 123 | + |
| 124 | +# all-wins / all-misses runs on a single cue, plus the 0.66 boundary case |
| 125 | +_ALL_WINS = ([5] * 8, [1] * 8) |
| 126 | +_ALL_MISS = ([2] * 8, [0] * 8) |
| 127 | +_BOUNDARY = ([4] * 51, [1] * 33 + [0] * 17 + [1]) # ratio hits exactly 33/50 = 0.66 |
| 128 | + |
| 129 | +_PARAMS = [(config.BASE_RT_S, config.RT_CHANGE_S), (config.BASE_RT_PRACTICE_S, 1 / 60.0)] |
| 130 | + |
| 131 | + |
| 132 | +# ── tests ───────────────────────────────────────────────────────────────────── |
| 133 | + |
| 134 | +@requires_engine |
| 135 | +@pytest.mark.parametrize("base_rt,rt_change", _PARAMS) |
| 136 | +@pytest.mark.parametrize("seed", [0, 1, 7, 42, 123]) |
| 137 | +def test_parity_random_sequences(base_rt, rt_change, seed): |
| 138 | + cues, wins = _random_sequence(seed, n=120) |
| 139 | + windows = _assert_parity(cues, wins, base_rt, rt_change) |
| 140 | + # guard against a vacuous pass: adaptation must actually have moved the window |
| 141 | + assert any(w != pytest.approx(base_rt) for w in windows) |
| 142 | + |
| 143 | + |
| 144 | +@requires_engine |
| 145 | +@pytest.mark.parametrize("base_rt,rt_change", _PARAMS) |
| 146 | +@pytest.mark.parametrize("cues,wins", [_ALL_WINS, _ALL_MISS, _BOUNDARY]) |
| 147 | +def test_parity_edge_cases(base_rt, rt_change, cues, wins): |
| 148 | + _assert_parity(list(cues), list(wins), base_rt, rt_change) |
| 149 | + |
| 150 | + |
| 151 | +@requires_engine |
| 152 | +def test_parity_check_is_discriminating(): |
| 153 | + """A deliberately wrong Python result must NOT match the reference — proves |
| 154 | + the comparison can actually fail, so a passing parity run is meaningful.""" |
| 155 | + cues, wins = _random_sequence(seed=99, n=40) |
| 156 | + ref = _octave_windows(cues, wins, config.BASE_RT_S, config.RT_CHANGE_S) |
| 157 | + perturbed = list(ref) |
| 158 | + perturbed[-1] += 0.020 # one trial off by a full rt_change step |
| 159 | + assert perturbed != pytest.approx(ref, abs=1e-12) |
0 commit comments