Binary Autoregressive Time Series

• Binary autoregressive time series are binary-valued stochastic processes whose current outcome depends on a finite history of past binary observations, permitting a finite-order Markov chain representation.
• The framework encompasses a range of models including logistic autoregressions, Bernoulli vector processes, de Bruijn models, dynamic multivariate probit, and generalized binary VAR, each with distinct parameterizations and inferential implications.
• These models are applied in diverse fields such as epidemiology, finance, and network analysis, offering insights into structural learning, mixing properties, and high-dimensional inference challenges.

Binary autoregressive time series are binary-valued stochastic processes in which the distribution of the current observation depends on lagged binary outcomes, and in some formulations also on exogenous covariates or contemporaneous latent shocks. In the literature, this umbrella includes at least five non-equivalent constructions: logistic autoregressions for univariate binary data, multivariate Bernoulli autoregressive processes on directed graphs, higher-order history-indexed Markov models based on de Bruijn graphs, dynamic multivariate probit-type threshold models, and generalized binary vector autoregressions built from random copy/flip/innovate mechanisms. What unifies them is that dependence is mediated by a finite lag history, so the process can be represented as a finite-order Markov chain after state augmentation, even though the inferential, computational, and interpretive consequences differ substantially across model classes (Gao et al., 2017, Katselis et al., 2016, Kimpton et al., 2022, Franchi et al., 2024, Dai et al., 29 Nov 2025).

1. Formal scope and principal model classes

Binary autoregressive modeling is not a single specification but a family of finite-memory dynamic models for data in $\{0,1\}$. In a univariate logistic autoregression of order $p$, the conditional success probability is

$$P(Y_t=1\mid \mathcal F_{t-1})=\frac{\exp(\beta_0+\sum_{j=1}^p \beta_j Y_{t-j})}{1+\exp(\beta_0+\sum_{j=1}^p \beta_j Y_{t-j})},$$

with inference conducted conditionally on the initial observations. In a multivariate Bernoulli autoregressive process, the next binary vector is generated coordinatewise from Bernoulli parameters formed by a sparse linear predictor and an innovation term,

$$X_{k+1}=\mathrm{Ber}\!\left(A f(X_k)+B W_{k+1}\right),$$

where the sign of influence is encoded by whether a parent contributes $X_j$ or $1-X_j$ inside $f_i$. In a binary de Bruijn process of order $m$, each lag history in $\{0,1\}^m$ has its own transition probability, so the model is equivalent to an $m$th-order binary Markov chain on letters or a first-order Markov chain on words. In the dynamic multivariate probit formulation, a latent Gaussian vector is thresholded,

$p$0

and in the generalized binary VAR, binary support is preserved through random selection matrices and Bernoulli innovations rather than a GLM link (Gao et al., 2017, Katselis et al., 2016, Kimpton et al., 2022, Franchi et al., 2024, Dai et al., 29 Nov 2025).

The main classes can be organized as follows.

Class	Representative form	Distinctive feature
LAR/LARX	$p$1	finite-sample Ex-FI
BAR	$p$2	sparse directed signed graph
DBP	$p$3	one parameter per history
Dynamic multivariate probit	$p$4	latent Gaussian dependence
gbVAR	$p$5	copy / flip / innovate

This taxonomy shows that “autoregressive” refers to dependence on lagged binary outcomes, not to a unique parameterization. A plausible implication is that model choice is driven less by the binary codomain itself than by whether the objective is conditional-probability interpretation, graph recovery, run-length analysis, composite-likelihood inference, or high-dimensional simultaneous inference.

2. Parameterization and dependence mechanisms

The simplest parametric mechanism is the logistic autoregression, where lagged binary responses enter additively on the logit scale. The linear predictor

$p$6

induces a conditional Bernoulli variance $p$7, and the first $p$8 observations are treated as initial conditions. The LARX($p$9) extension adds exogenous regressors $$P(Y_t=1\mid \mathcal F_{t-1})=\frac{\exp(\beta_0+\sum_{j=1}^p \beta_j Y_{t-j})}{1+\exp(\beta_0+\sum_{j=1}^p \beta_j Y_{t-j})},$$0 through $$P(Y_t=1\mid \mathcal F_{t-1})=\frac{\exp(\beta_0+\sum_{j=1}^p \beta_j Y_{t-j})}{1+\exp(\beta_0+\sum_{j=1}^p \beta_j Y_{t-j})},$$1 (Gao et al., 2017).

The BAR process replaces the logit link by a direct Bernoulli-parameter linear form. For node $$P(Y_t=1\mid \mathcal F_{t-1})=\frac{\exp(\beta_0+\sum_{j=1}^p \beta_j Y_{t-j})}{1+\exp(\beta_0+\sum_{j=1}^p \beta_j Y_{t-j})},$$2,

