Metadata-Version: 2.5
Name: regscan
Version: 0.1.2
Summary: Regression-based scan statistics for interval anomaly detection in smoothly varying 1D signals
Project-URL: Homepage, https://github.com/BeardyMan37/regscan
Project-URL: Repository, https://github.com/BeardyMan37/regscan
Project-URL: Paper, https://arxiv.org/abs/2608.22201
Author: Gazi Abdur Rakib
License-Expression: BSD-3-Clause
License-File: LICENSE
Keywords: anomaly detection,change point,kernel regression,radio astronomy,scan statistics
Classifier: Development Status :: 3 - Alpha
Classifier: Intended Audience :: Science/Research
Classifier: Programming Language :: Python :: 3
Classifier: Topic :: Scientific/Engineering :: Astronomy
Classifier: Topic :: Scientific/Engineering :: Information Analysis
Requires-Python: >=3.10
Requires-Dist: numba>=0.57
Requires-Dist: numpy>=1.22
Provides-Extra: dev
Requires-Dist: mypy; extra == 'dev'
Requires-Dist: pytest-cov; extra == 'dev'
Requires-Dist: pytest>=7; extra == 'dev'
Requires-Dist: ruff; extra == 'dev'
Description-Content-Type: text/markdown

# regscan

Regression-based scan statistics for detecting interval anomalies in smoothly
varying 1-D signals.

Given a signal and a function family `F`, the score of an interval `I = [a, b]`
is

```
S(I) = 1 - (SR_I + SR_O) / SR_A
```

where `SR_A`, `SR_I` and `SR_O` are the sums of squared residuals from fitting
`F` to the whole signal, to the inside of `I`, and to the outside. The score
approaches 0 when splitting the signal explains nothing and approaches 1 when
it explains everything. A scan returns the highest-scoring interval.

## Installation

```bash
pip install regscan
```

Requires Python 3.10 or later. The only dependencies are numpy and numba.

The kernel scan is a tight scalar loop, which numba compiles well and numpy
vectorizes poorly, so it is JIT-compiled. The first call in a process pays a
one-off compilation cost; numba caches the result on disk thereafter.

## Quick start

```python
import numpy as np
import regscan

t = np.linspace(0, 1, 300)
x = 4.7 + 0.10 * t + 0.05 * np.random.default_rng(0).normal(0, 1, 300)
x[140:170] -= 0.25                      # the anomaly

res = regscan.scan(x, method="nwkr_gaussian", w=12)
res.score, res.a, res.b                 # 0.248, 140, 169
res.width_frac, res.at_edge()           # 0.10, False
```

On the same signal with no anomaly planted, `nwkr_gaussian` scores 0.036 while
`mean` scores 0.185. The constant family is responding to the ramp, not to an
anomaly.

## Notation

| symbol | meaning |
|---|---|
| `n` | length of the signal, in samples |
| `w` | window width: the kernel bandwidth for `F_KR`, and the scale of local structure any family's fit can follow |
| `r` | range cap: the longest candidate interval considered, `3w` by default |
| `d` | polynomial degree for the `F_d` family (`poly_deg1` is `d = 1`) |
| `a`, `b` | inclusive start and end indices of a candidate interval |
| `I` | the candidate interval `[a, b]` |

`n` is fixed by the data. `w` and `r` are the parameters worth thinking about,
and both may be passed explicitly to `regscan.scan`.

## Methods

| method | family | model | cost |
|---|---|---|---|
| `mean` | `F_0` | constant | O(nr) |
| `poly_deg1` ... `poly_deg3` | `F_d` | degree-`d` polynomial | O(nr d³) |
| `nwkr_gaussian`, `nwkr_laplace` | `F_KR` | Nadaraya-Watson kernel regression | O(nrw) |
| `krr_gaussian`, `krr_laplace` | `F_KRR` | kernel ridge regression | O(n⁴w) |

`regscan.available()` lists them at run time.

`F_KR` is the method this package exists for. A weak family such as `F_0` or
`F_1` cannot represent a curved background, so it lowers the residual by
splitting the interval wherever the curvature is worst, which flags smooth
structure as an anomaly. A kernel fit follows that structure instead, so the
structure enters `SR_A`, `SR_I` and `SR_O` alike and cancels out of the score.

## Choosing `w` and `r`

`w` sets the scale of structure the fit can follow. Set it too small and the
kernel reproduces the anomaly itself, which cancels it from the score; set it
too large and the fit cannot follow the background, which is the failure mode
of the weaker families.

`r` is the longest interval the scan will consider. An anomaly wider than `r`
cannot be returned at all, and cost falls linearly with `r`, so it is the
parameter to reach for when the approximate width of the feature is known.

```python
regscan.scan(x, method="nwkr_gaussian", w=12, r=90)
```

When omitted, `w` defaults to `max(3, n // 16)` and `r` to `3 * w`. These
defaults exist so that a scan runs unattended; they are not recommendations.

Sweeping `w` across a range and checking whether `(a, b)` holds steady is a
cheap way to distinguish a resolved feature from an artifact of the bandwidth.
A real interval stays put, whereas one that tracks `w` is measuring the kernel.

