Skip to content

Array-API fallback for Combiner.sigma_clipping when the namespace is not numpy - #1001

Open
mwcraig wants to merge 4 commits into
astropy:mainfrom
mwcraig:sigma-clipping-array-api
Open

Array-API fallback for Combiner.sigma_clipping when the namespace is not numpy#1001
mwcraig wants to merge 4 commits into
astropy:mainfrom
mwcraig:sigma-clipping-array-api

Conversation

@mwcraig

@mwcraig mwcraig commented Aug 27, 2026

Copy link
Copy Markdown
Member

Part of #929 (the Combiner.sigma_clipping half). Companion PR: #1000 (the sigma_func half), built off the same main; the two touch different files and merge cleanly in either order (the only overlap is adjacent bullets in CHANGES.rst / docs/array_api.rst).

Combiner.sigma_clipping wrapped astropy.stats.sigma_clip unconditionally, which converts the data to numpy and fails outright on array-api-strict (3 of the 13 remaining strict failures).

Policy (numpy is untouched): when the namespace is numpy, the data still go to astropy.stats.sigma_clip exactly as before — masked=True, all extra keyword arguments forwarded, compiled fast path intact. Every other namespace uses a new private _sigma_clip_mask, an iterative clip written purely in terms of the array API that reproduces astropy's mask: the last iteration's bounds applied to the original data, non-finite values always masked, None/0 thresholds meaning 3 as in astropy. The 'median'/'mean'/'std'/'mad_std' string options map to the same tiered (bottleneck → native → fallback) NaN-aware reductions the combination methods use (new private _nanmadstd for 'mad_std'); astropy-only keyword arguments (grow, masked, return_bounds, …) raise TypeError off numpy. An integer maxiters runs with no host synchronisation, so on dask the clip stays one graph.

Latent bug fixed: axis, copy and maxiters were read with kwd.get and then forwarded again through **kwd, so passing any of them raised TypeError: got multiple values for keyword argument. They are popped now.

Docstring fixes: func defaults to 'mean', not 'median'; None thresholds are treated as 3 (as astropy does) rather than "no rejection" — the docstring still described the pre-2.4 behaviour. The None/0/negative threshold semantics are a behaviour question, not part of this PR — see #1002.

Tests: the three test_combiner_sigmaclip_* tests are parametrized over a backend-generic MAD callable and "mad_std" and no longer use numpy-only .all(). A differential test runs _sigma_clip_mask on every backend against astropy.stats.sigma_clip(..., masked=True, copy=False).mask over 5 data sets (outliers, NaN/inf incl. an all-NaN column, zero spread, a column that collapses mid-loop, ints) × 3 cenfuncs × 3 stdfuncs × 4 threshold pairs × maxiters ∈ {1, 3, None}; plus tests for silence on all-NaN slices, None/0 → 3, final-bounds-not-union semantics, numpy scalar thresholds, bad axes/strings, _nanmadstd vs mad_std, a spy proving astropy is called iff numpy, the axis/maxiters regression, astropy-only kwargs, the use_astropy deprecation, existing-mask preservation, and combine(sigma_clip=True) (previously not exercised at all). Design was prototype-verified against astropy over 864 structured cases (numpy + strict device1), 2400 random numpy trials and 243 cases each on jax and dask, 0 mismatches.

Docs: changelog entry; one bullet in docs/array_api.rst; docs/image_combination.rst now recommends func="median" over np.ma.median. The combiner.py sigma_clipping line leaves the escape baseline (confirmed by a full-suite dask regeneration; the diff is that one line).

Verified locally

🤖 Generated with Claude Code

https://claude.ai/code/session_01XGArZMnLPC1G7YP62vN18V