$$P(Y_t=1\mid \mathcal F_{t-1})=\frac{\exp(\beta_0+\sum_{j=1}^p \beta_j Y_{t-j})}{1+\exp(\beta_0+\sum_{j=1}^p \beta_j Y_{t-j})},$$3

Here $$P(Y_t=1\mid \mathcal F_{t-1})=\frac{\exp(\beta_0+\sum_{j=1}^p \beta_j Y_{t-j})}{1+\exp(\beta_0+\sum_{j=1}^p \beta_j Y_{t-j})},$$4 indicates an edge $$P(Y_t=1\mid \mathcal F_{t-1})=\frac{\exp(\beta_0+\sum_{j=1}^p \beta_j Y_{t-j})}{1+\exp(\beta_0+\sum_{j=1}^p \beta_j Y_{t-j})},$$5 in the directed graph, but positive and negative effects are not represented by positive and negative coefficients. Instead, all nonzero $$P(Y_t=1\mid \mathcal F_{t-1})=\frac{\exp(\beta_0+\sum_{j=1}^p \beta_j Y_{t-j})}{1+\exp(\beta_0+\sum_{j=1}^p \beta_j Y_{t-j})},$$6 are nonnegative, and polarity is carried by whether $$P(Y_t=1\mid \mathcal F_{t-1})=\frac{\exp(\beta_0+\sum_{j=1}^p \beta_j Y_{t-j})}{1+\exp(\beta_0+\sum_{j=1}^p \beta_j Y_{t-j})},$$7 or $$P(Y_t=1\mid \mathcal F_{t-1})=\frac{\exp(\beta_0+\sum_{j=1}^p \beta_j Y_{t-j})}{1+\exp(\beta_0+\sum_{j=1}^p \beta_j Y_{t-j})},$$8, meaning that either $$P(Y_t=1\mid \mathcal F_{t-1})=\frac{\exp(\beta_0+\sum_{j=1}^p \beta_j Y_{t-j})}{1+\exp(\beta_0+\sum_{j=1}^p \beta_j Y_{t-j})},$$9 or $$X_{k+1}=\mathrm{Ber}\!\left(A f(X_k)+B W_{k+1}\right),$$0 enters the predictor. Row normalization,

$$X_{k+1}=\mathrm{Ber}\!\left(A f(X_k)+B W_{k+1}\right),$$1

keeps all Bernoulli parameters in $$X_{k+1}=\mathrm{Ber}\!\left(A f(X_k)+B W_{k+1}\right),$$2 and yields a sparse description of a $$X_{k+1}=\mathrm{Ber}\!\left(A f(X_k)+B W_{k+1}\right),$$3 transition matrix by at most $$X_{k+1}=\mathrm{Ber}\!\left(A f(X_k)+B W_{k+1}\right),$$4 effective parameters when the maximum in-degree is $$X_{k+1}=\mathrm{Ber}\!\left(A f(X_k)+B W_{k+1}\right),$$5 (Katselis et al., 2016).

The binary de Bruijn process uses a different mechanism: it assigns a separate success probability to each binary lag history of length $$X_{k+1}=\mathrm{Ber}\!\left(A f(X_k)+B W_{k+1}\right),$$6. With state space $$X_{k+1}=\mathrm{Ber}\!\left(A f(X_k)+B W_{k+1}\right),$$7, each state has exactly two outgoing transitions, formed by deleting the oldest bit and appending either $$X_{k+1}=\mathrm{Ber}\!\left(A f(X_k)+B W_{k+1}\right),$$8 or $$X_{k+1}=\mathrm{Ber}\!\left(A f(X_k)+B W_{k+1}\right),$$9. This is equivalent to specifying an arbitrary map from the last $X_j$0 outcomes to the next-step success probability. Relative to logistic or probit autoregression, it automatically captures interactions among lags because each history pattern has its own parameter; the cost is $X_j$1 free transition probabilities (Kimpton et al., 2022).

The dynamic multivariate probit model places the dependence structure in a latent Gaussian layer. Lagged binary vectors enter through coefficient matrices $X_j$2, exogenous covariates through $X_j$3, and contemporaneous cross-sectional dependence through the covariance matrix $X_j$4 of $X_j$5. This separates dynamic cross-lag effects from instantaneous latent correlation. The generalized binary VAR again differs: it does not specify $X_j$6 via a link function, but builds $X_j$7 from row-wise multinomial selection among lagged coordinates or an innovation component, so each coordinate is obtained by copying a lagged bit, flipping it when the coefficient is negative, or replacing it by a Bernoulli innovation (Franchi et al., 2024, Dai et al., 29 Nov 2025).

3. Markov structure, stationarity, and mixing behavior

