Espartaco

“Is that to say we are against Free Trade? No, we are for Free Trade, because by Free Trade all economical laws, with their most astounding contradictions, will act upon a larger scale, upon the territory of the whole earth; and because from the uniting of all these contradictions in a single group, where they will stand face to face, will result the struggle which will itself eventuate in the emancipation of the proletariat.”

Karl Heinrich Marx · Marx-Engels Collected Works, Vol. VI, p. 290

EnglishEspañol

Tag: detección de rupturas estructurales

  • valueprhr: When Do Market Prices Reflect the Labor That Produced Them? A Modern R Toolkit for an Old Question

    valueprhr: When Do Market Prices Reflect the Labor That Produced Them? A Modern R Toolkit for an Old Question

    You can also find this library at CRAN and download it directly from R and RStudio.

    An introduction to an R package that brings Bayesian inference, panel data econometrics, and rigorous validation to one of political economy’s most enduring empirical debates.


    Why This Package Exists

    Here is a question that has occupied economists for over two centuries: when you pay for something, does the price you pay bear any systematic relationship to the labor required to make it?

    Adam Smith thought so. David Ricardo refined the idea. Karl Marx built an entire theory of exploitation on it. And since the mid-twentieth century, empirical researchers have been trying to measure the strength of this correspondence with real-world data.

    The challenge has always been methodological. The datasets are panel data — prices observed across many economic sectors over many time periods — and they demand techniques that respect both the cross-sectional structure (different industries behave differently) and the temporal dimension (relationships can shift over time). A simple scatterplot of values against prices, however illustrative, will not settle the question.

    valueprhr is an R package built to close this methodological gap. It provides a complete, reproducible pipeline: from raw price matrices to model estimation, from Bayesian inference to out-of-sample validation, from structural break detection to side-by-side model comparison. It was designed for political economy, but as we will see, its toolkit applies to any panel data problem where you need to assess the correspondence between two variables across entities and time.


    The Core Idea (in Plain Language)

    In the classical and Marxian tradition, the value of a commodity is determined by the total labor time — direct and indirect — required to produce it. If a table requires 10 hours of socially necessary labor and a chair requires 5, the table’s value is twice the chair’s.

    This gives rise to what economists call direct prices (denoted pd): prices that are strictly proportional to the labor embodied in each commodity. They represent what prices would be if they perfectly mirrored labor content.

    But capitalism does not work that way. Capital flows between sectors seeking the highest return, and competition tends to equalize the rate of profit across industries. The prices that emerge from this process are called prices of production (denoted pπ). They redistribute surplus value: sectors with higher organic composition of capital (more machinery relative to labor) tend to have prices of production above their direct prices, and vice versa.

    The central empirical question is: despite this redistribution, how closely do direct prices and prices of production correspond?

    The standard test is a log-linear regression:

    ln(pπit) = α + β · ln(pdit) + uit

    where i indexes sectors and t indexes time periods.

    Three hypotheses are at stake:

    • β ≈ 1: a one-percent increase in direct prices is associated with roughly a one-percent increase in production prices (proportionality).
    • R² ≈ 1: direct prices explain the vast majority of the variation in production prices.
    • Stability: the relationship holds consistently across time periods.

    If all three hold, the labor theory of value has strong empirical support. valueprhr gives you the tools to test each one rigorously.


    What’s Inside the Package

    valueprhr organizes its functionality into six modules. Here is what each does and why it matters.

    1. Data Preparation

    Real-world data rarely arrives in the format econometric methods require. The package accepts two data frames in wide format (rows = years, columns = sectors) — one for direct prices, one for production prices — and converts them into the long-format panel structure that econometric models expect.

    library(valueprhr)
    # Wide format: Year | Agriculture | Manufacturing | Mining | ...
    direct <- read.csv("direct_prices.csv")
    production <- read.csv("production_prices.csv")
    # Convert to long panel: Year, Sector, direct, production, log_direct, log_production
    panel <- prepare_panel_data(direct, production, log_transform = TRUE)
    head(panel)
    #> Year Sector direct production log_direct log_production
    #> 1 1960 Agriculture 45.2 48.1 3.81 3.87
    #> 2 1961 Agriculture 46.0 49.0 3.83 3.89
    #> ...

    The function prepare_log_matrices() does the same job but returns matrix format, which is what the multivariate methods (PLS, CCA) need.

    2. Panel Data Models

    This is where the core econometrics happens. The package implements two complementary specifications:

    Two-Way Fixed Effects (FE) controls for both sector-specific and time-specific unobserved heterogeneity:

    Yit = αi + γt + β · Xit + εit

    In plain terms: every sector has its own baseline (some sectors are systematically more expensive), every year has its own macroeconomic conditions (inflation, crises), and the model isolates the within variation to estimate the core relationship.

    fe <- fit_twoway_fe(panel, robust_se = TRUE, cluster_type = "group")
    print(fe)
    #> Two-Way Fixed Effects Model
    #> ============================
    #> Observations: 1200 | Sectors: 20 | Years: 60
    #> R-squared: 0.9876 | Adjusted R-squared: 0.9870
    #>
    #> log_direct coefficient:
    #> Estimate = 0.9754, SE = 0.0123, t = 79.30, p = 0.0000

    The cluster_type = "group" option computes cluster-robust standard errors at the sector level, which accounts for serial correlation within each sector’s time series.

    Mundlak Correlated Random Effects (CRE) takes a different route. Instead of dummy variables for every sector, it decomposes the predictor into a within-sector component (how Xit deviates from sector i‘s average) and a between-sector component (the sector average itself):

    Yit = α + βW · (Xiti) + βB · i + ui + εit

    In data science language: this is a way to control for group-level confounders without the computational cost of N dummy variables. If βW = βB, the within and between effects are the same, and a simpler Random Effects model suffices. If they differ, the relationship between values and prices operates differently within a sector over time than across sectors.

    # Add Mundlak terms
    panel_cre <- create_mundlak_data(panel, x_var = "log_direct")
    # Fit the model
    cre <- fit_mundlak_cre(panel_cre, include_time_fe = TRUE)
    print(cre)
    #> Mundlak Correlated Random Effects Model
    #> =========================================
    #> Within-sector effect (beta_W): 0.9680
    #> Between-sector effect (beta_B): 0.9912
    #>
    #> Mundlak test H0: beta_W = beta_B
    #> F-stat = 2.14, p-value = 0.1438
    #> -> Fail to reject H0: RE/CRE specification is consistent

    The function test_mundlak_specification() formalizes this check. A low p-value means you should stick with Fixed Effects; a high p-value means the simpler model is adequate.

    The package also includes a Panel Granger Causality test (the Dumitrescu-Hurlin procedure), which tests whether past values of direct prices help predict current production prices — and vice versa.

    panel_granger_test(panel, lags = c(1, 2))
    #> direction lag W_stat Z_stat p_value significant
    #> 1 direct -> production 1 8.432 3.126 0.0018 TRUE
    #> 2 direct -> production 2 6.215 2.441 0.0146 TRUE
    #> 3 production -> direct 1 5.890 2.103 0.0354 TRUE
    #> 4 production -> direct 2 4.012 1.332 0.1828 FALSE

    3. Bayesian Models

    Classical (frequentist) estimation gives you a single point estimate for β. Bayesian methods give you a full probability distribution over possible values, incorporating your prior beliefs and updating them with the data.

    In econometric language: instead of β̂ = 0.975 ± 0.012, you get a posterior distribution showing that β lies between 0.95 and 1.00 with 95% probability.

    The package offers two Bayesian approaches:

    Sector-by-Sector Bayesian GLM fits an independent Bayesian linear model for each sector, using weakly informative priors (the rstanarm package handles the MCMC sampling via Stan). Each sector gets its own slope and intercept, along with Leave-One-Out Cross-Validation (LOO-CV) scores.

    bayes <- fit_bayesian_glm_sectors(
    direct, production,
    chains = 4, iter = 4000
    )
    print(bayes$summary_table)
    #> Sector beta_mean beta_sd beta_lower beta_upper elpd looic n_obs
    #> 1 Agriculture 0.982 0.025 0.933 1.031 -42.3 84.6 60
    #> 2 Manufacturing 0.971 0.031 0.910 1.031 -38.7 77.4 60
    #> 3 Mining 0.958 0.042 0.876 1.041 -45.1 90.2 60
    #> ...

    In data science language: LOO-CV is a principled way to assess out-of-sample predictive performance without holding out data. The LOOIC (LOO Information Criterion) is the Bayesian analogue of AIC — lower is better.

    Bayesian Hierarchical Model goes further by pooling information across sectors. Instead of treating each sector in isolation, it assumes that sector-specific slopes are drawn from a common population distribution:

    βi ~ N(μβ, σβ2)

    Sectors with less data “borrow strength” from the population mean. This is especially valuable when some sectors have short time series.

    hier <- fit_bayesian_hierarchical(panel, include_time = TRUE)
    print(hier)
    #> Bayesian Hierarchical Model
    #> ============================
    #> Observations: 1200 | Sectors: 20
    #>
    #> LOO-CV:
    #> ELPD = -312.45
    #> LOOIC = 624.90
    #>
    #> Population-level effects:
    #> parameter mean sd 2.5% 97.5%
    #> 1 (Intercept) 0.1423 0.0892 -0.032 0.317
    #> 2 log_direct 0.9734 0.0145 0.945 1.002
    #> 3 Time_scaled 0.0031 0.0018 -0.0004 0.007

    4. Multivariate Analysis

    When the number of sectors (N) is large relative to the number of time periods (T), standard regression becomes unstable. This is the “small T, large N” problem common in panel data. The package offers three multivariate techniques to handle it:

    Partial Least Squares (PLS) extracts latent components that explain covariance between direct prices and production prices. It handles multicollinearity gracefully and is widely used in chemometrics, genomics, and now in value-price analysis.

    matrices <- prepare_log_matrices(direct, production)
    pls <- fit_pls_multivariate(
    matrices$X_clean, matrices$Y_clean,
    max_components = 8
    )
    print(pls)
    #> Partial Least Squares (PLS) Regression
    #> =======================================
    #> Optimal components: 3
    #>
    #> R-squared by component:
    #> n_components R2_train R2_cv
    #> 1 1 0.942 0.938
    #> 2 2 0.971 0.965
    #> 3 3 0.984 0.980

    Canonical Correlation Analysis (CCA) finds linear combinations of direct prices and production prices that are maximally correlated. In econometric language: CCA extracts the “shared economic signal” — the common factor driving both sets of prices.

    cca <- run_sparse_cca(matrices$X_clean, matrices$Y_clean, n_components = 3)
    print(cca)
    #> Canonical Correlation Analysis
    #> ===============================
    #> Components: 3
    #>
    #> Canonical correlations:
    #> CC1: r = 0.9987 (Var X: 92.3%, Var Y: 91.8%)
    #> CC2: r = 0.9841 (Var X: 5.1%, Var Y: 5.4%)
    #> CC3: r = 0.9523 (Var X: 1.8%, Var Y: 1.9%)

    The first canonical correlation above 0.99 indicates an extremely tight structural link between the two price systems.

    Panel VAR captures dynamic feedback: do lagged values of direct prices predict current production prices, and vice versa?

    pvar <- fit_panel_var(panel, lags = 2, transformation = "fd")

    5. Cross-Validation

    Standard k-fold cross-validation violates temporal ordering. If you train on 1960–1990 and test on 1985–1990, future information leaks into the training set. The package implements two time-aware approaches:

    Rolling Window CV trains on t₀ … tW, tests on tW+1tW+H, then rolls the window forward.

    cv <- rolling_window_cv(
    panel,
    window_sizes = c(20, 30),
    step_size = 2,
    test_horizon = 3
    )
    print(cv$summary)

    Leave-One-Sector-Out (LOSO) trains on all sectors except one and predicts the held-out sector. This tests cross-sectional generalization: does the value-price relationship estimated from other sectors hold for agriculture? For mining? For finance?

    loso <- leave_one_sector_out(panel)
    print(loso$summary)
    #> metric mean sd
    #> 1 RMSE 0.04521 0.01832
    #> 2 MAE 0.03587 0.01456
    #> 3 R_squared 0.96120 0.02340

    An average R² above 0.96 in LOSO-CV means the relationship generalizes robustly across sectors.

    6. Structural Break Tests

    Has the value-price relationship been stable over time? Or did it shift at some point — due to globalization, a methodological change in data construction, a technological revolution, or a regime shift in profit rate equalization?

    The package aggregates the panel to a time series and applies a battery of tests:

    breaks <- test_structural_breaks(panel, break_date = 1990)
    print(breaks)
    #> Structural Break Tests
    #> ========================
    #> Time-series observations: 60
    #>
    #> Chow Test:
    #> Break date: 1990
    #> F-stat = 1.8420, p = 0.1687
    #>
    #> supF / Bai-Perron Test:
    #> supF = 5.2130, p = 0.0842
    #> Breaks detected: 0

    A non-significant result is actually good news here: it means the value-price correspondence has been structurally stable across the entire sample period.


    The Full Pipeline in One Command

    If you want to run everything at once — data preparation, FE and CRE models, cross-validation, structural break tests, and model comparison — the package offers a single entry point:

    results <- run_full_analysis(
    direct,
    production,
    run_bayesian = FALSE, # Set TRUE if you have rstanarm installed
    run_cv = TRUE,
    run_breaks = TRUE,
    verbose = TRUE
    )
    # Access everything
    print(results$comparison)
    print(results$cv_summary)
    cat(format_break_results(results$breaks))

    You can then export the comparison table and CV results to CSV:

    export_results_csv(
    results$comparison,
    results$cv_summary,
    output_dir = "results/"
    )

    Beyond Political Economy: General Panel Data Applications

    Although valueprhr was built for the specific question of value-price correspondence, its methods are general-purpose panel data tools. Any research problem involving the relationship between two variables observed across entities and time can benefit from the package:

    • Health economics: Does out-of-pocket spending track underlying treatment costs across regions over time?
    • Environmental economics: Do carbon prices reflect the embodied emissions of goods across industries?
    • Education: Do standardized test scores correspond to instructional expenditure across school districts over decades?
    • Finance: Do book values predict market valuations across sectors?
    • Any two-variable panel regression where you need fixed effects, Mundlak decomposition, robust standard errors, time-aware cross-validation, or structural break detection.

    The key requirement is that your data has a panel structure (entities × time) and that the distributional assumptions of the models are reasonable for your context. The methods — two-way FE, Mundlak CRE, Bayesian hierarchical models, PLS, CCA, rolling-window CV, structural break tests — are econometric staples that transcend any particular application domain.


    Installation and Dependencies

    The package requires R ≥ 4.1.0. Core functionality depends only on base R and the Metrics package. Extended features (Bayesian models, panel data infrastructure, structural break tests) are handled through soft dependencies that are loaded on demand:

    # Install from GitHub
    install.packages("devtools")
    devtools::install_github("isadorenabi/valueprhr")
    # Optional: install all suggested packages at once
    suggested <- c(
    "rstanarm", "loo", "plm", "lme4", "pls", "vars",
    "panelvar", "strucchange", "lmtest", "sandwich",
    "dplyr", "tidyr", "tibble"
    )
    install.packages(suggested[!sapply(suggested, requireNamespace, quietly = TRUE)])

    Note for Bayesian models: rstanarm requires a working C++ toolchain — Rtools on Windows, Xcode Command Line Tools on macOS, or build-essential on Linux.


    What Makes This Package Methodologically Different

    Three features distinguish valueprhr from a hand-rolled analysis:

    1. Time-aware validation. Most applied work reports in-sample R² as evidence of fit. valueprhr pairs every model with rolling-window and leave-one-sector-out cross-validation, giving you out-of-sample performance that is honest about temporal dependence and cross-sectional generalization.
    2. The Mundlak decomposition. By splitting effects into within-sector and between-sector components, the package lets you test whether the value-price relationship operates at the sector level (structural), at the temporal level (cyclical), or both. This is a nuance that most empirical studies in this literature overlook.
    3. Bayesian hierarchical pooling. Sectors with short time series are a common headache. The hierarchical model lets small sectors borrow statistical strength from the population, producing more stable estimates than independent sector-by-sector regressions.

    A Note on the Underlying Data

    The wiki documentation mentions that market price indices used in this framework are constructed by temporally disaggregating the aggregate Consumer Price Index using the Input-Output matrix as a structural indicator, relying on closed-form Bayesian solutions from the BayesianDisaggregation library. This is a methodological detail worth understanding: the sectoral prices are not raw market quotes but statistically consistent decompositions of the macroeconomic aggregate. This ensures that the estimated sectoral price movements add up to the observed CPI, a property that many ad hoc sectoral price datasets lack.


    Citation

    If you use valueprhr in your research:

    @software{gomezjulian2025valueprhr,
    author = {Gómez Julián, José Mauricio},
    title = {valueprhr: Value-Price Analysis with Bayesian and Panel Data Methods},
    year = {2025},
    url = {https://github.com/isadorenabi/valueprhr},
    note = {R package version 0.1.0}
    }

    Author: José Mauricio Gómez Julián — ORCID — isadore.nabi@pm.me

    License: MIT

    Repository: github.com/IsadoreNabi/valueprhr


    The labor theory of value is either one of the most important ideas in the history of economics or one of the most contested. Either way, it deserves better tools than a spreadsheet and a prayer. valueprhr brings the full machinery of modern econometrics to the question — and lets the data speak for itself.

  • Detecting When the World Changes: A Deep Look at the RegimeChange R Package

    Detecting When the World Changes: A Deep Look at the RegimeChange R Package

    You can also find this library at CRAN and download it directly from R and RStudio.

    Imagine you’re monitoring a manufacturing line, and the sensor readings suddenly look… different. The numbers aren’t obviously wrong — they’re within range — but something has shifted. Or picture an economist staring at decades of GDP data, trying to pinpoint exactly when a country’s growth model fundamentally changed. Or a doctor watching a patient’s vitals, waiting for the moment “stable” tips into “critical.”

    All of these scenarios share the same underlying question: when did the system switch from one regime to another?

    This is the problem of regime change detection (also called changepoint detection), and it’s one of the most practically important — and theoretically rich — problems in time series analysis. A new open-source R package called RegimeChange tackles this problem with unusual breadth, combining classical statistics, Bayesian inference, and even deep learning under a single, coherent interface. In this post, I’ll walk you through what it does, how it works, and why it matters — in plain language, without cutting corners on the technical details.


    The Problem: Signal vs. Noise

    Every observable system fluctuates. The central question of changepoint detection is deceptively simple: is this fluctuation just noise (random variation within the same regime), or is it signal (evidence that the system has transitioned to a qualitatively different state)?

    This is harder than it sounds. There’s an irreducible tension between two types of error:

    • False alarm (Type I): You detect a change where none exists. In manufacturing, this means shutting down a line for no reason. In medicine, it’s a false positive that triggers unnecessary intervention.
    • Missed detection (Type II): You fail to detect a real change. In finance, this means holding a position through a regime shift. In epidemiology, it means missing the start of an outbreak.

    Every detection algorithm negotiates this trade-off through its threshold: how much evidence do you demand before declaring that the world has changed? Set it too low and you cry wolf; set it too high and you react too late. There is no universal answer — the right threshold depends on the relative cost of each type of error in your specific context.


    What Makes RegimeChange Different

    Several R packages already address changepoint detection — the venerable changepoint package, wbs, not, ecp, and the Python ruptures library, to name a few. Each has strengths and limitations. RegimeChange distinguishes itself by integrating approaches that are usually siloed:

    1. Frequentist and Bayesian methods together. Most packages pick one paradigm. RegimeChange gives you both, plus deep learning, all callable through the same detect_regimes() function.
    2. Offline and online modes. You can analyze a complete dataset retrospectively (“when did the changes occur?”) or process streaming data in real time (“is a change happening right now?”) using the same library.
    3. Native uncertainty quantification. Every changepoint estimate comes with confidence intervals (via bootstrap for frequentist methods) or full posterior distributions (for Bayesian methods). This isn’t just a point estimate — it’s a statement about how confident you should be.
    4. Robustness to messy data. Real-world data has outliers, heavy tails, and autocorrelation. RegimeChange includes configurable robust estimation and AR(1) pre-whitening — features that, as the benchmarks show, make a dramatic difference on contaminated data.
    5. An optional Julia backend. For large datasets, the package can transparently dispatch to high-performance Julia implementations while keeping the R interface you’re used to.

    Let’s look at each of these in more detail.


    The Three Method Families

    RegimeChange organizes its detection algorithms into three families. Here’s what each one does, explained simply.

    Frequentist Methods

    These are the classical, well-established workhorses of changepoint detection.

    CUSUM (Cumulative Sum) is the granddaddy of them all, dating to E. S. Page’s 1954 paper. The idea is beautifully intuitive: you accumulate deviations from a target (or “baseline”) value over time. If the system is in the same regime, positive and negative deviations cancel out, and the cumulative sum hovers near zero. When a change occurs, deviations start piling up in one direction, and the statistic grows. When it crosses a threshold, you sound the alarm. CUSUM runs in O(n) time and is ideal for detecting a single change in mean, especially in online monitoring scenarios.

    PELT (Pruned Exact Linear Time), introduced by Killick, Fearnhead, and Eckley in 2012, is the modern gold standard for detecting multiple changepoints. It uses dynamic programming to find the globally optimal segmentation of your data — the set of changepoints that minimizes a cost function (essentially: how well do the segments fit the data, penalized by the number of segments). The “pruning” part is clever: it discards candidate changepoint locations that provably can’t be optimal, which brings the average-case complexity down from O(n²) to O(n). RegimeChange supports several penalty criteria — BIC, AIC, MBIC, MDL — or you can specify a manual numeric penalty.

    Binary Segmentation is a simpler, greedy alternative: find the single best changepoint, split the data there, then recurse on each segment. It’s fast (O(n log n)) but doesn’t guarantee the global optimum. Wild Binary Segmentation (WBS), from Fryzlewicz (2014), improves on this by searching over many random sub-intervals rather than the full dataset, making it much more robust when changepoints are close together.

    The package also includes FPOP (Functional Pruning Optimal Partitioning), which maintains piecewise-quadratic cost functions for even better theoretical guarantees; E-Divisive, a nonparametric method using energy statistics that can detect changes in any aspect of the distribution; Kernel-based CPD, which maps data into a reproducing kernel Hilbert space to capture complex distributional shifts; and NOT (Narrowest-Over-Threshold), which excels at precise localization.

    Bayesian Methods

    Where frequentist methods give you a point estimate and (if you ask) a confidence interval, Bayesian methods give you something richer: a full probability distribution over when the change might have occurred.

    BOCPD (Bayesian Online Changepoint Detection), from Adams and MacKay (2007), is the flagship. At each time step, it maintains a posterior distribution over the “run length” — the number of observations since the last changepoint. When a new data point arrives, it updates this distribution: either the run length grows by one (no change), or it resets to zero (change detected). The probability of a reset at each step is governed by a “hazard function,” typically a geometric distribution representing a constant probability of change per time step.

    What makes BOCPD powerful is that it naturally produces, at every time point, the probability that a changepoint just occurred. This is a fundamentally different kind of output than a binary “change/no-change” flag — it’s a calibrated measure of confidence that updates with every new observation. RegimeChange implements BOCPD with several conjugate prior models: Normal-Gamma (for unknown mean and variance), Normal with known variance, Gamma-Poisson (for count data), and Normal-Wishart (for multivariate data).

    Shiryaev-Roberts is another Bayesian sequential method, based on work by Shiryaev from the 1960s. It accumulates likelihood ratios from all possible past changepoint times and is asymptotically optimal for minimizing detection delay — that is, for detecting changes as quickly as possible once they occur.

    Deep Learning Methods

    For complex, nonlinear patterns that classical methods struggle with, RegimeChange offers optional deep learning detectors (requiring keras and tensorflow):

    • Autoencoders are trained to reconstruct “normal” patterns. When the data changes regime, reconstruction error spikes — the model can’t rebuild what it hasn’t seen before.
    • Temporal Convolutional Networks (TCNs) use dilated causal convolutions to model long-range temporal dependencies, operating either in supervised mode (if you have labeled changepoints for training) or unsupervised mode (predicting the next value and flagging prediction errors).
    • Transformers bring self-attention to the problem, capturing long-range relationships in the time series.
    • Contrastive Predictive Coding (CPC) learns representations through self-supervised contrastive learning, detecting changes where learned embeddings shift significantly.

    An ensemble mode combines multiple deep learning methods, requiring a minimum number of them to agree before declaring a changepoint — a strategy that trades some sensitivity for robustness.


    One Interface to Rule Them All

    The genius of RegimeChange is that you don’t need to learn a different API for each method. Everything flows through a single function:

    library(RegimeChange)
    # Generate data with a changepoint at t=200
    set.seed(42)
    data <- c(rnorm(200, 0, 1), rnorm(200, 2, 1))
    # Detect using PELT (the default offline method)
    result <- detect_regimes(data, method = "pelt")
    print(result)
    plot(result)
    # Switch to Bayesian online detection
    result <- detect_regimes(data, method = "bocpd",
    prior = normal_gamma(mu0 = 0, kappa0 = 1,
    alpha0 = 1, beta0 = 1))
    # Use CUSUM for online monitoring
    result <- detect_regimes(data, method = "cusum",
    mode = "online", threshold = 5)

    You specify what kind of change you’re looking for (type = "mean", "variance", "both", "trend", or "distribution"), how many changes you expect ("single", "multiple", or a specific integer), and the penalty criterion for model complexity. The function returns a regime_result object containing the detected changepoints, their confidence intervals, segment-by-segment statistics (mean, variance, etc.), and — for Bayesian methods — the full posterior distribution.

    Online Detection

    For real-time monitoring, you create a persistent detector object and feed it data one observation at a time:

    detector <- regime_detector(method = "bocpd",
    prior = normal_gamma(),
    threshold = 0.5)
    for (x in data_stream) {
    detector <- update(detector, x)
    if (detector$last_result$alarm) {
    message("Changepoint detected at time ", detector$last_result$t)
    detector <- reset(detector)
    }
    }

    This is the pattern you’d use for industrial monitoring, fraud detection, or epidemiological surveillance — anywhere data arrives sequentially and you need to react in real time.


    Robustness: Where RegimeChange Really Shines

    Real data is messy. It has outliers, heavy tails, and autocorrelation. Standard changepoint methods, which typically assume independent, normally distributed observations, can fail dramatically under these conditions.

    RegimeChange addresses this with a configurable robustness layer in its PELT implementation. When you set robust = TRUE, the package:

    • Winsorizes the data (clips extreme values beyond a specified quantile) to limit the influence of outliers.
    • Uses Huber M-estimation or Tukey biweight loss functions instead of standard squared-error loss. These functions grow linearly (Huber) or flatten out entirely (Tukey) for large residuals, so a single outlier can’t dominate the cost calculation.
    • Employs the Qn scale estimator (from Rousseeuw and Croux, 1993), which is far more efficient than the median absolute deviation (MAD) — 82% efficiency vs. 37% — while maintaining the same 50% breakdown point (meaning up to half the data can be contaminated before the estimator breaks).

    You can also set robust = "auto", which analyzes the data’s contamination level — comparing the MAD to the standard deviation and counting outlier proportions — and automatically selects the appropriate robustness level ("none", "mild", "moderate", or "aggressive").

    For autocorrelated data, the correct_ar = TRUE option estimates the AR(1) coefficient and applies a pre-whitening transformation (subtracting the lagged value scaled by the autocorrelation coefficient) before running detection. This prevents autocorrelation from inflating false positive rates.

    The Benchmark Numbers

    The package’s wiki reports extensive benchmarks against established R packages (changepoint, wbs, not, ecp) across 17 scenarios with 50 replications each. The results are striking:

    • Overall: RegimeChange’s automatic selector achieved the highest mean F1 score (0.804) across all methods tested.
    • Contaminated data: The robust mode achieved F1 = 0.998 across all contamination scenarios, versus 0.479 for the standard changepoint PELT — a 108% improvement on messy data.
    • Subtle variance changes (2:1 ratio): RegimeChange scored 0.667 vs. 0.400 for reference packages — a 67% improvement.
    • Strong autocorrelation (AR coefficient 0.7): RegimeChange scored 0.900 vs. 0.720 — a 25% improvement.

    In 16 of 17 scenarios, RegimeChange either won or tied. The single loss was marginal (0.017 in one heavy-tailed scenario).


    Uncertainty Quantification: Beyond Point Estimates

    A changepoint location of “observation 200” is only half the story. The other half is: how confident are we? RegimeChange takes this seriously, providing two types of uncertainty:

    Location uncertainty — how precise is the estimate? For frequentist methods, the package uses block bootstrap resampling: it repeatedly resamples blocks of the data (preserving local dependence structure), reruns detection, and computes 95% confidence intervals from the distribution of bootstrap estimates. For Bayesian methods, uncertainty comes naturally from the posterior distribution over run lengths.

    Existence uncertainty — is there even a change at all? BOCPD produces, at each time point, the posterior probability that a changepoint just occurred. This is arguably more informative than a p-value: it’s a direct probabilistic statement that you can use for decision-making.


    The Julia Backend

    R is wonderful for data analysis, but it can be slow for computationally intensive tasks. RegimeChange ships with an optional Julia backend — a complete reimplementation of the core algorithms (PELT, FPOP, BOCPD, CUSUM, Kernel CPD, WBS, and multivariate PELT) in Julia, a language designed for numerical performance.

    The integration is seamless. You initialize Julia once:

    init_julia()
    julia_available() # Returns TRUE if ready

    After that, the package automatically dispatches to Julia for large datasets (n > 1000 by default) while using R for smaller ones. You can benchmark the difference with benchmark_backends(). The Julia implementation includes thoughtful numerical engineering: Welford’s algorithm for numerically stable running mean/variance, Kahan compensated summation for cumulative sums, log-domain computations to avoid underflow in BOCPD, and careful handling of ill-conditioned covariance matrices in the multivariate case.

    If Julia isn’t installed, the package simply falls back to the R implementation — no errors, no fuss.


    Evaluation and Comparison

    RegimeChange includes a comprehensive evaluation framework. If you know the true changepoints (e.g., in simulated data or a well-studied real dataset), you can compute a full battery of metrics:

    metrics <- evaluate(result, true_changepoints = c(100, 250), tolerance = 10)
    print(metrics)

    This gives you precision, recall, and F1 score (with a tolerance window for matching), Hausdorff distance (the worst-case localization error), Rand Index and Adjusted Rand Index (segmentation agreement), and the covering metric (a weighted IoU measure). You can also run compare_methods() to test multiple algorithms on the same data and see the results side by side.

    The package ships with several built-in datasets: well_log (geophysical measurements), industrial_sensor data, economic_cycles, and simulated_changepoints — giving you realistic examples to experiment with.


    When to Use What

    Here’s a practical guide based on the library’s design and benchmark results:

    Your SituationRecommendation
    Clean data, clear mean shiftsPELT with default settings
    Outliers or heavy tailsPELT with robust = TRUE (or "auto")
    Autocorrelated time seriesPELT with correct_ar = TRUE
    Need real-time detectionBOCPD or CUSUM via regime_detector()
    Need probability of change (not just yes/no)BOCPD
    Nonparametric, no distributional assumptionsE-Divisive or Kernel CPD
    Complex nonlinear patterns (and you have data to train on)Deep learning methods
    High-dimensional multivariate dataSparse Projection CPD
    Not sure which method to useEnsemble mode, or let method = "auto" choose

    The Bigger Picture

    What I find compelling about RegimeChange is its philosophical stance. The wiki frames changepoint detection not just as a technical problem but as an epistemological one: how do we distinguish real change from noise? The likelihood ratio — the ratio of how probable an observation is under the “changed” hypothesis vs. the “unchanged” hypothesis — is literally a quantification of surprise. Each observation casts a vote, evidence accumulates, and the decision threshold represents how much surprise we demand before acting.

    This framing connects the technical machinery to the real-world stakes. Whether you’re monitoring a factory, tracking a disease outbreak, or studying climate data, the fundamental structure is the same: the world was in one state, it transitioned to another, and you need to detect that transition as reliably and quickly as possible.

    RegimeChange gives you a unified toolkit for that task — one that spans the classical-to-modern spectrum, handles messy data gracefully, tells you not just where the change is but how confident you should be, and scales from your laptop to a Julia-accelerated backend when you need speed.


    Getting Started

    RegimeChange is available on GitHub and installs like any R package from source:

    # install.packages("devtools")
    devtools::install_github("IsadoreNabi/RegimeChange")

    It requires R ≥ 4.0.0 and depends on ggplot2, rlang, cli, and magrittr. The Julia backend is optional (requires Julia ≥ 1.6 and the JuliaCall R package), as are the deep learning methods (require keras and tensorflow). The package ships with three vignettes — an introduction, an offline detection guide, and a Bayesian methods tutorial — and the wiki is a comprehensive reference covering everything from mathematical foundations to detailed benchmarks.

    For researchers and practitioners who need rigorous, flexible, and robust changepoint detection in R, RegimeChange is well worth a close look.