mwcraig and others added 2 commits August 26, 2026 20:38
Combiner.sigma_clipping wrapped astropy.stats.sigma_clip unconditionally,
which converts the data to numpy and fails outright on array-api-strict
(ccdproc#929). NumPy users must not be slowed down, so numpy data still go
to astropy exactly as before (masked=True, all extra keyword arguments
forwarded, compiled fast path intact). Every other namespace now uses
_sigma_clip_mask, an iterative clip written purely in terms of the array
API standard that reproduces astropy's mask: the last iteration's bounds
applied to the original data, non-finite values always masked, None/0
thresholds meaning 3 as in astropy. The 'median'/'mean'/'std'/'mad_std'
string options map to the same tiered NaN-aware reductions the combination
methods use (new private _nanmadstd for 'mad_std'); astropy-only keyword
arguments raise TypeError off numpy. An integer maxiters runs without any
host synchronisation, so the clip stays a single graph on dask.

Also fixes a latent bug: axis, copy and maxiters were read with kwd.get and
then forwarded again through **kwd, so passing any of them raised
"TypeError: got multiple values for keyword argument". They are popped now.

Docstring fixes: func defaults to 'mean' (not 'median'); None thresholds
are treated as 3 rather than "no rejection"; kwd semantics per namespace.
Docs: image_combination.rst recommends func="median" over np.ma.median;
array_api.rst lists the limitations off numpy. The sigma_clipping escape
baseline entry is deleted (confirmed by a full-suite dask regeneration).

Verified: the design was checked against astropy.stats.sigma_clip(...,
masked=True, copy=False).mask over 864 structured cases (numpy and
array-api-strict on device1), 2400 random numpy trials, and 243 cases each
on jax and dask, with 0 mismatches. The new differential test
(5 data sets x 3 cenfunc x 3 stdfunc x 4 thresholds x 3 maxiters) runs on
every backend. Full suite: numpy 1123 passed; array-api-strict 10 failed
(down from 13; all remaining are ccdproc#929 sigma_func, ccdproc#936 and
ccdproc#983), 0 xpassed; jax 1122 passed; dask 1116 passed.

Co-Authored-By: Claude Fable 5 <noreply@anthropic.com>
Claude-Session: https://claude.ai/code/session_01XGArZMnLPC1G7YP62vN18V
Co-Authored-By: Claude Fable 5 <noreply@anthropic.com>
Claude-Session: https://claude.ai/code/session_01XGArZMnLPC1G7YP62vN18V
@codecov

codecov Bot commented Aug 27, 2026

Copy link
Copy Markdown

Codecov Report

✅ All modified and coverable lines are covered by tests.
✅ Project coverage is 97.92%. Comparing base (57fd985) to head (8fb4719).
⚠️ Report is 1 commits behind head on main.

Additional details and impacted files
@@            Coverage Diff             @@
##             main    #1001      +/-   ##
==========================================
+ Coverage   97.78%   97.92%   +0.13%     
==========================================
  Files           9        9              
  Lines        1808     1878      +70     
==========================================
+ Hits         1768     1839      +71     
+ Misses         40       39       -1     
Flag Coverage Δ
dask 97.11% <100.00%> (+0.16%) ⬆️
jax 97.27% <100.00%> (+0.15%) ⬆️
numpy 97.81% <100.00%> (+0.13%) ⬆️

Flags with carried forward coverage won't be shown. Click here to find out more.

☔ View full report in Codecov by Harness.
📢 Have feedback on the report? Share it here.

🚀 New features to boost your workflow:
  • ❄️ Test Analytics: Detect flaky tests, report on failures, and find test suite problems.

mwcraig and others added 2 commits August 26, 2026 20:59
Record in the test file why _sigma_clip_mask is tested by comparison
against astropy over a full cross product rather than a curated list
(the corners of astropy's behaviour only appear for particular
combinations of center, deviation, thresholds and maxiters), what each
data set is for, and what the grid costs. Comments only; no behaviour
change.

Co-Authored-By: Claude Fable 5 <noreply@anthropic.com>
Claude-Session: https://claude.ai/code/session_01XGArZMnLPC1G7YP62vN18V
astropy 8.1 (astropy#19858) changed the python loop that sigma_clip takes
for a callable cenfunc/stdfunc: with masked=True it now masks the union of
every iteration's rejections regardless of copy, where it used to apply the
final bounds to the data for copy=False. Its compiled path (string options)
still applies the final bounds, which is what _sigma_clip_mask reproduces.

Derive the reference mask from return_bounds=True, which both astropy paths
and versions agree on, and keep checking astropy's own mask against it on
the compiled path. Docstrings note the corner where astropy >= 8.1's numpy
path can differ for a callable.

Co-Authored-By: Claude Fable 5 <noreply@anthropic.com>
Claude-Session: https://claude.ai/code/session_01JnWCg95xE93SbhME52jxGJ
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

1 participant