All of these models admit a finite-state interpretation after suitable state augmentation. In the de Bruijn process, the word process $X_j$8 is first-order Markov on $X_j$9 states. In LAR($1-X_j$0), the lag vector $1-X_j$1 plays the role of the Markov state for the conditional likelihood. In multivariate threshold models, stacking

$1-X_j$2

produces a first-order representation in the augmented state (Kimpton et al., 2022, Gao et al., 2017, Franchi et al., 2024).

The BAR process yields one of the sharpest structural results. Because all $1-X_j$3, the chain on $1-X_j$4 is irreducible and aperiodic; since the state space is finite, it has a unique stationary distribution and is geometrically ergodic. Defining

$1-X_j$5

the paper proves

$1-X_j$6

hence $1-X_j$7. This logarithmic mixing in the dimension is derived by a synchronous coupling argument with shared innovations and shared uniforms, under which coordinate disagreements contract roughly by powers of $1-X_j$8. By contrast, single-site BAR random walks on the hypercube have $1-X_j$9 mixing, consistent with one-site-at-a-time updates (Katselis et al., 2016).

The 2024 multivariate probit theory establishes a different stationarity principle. For threshold-type autoregressive binary models with stationary exogenous covariates, existence of a stationary path is shown to be “almost automatic” when the innovation law is not compactly supported. Because the lagged endogenous state space is finite, sufficiently rich innovation support implies that threshold crossings in either direction remain possible from any state, so contraction restrictions of the sort common in continuous-state nonlinear autoregressions are not required. Under stronger irreducibility-type conditions, the stationary solution is unique in distribution and ergodic (Franchi et al., 2024).

The de Bruijn framework emphasizes stationary distributions over words. For irreducible, persistent, aperiodic chains,

$f_i$0

and for $f_i$1 the paper gives an explicit formula for $f_i$2. This makes clear that marginal success probability and serial dependence are distinct: the same one-dimensional marginal can coexist with markedly different clustering and run-length behavior (Kimpton et al., 2022).

4. Estimation and finite-sample inference

Inference in binary autoregressive models is complicated by the fact that lagged binary responses are endogenous covariates. In the logistic autoregressive setting, the exact conditional Fisher information is

$f_i$3

which averages the negative Hessian over the model-implied distribution of lag histories. The paper derives this matrix for finite $f_i$4, provides a recursive state-space computation through lag-state probabilities $f_i$5, and shows that Ex-FI approaches the asymptotic Fisher information as $f_i$6. In simulations, confidence intervals based on Ex-FI tend to be narrower than those based on the empirical Fisher information while maintaining type I error rates at or below nominal levels; in a respiratory-rate application, this narrower interval estimation changes whether some odds-based intervals exclude $f_i$7 (Gao et al., 2017).

The BAR paper addresses a different inferential target: structural learning of the directed signed graph from a single trajectory. Its algorithm, BARObs, has two stages. First, “supergraph selection” computes pairwise conditional influences

$f_i$8

and ranks $f_i$9 to construct an estimated super-neighborhood of size $m$0. Second, “supergraph trimming” evaluates all $m$1 conditional probabilities on the selected neighborhood, uses near-maximizers within $m$2 of the empirical maximum, and retains nodes that are fixed across these near-maximizers. Under a BAR identifiability condition controlling same-time correlations, exact graph recovery occurs with probability at least $m$3 when

$m$4

The computational complexity is $m$5 (Katselis et al., 2016).

In the de Bruijn process, likelihood factorization is simpler because conditioning is purely by lag history. With counts $m$6 of transitions between compatible words, the likelihood is a product of state-specific Bernoulli/binomial terms, so the MLE is an empirical transition frequency for each history. The paper also gives a Fisher-information expression, conjugate Beta-posterior updating for the transition probabilities, and two order-selection devices: $m$7 and Bayes factors $m$8 (Kimpton et al., 2022).

For the dynamic multivariate probit model, the full likelihood requires $m$9-dimensional Gaussian orthant probabilities and becomes computationally burdensome when $\{0,1\}^m$0 is moderate or large. The proposed alternative combines pseudo-likelihood, based on univariate conditional probabilities $\{0,1\}^m$1, and pairwise likelihood, based on bivariate Gaussian probabilities. Under stationarity, ergodicity, regularity, and an identifiability condition, the resulting estimators are consistent and asymptotically normal with sandwich covariance

$\{0,1\}^m$2

and analogous results hold in a panel-data setting via multi-indexed ergodic theorems and partial-sum CLTs (Franchi et al., 2024).

5. High-dimensional and network-oriented extensions

A major recent development is the move from low-dimensional binary time-series regression to sparse multivariate dependence graphs. In the BAR process, network structure is primary: each nonzero $\{0,1\}^m$3 corresponds to a directed edge $\{0,1\}^m$4, and the signed parent sets $\{0,1\}^m$5 encode whether a parent increases the next-step success probability through $\{0,1\}^m$6 or through $\{0,1\}^m$7. The sample complexity of BARObs is proportional to the mixing time and is described as nearly order-optimal because it is only a $\{0,1\}^m$8 factor away from a Fano-type lower bound,