## Configuration

```python
from regscan import ScanConfig

cfg = ScanConfig(
    kernel="laplace",     # "gaussian" or "laplace"
    buffer=24,            # exclude this many samples at each end
    min_width=0.01,       # interval length bounds, as a fraction of n
    max_width=0.25,
)
regscan.scan(x, method="nwkr_gaussian", w=16, config=cfg)
```

`ScanConfig` is immutable and passed explicitly. Nothing is stored in module
globals, so scanning several methods in one process cannot leak settings
between them.

`buffer` excises samples rather than merely forbidding interval placement
there, so the fit, `SR_A` and the outside residuals are all computed on the
trimmed signal. Because a non-zero buffer suppresses detections at the ends of
the signal, give every family the same value when comparing them.

`max_width` is worth capping. As an interval approaches `n/2` the inside and
outside become comparable in size, the statistic stops discriminating, and the
maximizer drifts toward whatever split best absorbs slow curvature.

## Super-resolution

`F_KR` can block-mean the signal, scan the shorter version, then search the
original samples around the winning blocks to recover exact endpoints.

```python
from regscan import ScanConfig

regscan.scan(x, method="nwkr_gaussian", w=100,
             config=ScanConfig(super_resolution=4))     # or "auto"
```

At `n = 1600` with `w = 100`:

| factor | time | interval |
|---|---|---|
| 1 (exact) | 3.11 s | (700, 819) |
| 2 | 0.44 s | (700, 819) |
| 4 | 0.08 s | (700, 819) |
| 8 | 0.02 s | (704, 815) |

A factor of 4 runs 39 times faster and returns the same interval. A factor of 8
does not, which is the trade-off: the coarse pass locates each endpoint only to
within a block of `factor` samples, and a wrong block gives a wrong answer.

Refinement therefore searches one block either side of each coarse endpoint.
This matters more than it may appear. Block means smooth an anomaly's edges, so
the coarse pass selects a neighboring block often enough that confining the
search to the winning block alone recovered the exact interval in only 20 of 40
trials at factor 4. Including the neighbors recovered all 40.

It remains an approximation. Verify against `super_resolution=1` on a sample of
your own data before relying on it.

Passing `"auto"` derives the factor from the signal length: 1 below 450
samples, then doubling at 900, 1800 and so on. `sr_cap` bounds it, which
matters for narrow features, since an interval must survive decimation to be
found. The default is `1`, which scans every sample exactly.

## Performance

One full scan with `w = n // 16`, after JIT warm-up:

| n | `mean` | `poly_deg1` | `nwkr_gaussian` |
|---|---|---|---|
| 100 | 3 ms | 44 ms | 4 ms |
| 200 | 12 ms | 181 ms | 18 ms |
| 400 | 51 ms | 754 ms | 83 ms |
| 800 | 187 ms | 2.7 s | 0.45 s |
| 1600 | 738 ms | 11.0 s | 3.2 s |

At `n = 800`, `F_KR` is six times faster than `F_1` despite fitting a far
richer model. That is the practical case for it: the polynomial family pays
O(d³) for every candidate interval, whereas the kernel family pays O(r) to
extend one.

Each family reaches its stated complexity by carrying state rather than
refitting. `mean` scores an interval in O(1) from prefix sums of `y` and `y²`.
`poly_deg1` uses prefix sums of the moments `t^p` and `t^p y`, so a fit costs
O(d³) to solve regardless of how many samples the interval spans. `F_KR` grows
an interval one sample at a time and updates in O(r): the inside buffer and
`sse_in`, the `nin` and `din` arrays holding the inside points' kernel
contribution to every index, and `sse_out`, obtained from those by subtraction
from the all-points totals. Since `sse_out` is adjusted rather than recomputed,
it is refreshed exactly on a fixed cadence to keep floating-point error from
accumulating.

Super-resolution reduces these times substantially again on longer signals.

## Comparing scores

Scores are comparable within a family but not across families. Each family
divides by its own `SR_A`, and a kernel fit has a smaller `SR_A` than a
constant fit before any interval is chosen, so the same interval scores
differently under `F_0` and `F_KR`. Compare families by rank, by whether they
agree on the interval, or by the contrast between anomalous and clean signals,
rather than by absolute score.

## Citation

Rakib et al., *Efficient Regression Models for Scan Statistics*,
[arXiv:2608.22201](https://arxiv.org/abs/2608.22201) (2026).

```bibtex
@misc{rakib2026efficientregressionmodelsscan,
      title={Efficient Regression Models for Scan Statistics},
      author={Gazi Abdur Rakib and Tristan Ashton and Ryan A. Loomis and Brian S. Mason and Eric J. Murphy and Ci Xue and Jeff M. Phillips},
      year={2026},
      eprint={2608.22201},
      archivePrefix={arXiv},
      primaryClass={stat.ME},
      url={https://arxiv.org/abs/2608.22201},
}
```

## License

BSD 3-Clause. See [LICENSE](LICENSE).