US Time-Series Example

Audience: Maintainers and reviewers. Stata-parity audit document.

This example documents the current audited reproduction slice for the Stata xtbreak COVID-19 time-series example. It covers the paper workflow from xtbreak estimate through Python-native post-estimation payloads for regime indicators, split designs, scatter data, and the single-break SSR path.

Source Anchors

  • Software paper: paper/arXiv-2110.14550v3/xtbreak_article.tex lines 444-470 define the US weekly COVID-19 model, the D.deaths and D.L(1/3).cases specification, the default automatic xtbreak workflow, and the conclusion that two breaks are selected at the 5% level.

  • Software paper: paper/arXiv-2110.14550v3/xtbreak_article.tex lines 481-487 define the estat split follow-up regression and long-run multiplier interpretation.

  • Software paper: paper/arXiv-2110.14550v3/xtbreak_article.tex lines 537-549 give the explicit xtbreak estimate d.deaths d.L(1/3).cases, breaks(2) and estat scatter estimation path.

  • Stata reference: xtbreak-2.2/examples/EmpiricalExample_1.do lines 20-45 run the corresponding automatic command, estat split, hypothesis tests, two-break estimate, and scatter post-estimation path.

  • Stata help: xtbreak-2.2/ado/xtbreak_estimate.sthlp lines 62-79 document SSR-minimizing dynamic programming and the partial dynamic-program path for variables without breaks.

  • Stata help: xtbreak-2.2/ado/xtbreak_estimate.sthlp lines 100-141 document e(breaks), e(CI), e(SSRvec), estat split, estat scatter, and estat ssr.

Audited Fixture Path

The current Python example hydrates EstimateResult from Stata oracle JSON files:

  • artifacts/verification/stata/phase1_us_ts_estimate_breaks2.json

  • artifacts/verification/stata/phase1_us_ts_estimate_breaks1_ssr_path.json

The empirical reproduction test also reads xtbreak-2.2/data/US.dta and reconstructs the Stata sample with 79 observations.

Two-Break Estimate

import json
from pathlib import Path

import xtbreak

payload = json.loads(
    Path("artifacts/verification/stata/phase1_us_ts_estimate_breaks2.json").read_text(
        encoding="utf-8"
    )
)
results = payload["results"]

estimate = xtbreak.EstimateResult.from_break_matrices(
    command=payload["command"],
    break_matrix=results["breaks"],
    ci_matrix=results["ci"],
    ssr=results["ssr"],
    metadata={"scenario_id": payload["scenario_id"]},
)

assert estimate.break_table() == [
    {
        "break": 1,
        "index": 15,
        "time_value": 3141,
        "ci_index_low": 14,
        "ci_index_high": 16,
        "ci_time_low": 3140,
        "ci_time_high": 3142,
    },
    {
        "break": 2,
        "index": 45,
        "time_value": 3171,
        "ci_index_low": 44,
        "ci_index_high": 46,
        "ci_time_low": 3170,
        "ci_time_high": 3172,
    },
]
assert estimate.ssr == 49.0738641567431

The paper output rounds the SSR display; the automated reproduction test compares the Python partial-break search to the Stata oracle at abs=1e-10 and to the paper display at abs=0.005.

Post-Estimation Payload Crosswalk

Phase 23 adds a direct EstimateResult payload crosswalk for the same 79-observation US sample. The fitted sample is hydrated with dependent name D.deaths, break variables D.L.cases, D.L2.cases, and D.L3.cases, and the Stata weekly time axis 3127..3205. The scatter step is scatter_plot_data("D.L.cases").

indicator = estimate.regime_indicators()
split = estimate.split_design()
scatter = estimate.scatter_plot_data("D.L.cases")

assert [indicator["values"].count(regime) for regime in (1, 2, 3)] == [15, 30, 34]
assert split["column_names"][:3] == [
    "D.L.cases_regime_1",
    "D.L.cases_regime_2",
    "D.L.cases_regime_3",
]
assert [len(trace["row_ids"]) for trace in scatter["traces"]] == [15, 30, 34]

The locked regime counts [15, 30, 34] correspond to the Stata break indexes [15, 45]. The scatter helper returns data only; callers can render it with any plotting library or with no plotting dependency at all.

Single-Break SSR Path

payload = json.loads(
    Path("artifacts/verification/stata/phase1_us_ts_estimate_breaks1_ssr_path.json").read_text(
        encoding="utf-8"
    )
)
results = payload["results"]

estimate = xtbreak.EstimateResult.from_break_matrices(
    command=payload["command"],
    break_matrix=results["breaks"],
    ci_matrix=results["ci"],
    ssr=results["ssr"],
    ssr_path=results["ssr_path"],
    sample_time_values=results["sample_time_values"],
    metadata={"scenario_id": payload["scenario_id"]},
)

ssr_plot = estimate.ssr_path_plot_data()
assert ssr_plot["selected_break"] == {"break": 1, "index": 16, "time_value": 3142}

The selected single break is the minimum of the Stata SSRvec path and is used for estat ssr style plotting.

Verification Command

PYTHONDONTWRITEBYTECODE=1 PYTHONPATH=xtbreak-py/src \
  python -m pytest test/empirical/test_phase7_paper_output_reproduction.py \
  test/parity/test_estimate_result_example_smoke.py -q

This covers the Stata oracle files, the paper logs, the Python partial structural-change search, the EstimateResult table and plot payloads, the two-break SSR 49.0738641567431, and the single-break SSR-path selected break {index: 16, time_value: 3142}.