$\{0,1\}^m$9

For bounded in-degree $m$0, the lower bound scales as $m$1, while the upper bound inherits the logarithmic mixing factor (Katselis et al., 2016).

The generalized binary VAR pushes this high-dimensional viewpoint further. For gbVAR($m$2), the parameter block matrix is

$m$3

with row-sum constraint

$m$4

The main theory is developed for gbVAR(1), where the Yule–Walker identity

$m$5

drives estimation. Because $m$6 is singular when $m$7, the procedure is row-wise: Lasso selection,

$m$8

followed by hard-thresholded support recovery and post-selection least squares. Under sparsity assumption $m$9 and further regularity conditions, the paper proves support recovery, Gaussian approximation for the max norm of the post-selection estimator, and the rate

$p$00

To avoid estimating the complicated limiting covariance directly, it proposes a second-order wild bootstrap based on perturbed lag-1 covariance errors and establishes bootstrap validity for simultaneous confidence regions (Dai et al., 29 Nov 2025).

This high-dimensional strand is methodologically distinct from logistic or probit autoregression. The 2025 paper explicitly states that gbVAR is not a GLM-based model with a logistic, probit, or complementary log-log link, but a structural binary vector autoregression that preserves binary support through random selection matrices. A common misconception is therefore that all high-dimensional binary autoregression reduces to sparse generalized linear modeling; the gbVAR construction shows that structurally different, moment-based autoregressive theories are also possible (Dai et al., 29 Nov 2025).

6. Applications, neighboring models, and recurrent limitations

The empirical domains motivating binary autoregressive models are diverse. The BAR process is proposed for opinion dynamics, epidemics, financial and biological time series, and the paper includes a pseudoreal abscisic acid signaling network with $p$01 nodes. The Ex-FI paper analyzes hourly respiratory-rate data from $p$02 expectant mothers after dichotomizing whether respiratory rate exceeds $p$03 breaths per minute. The de Bruijn paper studies daily precipitation occurrence at Eskdalemuir, UK, and annual Oxford–Cambridge boat-race winners, while the 2024 multivariate probit framework is motivated by absence-presence data in ecology collected across sites and time. The gbVAR paper applies its methods to seven NASDAQ stocks and to annual trade relations among seven countries (Katselis et al., 2016, Gao et al., 2017, Kimpton et al., 2022, Franchi et al., 2024, Dai et al., 29 Nov 2025).

Several recurring limitations cut across the literature. First, flexibility often brings exponential state growth. In the de Bruijn process, the number of free transition parameters is $p$04, so large $p$05 is data hungry and some histories may be rarely observed. In exact Fisher-information calculations for LAR($p$06), the lag-state space is also $p$07, making the recursive exact computation suitable for modest $p$08 rather than very large lag order (Kimpton et al., 2022, Gao et al., 2017).

Second, inferential convenience differs sharply by model class. Logistic autoregression offers conditional-likelihood-based parameter interpretation, but standard GLM variance calculations based on empirical Fisher information can be inadequate because the lagged regressors are endogenous. Dynamic multivariate probit handles contemporaneous dependence through a latent Gaussian covariance matrix, but full likelihood is expensive and composite likelihood sacrifices efficiency relative to full MLE. BAR yields interpretable directed signed influence structure and rapid mixing, but finite-sample recovery depends on an identifiability condition controlling harmful same-time correlations. gbVAR supports high-dimensional simultaneous inference, but its theory relies on strong row-sparsity and on asymptotically correct support recovery (Gao et al., 2017, Franchi et al., 2024, Katselis et al., 2016, Dai et al., 29 Nov 2025).

Third, binary autoregressive models should not be conflated with neighboring literatures. The BAR paper distinguishes its directed simultaneous-update dynamics from Ising/Glauber dynamics and from autologistic/logistic autoregressions. The de Bruijn paper shows that a fully history-indexed finite-state transition rule is different from a low-dimensional link-function model, even though both are autoregressive in the sense of depending on lagged outcomes. The 2024 framework further suggests that latent-threshold and observation-driven perspectives coexist within the same broad subject (Katselis et al., 2016, Kimpton et al., 2022, Franchi et al., 2024).

Taken together, these results indicate that binary autoregressive time series should be understood not as a single canonical model but as a field organized around finite-memory dependence on past binary states. The central design choices are whether dependence is parameterized through a link function, a sparse directed graph, a complete history table, a latent Gaussian threshold system, or a structural binary VAR recursion. Each choice determines the relevant ergodic theory, inferential machinery, and substantive interpretation.