stireg

stireg is a Python framework for retrospective analysis of changes in localized, scale-resolved second-order summaries of spatio-temporal event patterns. It separates pair-score construction, information diagnostics, single-change inference, exploratory multiple-change segmentation, and visualization. It is independent of ST-SCKM and does not use a clustering objective, centroid update, or graph penalty.

The statistical target is precise: changes in intensity-reweighted event-pair projections indexed by prespecified spatial cells and spatial-temporal lag bins. The package detects changes in that discretized summary, not in every aspect of an underlying point-process mechanism.

Current statistical scope

The implemented software supports one unmarked continuous event pattern in a rectangular planar window. Both events in each pair must lie in the same temporal block, cells must be inside the window eroded by the maximum spatial lag, and intensities must be known or conditionally estimated at event locations. The confirmatory procedure tests one retrospective mean change in the complete score vector. Recursive binary segmentation can return several candidate boundaries, with explicitly segmentwise rather than family-wise error semantics.

The following are not yet validated: arbitrary polygonal windows, mark-specific scores, simultaneous affected-support confidence sets, automatic scientific scale selection, and unrestricted cross-block dependence. These are research tasks, not hidden package options.

Installation

Install the reviewed release from PyPI:

python -m pip install stireg==0.2.1

For a submission review without using PyPI, install the attached source distribution from the directory containing it:

python -m pip install stireg-0.2.1.tar.gz

Create an isolated environment and install the package from this checkout:

python3.12 -m venv .venv
source .venv/bin/activate
python -m pip install --upgrade pip
python -m pip install -e ".[test,applications]"

Alternatively, create the recorded Conda environment:

conda env create -f environment.yml
conda activate stireg-research

All commands in the replication section below are run from the repository root, the directory containing run_all.py, pyproject.toml, and requirements.txt.

Quick start

import numpy as np

from stireg import (
    RectWindow,
    block_pair_scores,
    information_report,
    regular_cells,
    simulate_snapshot_cluster_change,
    multiple_change_test,
    single_change_test,
)

window = RectWindow(0, 1, 0, 1, 0, 41)
simulation_blocks = np.arange(0, 42)
analysis_blocks = np.arange(1, 42)

events = simulate_snapshot_cluster_change(
    window,
    block_edges=simulation_blocks,
    change_block=21,
    target_rate=200,
    parent_rate_after=20,
    cluster_sigma_after=0.02,
    mean_offspring_after=10,
    random_state=20260928,
)

interior = window.eroded(0.15)
x_edges, y_edges = regular_cells(interior, nx=2, ny=2)
scores = block_pair_scores(
    events,
    block_edges=analysis_blocks,
    x_edges=x_edges,
    y_edges=y_edges,
    spatial_lags=[0, 0.05, 0.10, 0.15],
    temporal_lags=[0, 0.5, 1.0],
)

diagnostics = information_report(scores)
result = single_change_test(
    scores.matrix(),
    min_segment=8,
    n_bootstrap=499,
    bootstrap_block_length=1,
    random_state=20260928,
)

print(result.change_index, result.p_value)
print(diagnostics.n_low_information, diagnostics.n_components)
print(result.summary())
print(scores.component(result.component_index))

multiple = multiple_change_test(
    scores.matrix(), min_segment=8, max_changes=3,
    n_bootstrap=499, random_state=20260928,
)
print(multiple.summary())

# Optional: python -m pip install "stireg[plot]"
result.plot(scores.matrix())

The result reports the estimated split in units of temporal blocks, the maximum statistic, its bootstrap p-value, and the component with the largest standardized contrast. The component is diagnostic and is not a multiplicity-adjusted affected-support estimate.

Mathematical overview

For events (z_i=(s_i,t_i)) with intensity (\lambda(z_i)), each block-cell-lag score has the form

[ Y_{bp}=\frac{1}{E_{bp}} \sum_{i\ne j} \frac{f_{bp}(z_i,z_j)}{\lambda(z_i)\lambda(z_j)}-1, ]

where (f_{bp}) selects one temporal block, interior spatial cell, spatial annulus, and positive temporal-lag bin. The denominator (E_{bp}) is an exact geometric exposure for the restricted rectangular design. The second-order Campbell formula shows that the population mean is a bin average of the pair correlation minus one, so first-order intensity cancels when the supplied intensity is correct.

The retrospective test computes standardized CUSUM contrasts over candidate split times and components, then calibrates their maximum with whole-block multipliers. Under a Poisson reference, within-block pair construction prevents event sharing across temporal blocks. It does not guarantee independence for self-exciting or Cox processes.

API overview

• RectWindow and EventPattern: validated event data and observation geometry.
• regular_cells: regular rectangular cell boundaries.
• block_pair_scores: indexed pair enumeration, exact exposure, scores, and raw pair counts.
• ScoreCube.component: physical cell and lag intervals for any flattened component index.
• information_report: transparent sparsity diagnostics based only on pair counts.
• single_change_test: one-change maximum-CUSUM test and multiplier calibration.
• multiple_change_test: transparent recursive segmentation for several candidate changes.
• SingleChangeResult.summary() and .plot(): practitioner-facing reporting and diagnostics.
• simulate_piecewise_poisson: known-intensity piecewise Poisson simulation.
• simulate_snapshot_cluster and simulate_snapshot_cluster_change: independent clustered-snapshot null and interaction-change generators.
• thin_pattern: independent thinning with retained-process intensity updates.

All public functions include type hints, validation, and deterministic behavior when a random seed is supplied.

Validation status

The initial 200-replicate known-intensity experiment produced rejection rates of 0.020 under stationary Poisson data and 0.025 under a pure intensity change, with detection 0.890 for a Poisson-to-cluster change. A separate 200-replicate robustness experiment produced rejection 0.045 under stationary clustered snapshots and 0.010 under an independent retention change. Detection was 0.115 at low event density and 0.920 at high density. These values describe the recorded designs only.

Real-data acquisition

The USGS application is reproducible with:

python scripts/acquire_usgs_ridgecrest.py
python scripts/preprocess_usgs_ridgecrest.py

The query, retrieval time, checksums, coordinate reference systems, and attrition are recorded under data/. Global Fire Atlas acquisition uses the versioned Zenodo record DOI 10.5281/zenodo.11400062. The large archive remains external to any JSS upload.

Reproducing the article

Create a clean environment from the repository root:

python3.12 -m venv .venv
source .venv/bin/activate
python -m pip install --upgrade pip
python -m pip install -r requirements.txt
python -m pip install .

Run unit tests first:

python -m pytest

Regenerate the tests, external benchmark, targeted sensitivity analysis, real-data analyses, manuscript tables, and figures with the single replication entry point:

python run_all.py

To rerun all four 200-replicate calibrated simulation studies before rebuilding the downstream artifacts, use python run_all.py --full. The default mode uses the committed calibrated simulation outputs so that routine replication does not silently overwrite a long-running experiment.

More complete examples are in docs/quickstart.md, with the public interface summarized in docs/api.md.

Citation and license

Author information, ORCID, version, and citation guidance are in CITATION.cff. The software is distributed under the MIT License. MIT is GPL-compatible and therefore satisfies the JSS software licensing requirement. Dataset licenses and citations remain separate and are recorded with each acquisition script.