Skip to content

Commit baace01

Browse files
yilibinbinAaronhaohaoclaude
authored
fix: precision-correctness audit findings (F6/F11/F2/F1) (#74)
* fix: precision-correctness audit findings F6/F11/F2/F1 F6 (cli/main.py) — CLI batch `calc` serialized the extrapolation limit with bare str(mpf) outside any precision_guard, truncating a requested high-precision (e.g. 50-digit) result to the ambient mp.dps (~15). Wrap compute + serialization in precision_guard(job.precision) and render via mp.nstr at that precision. F11 (app_desktop/workspace_controller.py) — workspace save persisted config.common (mpmath_precision, uncertainty/display digits, scientific) and config.latex (input digits, group size, dcolumn) but restore never read them back, silently resetting compute-affecting precision to defaults on reload. Add _restore_common_config, called from _restore_workspace_contents. F2 (datalab_core/mcmc_refine.py) — the Gaussian log-probability divided the already sigma-weighted residual sum (Σ rᵢ²/σᵢ²) by rmse² again, double-scaling the likelihood and corrupting credible-interval WIDTHS for any weighted fit. Extract _gaussian_log_probability: -0.5·residuals_sq when weighted (weights already 1/σ²), keep the /rmse² plug-in variance only in the unweighted branch. F1 (fitting/hp_fitter.py) — _compute_covariance took mp.sqrt of a covariance diagonal guarded only for NaN, so a negative diagonal (finite-precision indefinite inverse) yielded a complex mpc that crashed combine_error_components. Add _error_from_variance (negative/NaN/Inf → NaN) and flag negative variances in the covariance warning, so the fit degrades gracefully instead of crashing. Each fix is TDD'd (pure-function unit tests + workspace round-trip + CLI calc high-precision). Neighborhood regression green (218 tests). * fix(workspace): restore LaTeX caption + engine on reload; stop mp.dps leak in test Review follow-ups on the F11 precision/round-trip fix: - workspace_controller._restore_common_config captured config.latex.caption and config.latex.engine at save (lines 1125-1126) but never read them back, so a saved caption reset to empty and the TeX engine reverted to the "tectonic" default on reload — the same "captured but never restored" class F11 set out to fix, just incomplete for these two fields. Restore both, and extend the round-trip test to set and assert them (it previously only covered the precision/digit fields, so the gap was uncaught). - test_cli_batch set mp.dps = 15 to simulate ambient CLI precision but never restored it; since mp.dps is process-global that leaked into later tests. Wrap the mutation in try/finally. ruff clean; 148 workspace/cli tests pass. Round-trip test is RED without the restore fix (caption not preserved), GREEN with it. Co-Authored-By: Claude Opus 4.8 <noreply@anthropic.com> * fix(mcmc): index likelihood_weights via the None-check so mypy --strict narrows The F2 refactor hoisted the None-check into a `weighted` boolean and then indexed likelihood_weights[index] guarded by `if weighted`. mypy --strict cannot use a separate local flag to narrow `Sequence[Any] | None`, so it reported "Value of type 'Sequence[Any] | None' is not indexable", failing the core-layer strict gate (test_mypy_strict_zero_errors_on_full_core_layer) — main passed it. Index against the direct `if likelihood_weights is not None` instead (mypy narrows it); keep the `weighted` flag for the _gaussian_log_probability call. Behavior is identical (weighted is True exactly when likelihood_weights is not None). mypy strict gate + MCMC log-prob tests pass; ruff clean. Co-Authored-By: Claude Opus 4.8 <noreply@anthropic.com> --------- Co-authored-by: fanghao <fanghaoaa@gmail.com> Co-authored-by: Claude Opus 4.8 <noreply@anthropic.com>
1 parent 75af565 commit baace01

8 files changed

Lines changed: 252 additions & 12 deletions

File tree

app_desktop/workspace_controller.py

Lines changed: 31 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -730,6 +730,36 @@ def _restore_root_config(window: Any, config: Any, *, restore_constants: bool =
730730
_restore_display_units_controls(window, "root", units, include_constants=True)
731731

732732

733+
def _set_checked_if(window: Any, attr: str, value: Any) -> None:
734+
widget = getattr(window, attr, None)
735+
if widget is not None and value is not None and hasattr(widget, "setChecked"):
736+
widget.setChecked(bool(value))
737+
738+
739+
def _restore_common_config(window: Any, common: Any, latex: Any) -> None:
740+
"""Restore the mode-independent common + latex settings. These are captured
741+
at save (config.common / config.latex) but were never read back, silently
742+
resetting compute-affecting precision and display options to defaults on
743+
reload (audit F11). mpmath_precision in particular changes the numerical
744+
result, so it must survive a round-trip."""
745+
if isinstance(common, dict):
746+
_set_value(getattr(window, "mpmath_precision_spin", None), common.get("mpmath_precision"))
747+
_set_value(getattr(window, "uncertainty_digits_spin", None), common.get("uncertainty_digits"))
748+
_set_value(getattr(window, "display_digits_spin", None), common.get("display_digits"))
749+
_set_checked_if(window, "generate_latex_checkbox", common.get("generate_latex"))
750+
_set_checked_if(window, "generate_plots_checkbox", common.get("generate_plots"))
751+
_set_checked_if(window, "verbose_checkbox", common.get("verbose"))
752+
_set_checked_if(window, "scientific_checkbox", common.get("display_scientific"))
753+
if isinstance(latex, dict):
754+
_set_value(getattr(window, "latex_input_precision_spin", None), latex.get("input_digits"))
755+
_set_value(getattr(window, "latex_group_size_spin", None), latex.get("group_size"))
756+
_set_checked_if(window, "dcolumn_checkbox", latex.get("use_dcolumn"))
757+
_set_checked_if(window, "caption_checkbox", latex.get("use_caption"))
758+
_set_text(getattr(window, "output_file_edit", None), str(latex.get("output_path") or ""))
759+
_set_text(getattr(window, "caption_edit", None), str(latex.get("caption") or ""))
760+
_set_combo_data(getattr(window, "latex_engine_combo", None), str(latex.get("engine") or "tectonic"))
761+
762+
733763
def _restore_extrapolation_config(window: Any, config: Any) -> None:
734764
if not isinstance(config, dict):
735765
return
@@ -1910,6 +1940,7 @@ def _restore_workspace_contents(window: Any, manifest: dict[str, Any], attachmen
19101940
_restore_data_section(window, unified_constants, constants=True)
19111941

19121942
config = workspace.get("config") or {}
1943+
_restore_common_config(window, config.get("common"), config.get("latex"))
19131944
_restore_extrapolation_config(window, config.get("extrapolation"))
19141945
_restore_error_config(window, config.get("error"))
19151946
_restore_statistics_config(window, config.get("statistics"))

cli/main.py

Lines changed: 20 additions & 7 deletions
Original file line numberDiff line numberDiff line change
@@ -158,29 +158,42 @@ def _run_calc(job: BatchJob) -> dict[str, object]:
158158
method hard-defaults to ``wynn_epsilon``; future iterations can
159159
expose ``method`` via the YAML schema.
160160
"""
161+
from mpmath import mp
162+
161163
from extrapolation_methods.accelerators import (
162164
SequenceAcceleratorConfig,
163165
apply_sequence_accelerator,
164166
)
167+
from shared.precision import precision_guard
165168

166169
xs, ys = _read_xy_csv(job.data_path)
167170
if len(ys) < 3:
168171
raise ValueError(
169172
f"Job {job.name!r}: calc needs at least 3 points, got {len(ys)}"
170173
)
171-
config = SequenceAcceleratorConfig(precision=job.precision)
172-
result = apply_sequence_accelerator("wynn_epsilon", ys, config)
174+
175+
def _serialize_mpf(value: object) -> object:
176+
# str(mpf) emits only ambient-mp.dps digits, silently truncating a
177+
# high-precision result; nstr at the requested precision preserves it.
178+
return mp.nstr(value, job.precision) if hasattr(value, "_mpf_") else value
179+
180+
# Wrap compute AND serialization in the precision guard so the accelerator
181+
# runs at, and the limit/metadata are rendered at, the requested precision
182+
# regardless of the ambient mp.dps the CLI process happens to be at.
183+
with precision_guard(job.precision):
184+
config = SequenceAcceleratorConfig(precision=job.precision)
185+
result = apply_sequence_accelerator("wynn_epsilon", ys, config)
186+
limit_text = _serialize_mpf(result.value)
187+
metadata = {k: _serialize_mpf(v) for k, v in (result.metadata or {}).items()}
188+
173189
return {
174190
"job_name": job.name,
175191
"operation": "calc",
176192
"method": "wynn_epsilon",
177193
"data_path": str(job.data_path),
178194
"precision": job.precision,
179-
"limit": str(result.value),
180-
"metadata": {
181-
k: (str(v) if hasattr(v, "_mpf_") else v)
182-
for k, v in (result.metadata or {}).items()
183-
},
195+
"limit": limit_text,
196+
"metadata": metadata,
184197
}
185198

186199

datalab_core/mcmc_refine.py

Lines changed: 24 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -84,6 +84,25 @@ def _evaluate_prediction(evaluator: Any, params: Mapping[str, mp.mpf], observati
8484
raise TypeError("MCMC evaluator does not accept supported argument shapes")
8585

8686

87+
def _gaussian_log_probability(residuals_sq: float, rmse: float, *, weighted: bool) -> float:
88+
"""Gaussian log-likelihood (up to an additive constant) for MCMC.
89+
90+
``residuals_sq`` is Σ wᵢ·rᵢ². When the fit is WEIGHTED the weights are already
91+
1/σᵢ² (see ``likelihood_weights_from_series``), so ``residuals_sq`` is the
92+
correct χ² = Σ rᵢ²/σᵢ² and the log-likelihood is simply ``-0.5·residuals_sq``.
93+
Dividing again by rmse² (the old code path) double-scaled the likelihood,
94+
acting as a constant "temperature" that rescaled every posterior variance and
95+
made the reported credible intervals wrong whenever σ/weights were supplied.
96+
97+
When UNWEIGHTED (weights all 1, data noise unknown), rmse² is the plug-in
98+
variance estimate, so the Gaussian exponent is ``-0.5·Σrᵢ²/rmse²`` — keep the
99+
division there.
100+
"""
101+
if weighted:
102+
return -0.5 * residuals_sq
103+
return -0.5 * residuals_sq / (rmse**2)
104+
105+
87106
def _estimate_rmse(targets: Sequence[Any], fit_result: Any) -> float:
88107
"""Estimate a strictly positive RMSE for MCMC proposal scaling."""
89108
residuals = fit_result.details.get("residuals") if fit_result.details else None
@@ -141,6 +160,8 @@ def refine_fit_with_mcmc(
141160
initial_guess = [float(fit_result.params[name]) for name in param_names] # float-bridge: emcee walker seeds
142161
rmse = _estimate_rmse(targets, fit_result)
143162

163+
weighted = likelihood_weights is not None
164+
144165
def _log_probability(theta: Sequence[object]) -> float:
145166
if not param_names or rmse <= 0:
146167
return float("-inf") # float-bridge: emcee log-prob sentinel
@@ -153,11 +174,13 @@ def _log_probability(theta: Sequence[object]) -> float:
153174
if not math.isfinite(pred):
154175
return float("-inf") # float-bridge: emcee log-prob sentinel
155176
residual = float(target) - pred # float-bridge: emcee log-prob
177+
# Index against the direct None-check (not the hoisted `weighted`
178+
# flag) so mypy --strict narrows Sequence | None to indexable.
156179
weight = float(likelihood_weights[index]) if likelihood_weights is not None else 1.0 # float-bridge: emcee log-prob
157180
residuals_sq += weight * (residual**2)
158181
if not math.isfinite(residuals_sq):
159182
return float("-inf") # float-bridge: emcee log-prob sentinel
160-
return -0.5 * residuals_sq / (rmse**2)
183+
return _gaussian_log_probability(residuals_sq, rmse, weighted=weighted)
161184
except (TypeError, ValueError, ArithmeticError, OverflowError, KeyError):
162185
return float("-inf") # float-bridge: emcee log-prob sentinel
163186

fitting/hp_fitter.py

Lines changed: 24 additions & 4 deletions
Original file line numberDiff line numberDiff line change
@@ -293,6 +293,21 @@ def _compute_statistics(
293293
)
294294

295295

296+
def _error_from_variance(variance: mp.mpf) -> mp.mpf:
297+
"""Parameter standard error = sqrt(variance), degrading to NaN for a
298+
non-finite or NEGATIVE variance. A negative covariance diagonal (finite-
299+
precision indefiniteness of the inverse) would otherwise make mp.sqrt return
300+
a complex mpc that later crashes combine_error_components with
301+
'cannot create mpf from mpc' — degrade to NaN so the fit surfaces a
302+
covariance warning instead of a hard failure (audit finding F1)."""
303+
try:
304+
if mp.isnan(variance) or mp.isinf(variance) or variance < 0:
305+
return mp.nan
306+
return mp.sqrt(variance)
307+
except (TypeError, ValueError):
308+
return mp.nan
309+
310+
296311
def _compute_covariance(
297312
model: ModelSpecification,
298313
params: dict[str, mp.mpf],
@@ -343,12 +358,17 @@ def _compute_covariance(
343358
noise = chi2 / dof if dof > 0 else mp.nan
344359
covariance = [[inv[i, j] * noise for j in range(k)] for i in range(k)]
345360
errors = {
346-
name: mp.sqrt(covariance[idx][idx])
347-
if not mp.isnan(covariance[idx][idx])
348-
else mp.nan
361+
name: _error_from_variance(covariance[idx][idx])
349362
for idx, name in enumerate(free_params)
350363
}
351-
if any(mp.isnan(val) or mp.isinf(val) for row in covariance for val in row):
364+
# Flag ill-conditioning: non-finite cells OR a negative variance on the
365+
# diagonal. A negative diagonal (finite-precision indefiniteness of the
366+
# inverse) is neither NaN nor Inf, so it would otherwise pass silently and
367+
# yield a complex sqrt that crashes downstream error combination.
368+
has_negative_variance = any(covariance[i][i] < 0 for i in range(k))
369+
if has_negative_variance or any(
370+
mp.isnan(val) or mp.isinf(val) for row in covariance for val in row
371+
):
352372
cov_warning = "协方差矩阵病态或奇异,参数不确定度可能不可靠。 / Covariance matrix is ill-conditioned or singular; parameter uncertainties may be unreliable."
353373
return covariance, errors, cov_warning
354374

tests/test_cli_batch.py

Lines changed: 41 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -57,6 +57,47 @@ def _batch_config(_data_csv: Path, tmp_path: Path) -> Path:
5757
return path
5858

5959

60+
def test_run_calc_preserves_requested_high_precision_in_limit(tmp_path: Path):
61+
"""CLI batch `calc` must serialize the extrapolated limit at the requested
62+
precision, not truncate it to the ambient mp.dps via bare str(mpf) outside a
63+
precision_guard (audit finding F6)."""
64+
from mpmath import mp
65+
66+
from cli.batch_config import BatchJob
67+
from cli.main import _run_calc
68+
69+
# A geometric series 1 + 1/2 + 1/4 + ... whose Wynn-eps limit is exactly 2,
70+
# but we assert the SIGNIFICANT-DIGIT COUNT of the serialized limit, which is
71+
# what truncation to dps=15 would cap. Use a slowly-varying convergent whose
72+
# accelerated value carries many digits: partial sums of sum 6/k^2 -> pi^2.
73+
data = tmp_path / "seq.csv"
74+
rows = ["n,y"]
75+
total = mp.mpf(0)
76+
with mp.workdps(80):
77+
for n in range(1, 40):
78+
total += mp.mpf(6) / mp.mpf(n) ** 2
79+
rows.append(f"{n},{mp.nstr(total, 60)}")
80+
data.write_text("\n".join(rows), encoding="utf-8")
81+
82+
# mp.dps is process-global; set a low ambient precision to simulate the CLI
83+
# before the per-job guard, then restore it so this mutation can't leak into
84+
# later tests (the per-job precision_guard is exactly what we're asserting).
85+
original_dps = mp.dps
86+
try:
87+
mp.dps = 15 # simulate the ambient CLI precision before per-job guard
88+
job = BatchJob(
89+
name="hp-calc", operation="calc", data_path=data,
90+
output_dir=tmp_path, model=None, precision=50,
91+
)
92+
result = _run_calc(job)
93+
finally:
94+
mp.dps = original_dps
95+
96+
limit = str(result["limit"])
97+
digits = len(limit.replace(".", "").replace("-", "").lstrip("0"))
98+
assert digits > 20, f"limit truncated to ambient dps: {limit!r} ({digits} digits)"
99+
100+
60101
def test_load_batch_config_round_trips(_batch_config: Path):
61102
"""Parse a batch YAML into a BatchConfig dataclass."""
62103
from cli.batch_config import load_batch_config
Lines changed: 40 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,40 @@
1+
"""Covariance error guard against a negative variance diagonal (audit F1).
2+
3+
A near-rank-deficient inverse can have a NEGATIVE diagonal (finite-precision
4+
indefiniteness). sqrt of that is a complex mpc, which passes an isnan-only guard
5+
and then crashes combine_error_components ('cannot create mpf from mpc'). The
6+
error must degrade to NaN so the fit surfaces a covariance warning instead.
7+
"""
8+
9+
from __future__ import annotations
10+
11+
from mpmath import mp
12+
13+
from fitting.hp_fitter import _error_from_variance, combine_error_components
14+
15+
16+
def test_error_from_negative_variance_is_nan_not_complex():
17+
result = _error_from_variance(mp.mpf("-1e-5"))
18+
assert mp.isnan(result)
19+
# Must NOT be a complex mpc (the pre-fix behavior).
20+
assert not isinstance(result, mp.mpc)
21+
22+
23+
def test_error_from_valid_variance_is_sqrt():
24+
assert _error_from_variance(mp.mpf("4")) == mp.mpf("2")
25+
26+
27+
def test_error_from_nan_and_inf_variance_is_nan():
28+
assert mp.isnan(_error_from_variance(mp.nan))
29+
assert mp.isnan(_error_from_variance(mp.inf))
30+
31+
32+
def test_nan_error_does_not_crash_combine_error_components():
33+
"""The NaN degradation flows through error combination without the
34+
'cannot create mpf from mpc' TypeError the negative-variance path caused."""
35+
stat = _error_from_variance(mp.mpf("-1e-5")) # NaN, not mpc
36+
stat_map, _sys_map, total_map = combine_error_components(
37+
{"a": mp.mpf("1")}, {"a": stat}, {}
38+
)
39+
assert mp.isnan(stat_map["a"])
40+
assert mp.isnan(total_map["a"])

tests/test_mcmc_log_probability.py

Lines changed: 35 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,35 @@
1+
"""MCMC Gaussian log-probability scaling (audit finding F2).
2+
3+
The weighted branch must NOT divide the already-weighted chi-square (Σ rᵢ²/σᵢ²)
4+
by rmse² again — that double-scaling corrupted the posterior width and the
5+
reported credible intervals whenever sigma/weights were supplied. The unweighted
6+
branch keeps the /rmse² plug-in variance. This pure function needs no emcee.
7+
"""
8+
9+
from __future__ import annotations
10+
11+
from datalab_core.mcmc_refine import _gaussian_log_probability
12+
13+
14+
def test_weighted_log_probability_is_not_divided_by_rmse_squared():
15+
residuals_sq = 4.0 # already Σ rᵢ²/σᵢ² when weighted
16+
rmse = 0.01 # small rmse would blow up the exponent if wrongly divided
17+
lp = _gaussian_log_probability(residuals_sq, rmse, weighted=True)
18+
# Correct weighted Gaussian: -0.5 * chi-square, independent of rmse.
19+
assert lp == -0.5 * residuals_sq
20+
# And crucially NOT the double-scaled value.
21+
assert lp != -0.5 * residuals_sq / (rmse**2)
22+
23+
24+
def test_weighted_log_probability_ignores_rmse():
25+
# Same chi-square with wildly different rmse must give the SAME log-prob.
26+
assert _gaussian_log_probability(4.0, 0.01, weighted=True) == _gaussian_log_probability(
27+
4.0, 1000.0, weighted=True
28+
)
29+
30+
31+
def test_unweighted_log_probability_keeps_rmse_normalization():
32+
residuals_sq = 4.0 # Σ rᵢ² (weights all 1)
33+
rmse = 2.0
34+
lp = _gaussian_log_probability(residuals_sq, rmse, weighted=False)
35+
assert lp == -0.5 * residuals_sq / (rmse**2)

tests/test_workspace_controller.py

Lines changed: 37 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -206,6 +206,43 @@ def test_workspace_controller_captures_manual_state(qtbot) -> None:
206206
assert bundle.attachments["attachments/plots/plot-001.png"] == PNG_1X1
207207

208208

209+
def test_workspace_round_trips_common_and_latex_precision_settings(qtbot) -> None:
210+
"""The common (mpmath_precision, uncertainty/display digits, scientific) and
211+
latex (input digits, group size) settings are saved but were never restored,
212+
silently resetting compute-affecting precision to defaults on reload
213+
(audit finding F11)."""
214+
from app_desktop.window import ExtrapolationWindow
215+
from app_desktop.workspace_controller import capture_workspace, restore_workspace
216+
217+
source = ExtrapolationWindow()
218+
qtbot.addWidget(source)
219+
source.mpmath_precision_spin.setValue(50)
220+
source.uncertainty_digits_spin.setValue(3)
221+
source.display_digits_spin.setValue(25)
222+
source.scientific_checkbox.setChecked(True)
223+
source.latex_input_precision_spin.setValue(40)
224+
source.latex_group_size_spin.setValue(5)
225+
# caption text + TeX engine are also captured at save; they must round-trip
226+
# too, or the F11 "captured but never restored" fix is incomplete.
227+
source.caption_edit.setText("Table 1: extrapolated limits")
228+
source.latex_engine_combo.setCurrentText("pdflatex")
229+
230+
bundle = capture_workspace(source, title="precision settings")
231+
232+
target = ExtrapolationWindow()
233+
qtbot.addWidget(target)
234+
restore_workspace(target, bundle.manifest, bundle.attachments)
235+
236+
assert target.mpmath_precision_spin.value() == 50
237+
assert target.uncertainty_digits_spin.value() == 3
238+
assert target.display_digits_spin.value() == 25
239+
assert target.scientific_checkbox.isChecked() is True
240+
assert target.latex_input_precision_spin.value() == 40
241+
assert target.latex_group_size_spin.value() == 5
242+
assert target.caption_edit.text() == "Table 1: extrapolated limits"
243+
assert target.latex_engine_combo.currentText() == "pdflatex"
244+
245+
209246
def test_workspace_preserves_raw_constants_text_view_draft(qtbot) -> None:
210247
from app_desktop.window import ExtrapolationWindow
211248
from app_desktop.workspace_controller import capture_workspace, restore_workspace

0 commit comments

Comments
 (0)