Package 'RSDC'

Augment data with fitted RSDC regime information

Description

broom-style augmentation: returns the estimation residuals together with the smoothed regime probabilities and the most likely (Viterbi) regime path.

Usage

## S3 method for class 'rsdc_fit'
augment(x, ...)

Arguments

x	An "rsdc_fit" object (must carry its estimation residuals).
...	Unused; for generic compatibility.

Value

A data frame with one row per observation: the residual columns (.resid1, ...), smoothed regime probabilities (.smoothed_p1, ...) and the Viterbi state .state.

Examples

y <- scale(matrix(rnorm(200 * 2), 200, 2))
fit <- rsdc_estimate("noX", residuals = y, N = 2)
head(generics::augment(fit))

---

Plot a fitted RSDC model with ggplot2

Description

Plots the smoothed regime probabilities over the sample as stacked areas. Requires the suggested ggplot2 package.

Usage

## S3 method for class 'rsdc_fit'
autoplot(object, ...)

Arguments

object	An "rsdc_fit" object.
...	Unused; for generic compatibility.

Value

A ggplot object.

Examples

if (requireNamespace("ggplot2", quietly = TRUE)) {
  y <- scale(matrix(rnorm(200 * 2), 200, 2))
  fit <- rsdc_estimate("noX", residuals = y, N = 2)
  ggplot2::autoplot(fit)
}

---

Glance at a fitted RSDC model

Description

broom-style one-row model summary.

Usage

## S3 method for class 'rsdc_fit'
glance(x, ...)

Arguments

x	An "rsdc_fit" object.
...	Unused; for generic compatibility.

Value

A one-row data frame with method, n_regimes, logLik, AIC, BIC, df and nobs.

Examples

y <- scale(matrix(rnorm(200 * 2), 200, 2))
fit <- rsdc_estimate("noX", residuals = y, N = 2)
generics::glance(fit)

---

Green vs Brown portfolio dataset

Description

Daily returns for a green and a brown portfolios constructed following the equal-weighted 10-90 percentile approach.

Usage

data(greenbrown)

Format

A data frame with 2266 rows and three columns:

DATE

Dates ranging from 2014-01-02 to 2022-12-30.

return_green

Numeric returns for the green portfolio.

return_brown

Numeric returns for the brown portfolio.

Source

Built from the source workbook in data-raw/green-brown-ptf.xlsx (not shipped on CRAN); see data-raw/greenbrown.R.

Examples

data("greenbrown")
str(greenbrown)
head(greenbrown)

---

Media Climate Change Concerns (MCCC) index

Description

Daily Media Climate Change Concerns (MCCC) index used as the exogenous covariate driving the time-varying transition probabilities (TVTP) in the empirical illustration. Only the headline Aggregate series is retained. The index is published at the monthly frequency; here it is forward-filled across the calendar days of each month and aligned to the greenbrown trading calendar, so that mccc is row-for-row mergeable with greenbrown on DATE.

Usage

data(mccc)

Format

A data frame with 2266 rows and two columns:

DATE

Dates ranging from 2014-01-02 to 2022-12-30, matching greenbrown.

mccc

Numeric Media Climate Change Concerns (Aggregate) index value for the month containing DATE. Provided in raw (unstandardized) form; standardize over your estimation window before use as a covariate.

Source

Monthly Aggregate series in inst/extdata/mccc-monthly.csv, extracted from the Media Climate Change Concerns workbook of (Ardia et al. 2023).

References

Ardia D, Bluteau K, Boudt K, Inghelbrecht K (2023). “Climate Change Concerns and the Performance of Green vs. Brown Stocks.” Management Science, 69(12), 7607–7632. doi:10.1287/mnsc.2022.4636.

Examples

data("mccc")
str(mccc)
head(mccc)

# Row-aligned with greenbrown; build a TVTP covariate matrix (intercept + MCCC)
data("greenbrown")
stopifnot(identical(mccc$DATE, greenbrown$DATE))
X <- cbind(intercept = 1, MCCC = as.numeric(scale(mccc$mccc)))

---

Estimate Regime-Switching or Constant Correlation Model (Wrapper)

Description

Unified front-end that dispatches to one of three estimators:

• f_optim() — TVTP specification (method = "tvtp").

• f_optim_noX() — fixed transition matrix (method = "noX").

• f_optim_const() — constant correlation, single regime (method = "const").

Usage

## S3 method for class 'rsdc_fit'
print(x, digits = 4, ...)

## S3 method for class 'rsdc_fit'
coef(object, ...)

## S3 method for class 'rsdc_fit'
nobs(object, ...)

## S3 method for class 'rsdc_fit'
logLik(object, ...)

## S3 method for class 'rsdc_fit'
vcov(object, type = c("hessian", "opg", "sandwich", "bootstrap"), ...)

## S3 method for class 'rsdc_fit'
confint(
  object,
  parm,
  level = 0.95,
  type = c("hessian", "opg", "sandwich", "bootstrap"),
  ...
)

## S3 method for class 'rsdc_fit'
summary(object, ...)

## S3 method for class 'rsdc_fit'
predict(
  object,
  residuals,
  sigma_matrix,
  value_cols,
  X = NULL,
  out_of_sample = FALSE,
  ...
)

## S3 method for class 'rsdc_fit'
simulate(object, nsim = 1, seed = NULL, X = NULL, n = NULL, ...)

## S3 method for class 'rsdc_fit'
plot(x, which = c("smoothed", "filtered"), ...)

rsdc_estimate(
  method = c("tvtp", "noX", "const"),
  residuals,
  N = 2,
  X = NULL,
  out_of_sample = FALSE,
  control = list()
)

Arguments

x	An object of class "rsdc_fit".
digits	Number of significant digits for printing.
...	Further arguments passed to methods (currently unused).
object	An object of class "rsdc_fit".
type	Covariance estimator: one of "hessian", "opg", "sandwich" (numerical) or "bootstrap" (simulation-based).
parm	Vector of parameter names/indices (default: all).
level	Confidence level.
residuals	Numeric matrix $T \times K$. Typically standardized residuals/returns.
sigma_matrix	Numeric matrix of conditional standard deviations.
value_cols	Columns of sigma_matrix giving the asset order.
X	Numeric matrix $T \times p$ of exogenous covariates (required for "tvtp").
out_of_sample	Logical. If TRUE, a fixed 70/30 split is applied prior to estimation.
nsim	Unused (kept for generic compatibility; one path is returned).
seed	Optional RNG seed.
n	Series length for "noX"/"const" (defaults to the fitted T).
which	Either "smoothed" (default) or "filtered".
method	Character. One of "tvtp", "noX", "const".
N	Integer. Number of regimes. Ignored when method = "const".
control	Optional list forwarded to the backends and optimizers: seed (default 123) and do_trace (default FALSE); optimizer settings itermax, NP, parallelType, steptol (DEoptim) and maxit (optim); compute_se (default TRUE) to toggle the observed-information standard errors; and n_starts (default 1) to repeat the global+local search from several seeds (seed, seed+1, ...) and keep the highest-likelihood fit. Because DEoptim is a stochastic global search, its result depends on the seed; n_starts > 1 is therefore primarily a stability diagnostic (not a replacement for the global search): the returned object carries start_logliks (one log-likelihood per start) so you can check whether the optimum is reproducible across seeds. A warm start (start) disables multi-start.

Details

• Method selection: match.arg() validates method.

• Inputs: "tvtp" requires non-NULL X; N is ignored for "const".

• Split: If out_of_sample = TRUE, the first 70% is used for fitting.

Value

An object of class "rsdc_fit": a list with components

transition_matrix

Estimated transition matrix ($1 \times 1$ for "const"; for "tvtp", evaluated at the in-sample covariate means).

correlations

Regime lower-triangular correlations (rows ordered by ascending mean).

covariances

Array $K \times K \times N$ of regime correlation matrices.

log_likelihood

In-sample log-likelihood; full-sample when out_of_sample = FALSE.

log_likelihood_oos

OOS log-likelihood on held-out 30 percent, or NULL.

beta

TVTP coefficients (only for "tvtp"; NULL otherwise).

par, coefficients

Named vector of estimated parameters in objective order.

vcov

Observed-information variance-covariance at the MLE, or NULL if the Hessian is singular/unavailable (e.g. at a bound).

se

Named standard errors (sqrt(diag(vcov))), or NULL.

convergence

optim convergence code (0 = converged).

npar, nobs

Number of free parameters and observations used in estimation.

method, N, K, p, call

Fit metadata.

Standard methods are provided: print, summary, coef, logLik (so AIC/BIC work), nobs, vcov, confint, predict, and simulate.

Functions

• print(rsdc_fit): Compact printer for a fitted model.

• coef(rsdc_fit): Coefficient vector (named, in objective order).

• nobs(rsdc_fit): Number of observations used in estimation.

• logLik(rsdc_fit): Log-likelihood (carries df and nobs so stats::AIC/BIC work out of the box).

• vcov(rsdc_fit): Variance-covariance matrix of the estimates. The numerical estimators are type = "hessian" (default, inverse observed information), "opg" (outer product of gradients) and "sandwich" (QML/robust $H^{-1} (\sum_t s_t s_t') H^{-1}$); the simulation-based estimator is "bootstrap", which calls rsdc_bootstrap (pass B and seed through ...). The bootstrap is recomputed on each call.

• confint(rsdc_fit): Wald confidence intervals from the chosen covariance (numerical or "bootstrap"). For bootstrap percentile intervals instead of Wald, use rsdc_bootstrap directly.

• summary(rsdc_fit): Summary with a coefficient table (estimate, SE, z, p).

• predict(rsdc_fit): Forecast from a fitted model (wraps rsdc_forecast); supply the residuals, conditional volatilities, and (for "tvtp") the covariate matrix X.

• simulate(rsdc_fit): Simulate from a fitted model. For "tvtp" supply a covariate matrix X (its row count sets the length); for "noX" the fixed transition matrix is used; for "const" a single regime is drawn.

• plot(rsdc_fit): Plot the smoothed regime probabilities (one panel per regime) for a fitted "noX"/"tvtp" model. which = "filtered" plots filtered probabilities instead.

References

Mullen K, Ardia D, Gil D, Windover D, Ulrich J (2011). “DEoptim: An R Package for Global Optimization by Differential Evolution.” Journal of Statistical Software, 40(6), 1–26. doi:10.18637/jss.v040.i06.

Hamilton JD (1989). “A New Approach to the Economic Analysis of Nonstationary Time Series and the Business Cycle.” Econometrica, 57(2), 357–384. doi:10.2307/1912559.

Pelletier D (2006). “Regime switching for dynamic correlations.” Journal of Econometrics, 131(1-2), 445–473. doi:10.1016/j.jeconom.2005.01.013.

See Also

rsdc_hamilton and rsdc_likelihood.

Examples

y <- scale(matrix(rnorm(100 * 3), 100, 3))
fit <- rsdc_estimate("noX", residuals = y, N = 2)
print(fit)
summary(fit)
c(AIC = AIC(fit), BIC = BIC(fit))
X <- cbind(1, scale(seq_len(nrow(y))))
fit_tvtp <- rsdc_estimate("tvtp", residuals = y, N = 2, X = X)
coef(fit_tvtp)

---

RSDC: Regime-Switching Dynamic Correlation Models

Description

The RSDC package provides a comprehensive framework for modeling, estimating, and forecasting correlation structures in multivariate time series under regime-switching dynamics. It supports both fixed transition probabilities and time-varying transition probabilities (TVTP) driven by exogenous variables.

The methodology is particularly suited to empirical asset pricing and portfolio management applications, enabling users to incorporate macroeconomic, financial, or climate-related predictors into the regime dynamics. The package integrates the full workflow — from model estimation to covariance matrix reconstruction and portfolio optimization — in a single, reproducible pipeline.

Main Features

• Model estimation and filtering: rsdc_hamilton (Hamilton filter), rsdc_likelihood (likelihood computation), rsdc_estimate (parameter estimation).

• Correlation and covariance forecasting: rsdc_forecast (in-sample / 70-30 path), rsdc_forecast_ahead (multi-step-ahead), rsdc_corr_bands (parameter-uncertainty bands), rsdc_viterbi (most likely regime path).

• Portfolio construction: rsdc_minvar (minimum-variance portfolios), rsdc_maxdiv (maximum-diversification portfolios).

• Simulation: rsdc_simulate (simulate TVTP regime-switching series).

• Fitted-model methods: rsdc_estimate() returns an object of class "rsdc_fit" with print, summary, coef, logLik (so AIC/BIC work), nobs, vcov, confint, predict, and simulate methods.

Authors

David Ardia, Benjamin Seguin and Roosevelt Ymele Nguemo. David Ardia is the maintainer (also copyright holder and funder).

References

Engle RF (2002). “Dynamic conditional correlation: A simple class of multivariate generalized autoregressive conditional heteroskedasticity models.” Journal of Business & Economic Statistics, 20(3), 339–350. doi:10.1198/073500102288618487.

Hamilton JD (1989). “A New Approach to the Economic Analysis of Nonstationary Time Series and the Business Cycle.” Econometrica, 57(2), 357–384. doi:10.2307/1912559.

Pelletier D (2006). “Regime switching for dynamic correlations.” Journal of Econometrics, 131(1-2), 445–473. doi:10.1016/j.jeconom.2005.01.013.

---

Parametric bootstrap standard errors for a fitted RSDC model

Description

Simulates B data sets from the fitted data-generating process, re-estimates each one — warm-started at the original MLE so the global search is skipped — and returns the bootstrap covariance, standard errors and percentile intervals of the parameters. This complements the analytic (Hessian/OPG/sandwich) covariances from vcov.rsdc_fit and, being simulation-based, respects the parameter bounds ($\rho \in (-1,1)$, transition probabilities in $[0,1]$).

Usage

rsdc_bootstrap(object, B = 199, X = NULL, seed = NULL, level = 0.95)

Arguments

object	An object of class "rsdc_fit".
B	Integer number of bootstrap replications (default 199).
X	Covariate matrix for "tvtp" (defaults to the stored estimation covariates).
seed	Optional integer seed for reproducibility.
level	Confidence level for the percentile intervals (default 0.95).

Details

For method = "tvtp" the data are simulated with rsdc_simulate from the fitted beta and covariate matrix X; for "noX" from the fitted (fixed) transition matrix; for "const" as i.i.d. draws from the single regime. Each replicate is re-estimated with the same method/N and control = list(start = coef(object), compute_se = FALSE), which warm-starts the local optimizer at the original estimates (fast, deterministic given a seed). Replicates whose re-estimation fails are dropped.

Value

A list with replicates ($B' \times$ npar matrix of successful re-estimates), vcov, se, ci (percentile intervals) and the effective number of replicates B.

See Also

vcov.rsdc_fit, rsdc_simulate, rsdc_estimate

Examples

set.seed(1)
X <- cbind(1, as.numeric(scale(seq_len(300))))
y <- scale(matrix(rnorm(300 * 2), 300, 2))
fit <- rsdc_estimate("tvtp", residuals = y, N = 2, X = X)
bs <- rsdc_bootstrap(fit, B = 50, seed = 1)
bs$se
bs$ci

---

Uncertainty bands for the predicted correlation path

Description

Propagates parameter uncertainty into the (smoothed-probability-weighted) predicted correlation path of a fitted "rsdc_fit" model. Parameter vectors are drawn from the asymptotic sampling distribution $\mathcal{N}(\hat\theta, \widehat{\mathrm{Var}}(\hat\theta))$ (using the covariance returned by vcov.rsdc_fit); for each draw the Hamilton filter/smoother is rerun and the predicted correlation path recomputed. The pointwise quantiles of these paths form the band.

Usage

rsdc_corr_bands(
  object,
  residuals = NULL,
  X = NULL,
  B = 500L,
  level = 0.95,
  seed = NULL
)

Arguments

object	An "rsdc_fit" object (from rsdc_estimate) estimated with control = list(compute_se = TRUE) so a covariance is available.
residuals	Optional $T \times K$ residual matrix; defaults to the matrix stored on the fit.
X	Optional $T \times p$ covariate matrix ("tvtp" only); defaults to the matrix stored on the fit.
B	Integer. Number of parameter draws (default 500).
level	Confidence level for the pointwise band (default 0.95).
seed	Optional integer seed for reproducibility.

Details

This reflects estimation (parameter) uncertainty in the correlation parameters and transition dynamics; it does not add the irreducible regime-classification uncertainty already summarized by the smoothed probabilities.

Draws are taken from the unconstrained Gaussian approximation; draws that yield an invalid model (e.g. $|\rho| \ge 1$, a non-positive-definite correlation matrix, or an invalid transition matrix) are dropped with a warning and the band is computed from the remaining draws. When the estimate sits close to a parameter bound this truncates the sampling distribution, so the reported band can be somewhat too narrow there; the attribute B_used reports how many draws survived.

Value

A list with one element per correlation pair ($C = K(K-1)/2$ of them), each a $T \times 3$ matrix with columns fit, lower, upper; plus attributes level and B_used (draws that yielded a valid path).

See Also

rsdc_forecast, rsdc_bootstrap

Examples

y <- scale(matrix(rnorm(250 * 2), 250, 2))
fit <- rsdc_estimate("noX", residuals = y, N = 2)
bands <- rsdc_corr_bands(fit, B = 100)
head(bands[[1]])

---

Forecast Covariance/Correlation Paths from an RSDC Model

Description

Generates per-period correlation and covariance matrices from a fitted model: "const" (constant correlation), "noX" (fixed transition matrix), or "tvtp" (time-varying transition probabilities).

Usage

rsdc_forecast(
  method = c("tvtp", "noX", "const"),
  N,
  residuals,
  X = NULL,
  final_params,
  sigma_matrix,
  value_cols,
  out_of_sample = FALSE,
  control = list()
)

Arguments

method	Character. One of "tvtp", "noX", "const".
N	Integer. Number of regimes (ignored for "const").
residuals	Numeric matrix $T \times K$ used to compute correlations or run the filter.
X	Optional numeric matrix $T \times p$ (required for "tvtp").
final_params	List with fitted parameters (e.g., from rsdc_estimate): must include correlations, and either transition_matrix ("noX") or beta ("tvtp"); include log_likelihood for BIC computation.
sigma_matrix	Numeric matrix $T \times K$ of forecasted standard deviations.
value_cols	Character/integer vector of columns in sigma_matrix that define asset order.
out_of_sample	Logical. If TRUE, use a fixed 70/30 split; otherwise use the whole sample.
control	Optional list; supports threshold (in (0,1), default 0.7).

Details

• Forecast horizon: If out_of_sample = TRUE, filter on the first threshold fraction and forecast on the remainder.

• Correlation paths:

  • "const" — constant correlation from final_params$correlations[1,], repeated across time.

  • "noX"/"tvtp" — regime-probability weighted average of regime correlations. Out-of-sample (out_of_sample = TRUE) the weights are the filtered probabilities $\Pr(S_t \mid \Omega_t)$, which use no future information beyond time $t$ (the OOS window is initialized from the in-sample terminal filtered state). Note these condition on the contemporaneous $U_t$, so they are a real-time nowcast, not a strictly pre-observation one-step-ahead forecast; for the latter use rsdc_forecast_ahead. The full-sample case uses smoothed probabilities $\Pr(S_t \mid \Omega_T)$. The returned smoothed_probs field always contains the smoothed probabilities.

• Covariance build: Reconstruct $R_t$ from the pairwise vector (columns ordered by combn(K, 2)), set $D_t = \mathrm{diag}(\sigma_{t,1},\dots,\sigma_{t,K})$, and $\Sigma_t = D_t R_t D_t$.

• BIC: Parameter count $k$ is N * ncol(X) + N * K * (K - 1) / 2 for "tvtp" with $N=2$, N * (N-1) * ncol(X) + N * K * (K - 1) / 2 for "tvtp" with $N=3$, N * (N - 1) + N * K * (K - 1) / 2 for "noX", and K * (K - 1) / 2 for "const".

Value

smoothed_probs

$N \times T^\ast$ smoothed probabilities ("noX"/"tvtp" only).

sigma_matrix

$T^\ast \times K$ slice aligned to the forecast horizon.

cov_matrices

List of $K \times K$ covariance matrices $\Sigma_t = D_t R_t D_t$.

predicted_correlations

$T^\ast \times \binom{K}{2}$ pairwise correlations in combn(K, 2) order.

BIC

Information criterion $\log(n)\,k - 2\,\ell$. When out_of_sample = TRUE and final_params$log_likelihood_oos is non-NULL, uses the OOS log-likelihood and OOS sample size; otherwise uses the IS values. Caveat: with out_of_sample = TRUE this is an out-of-sample predictive score (held-out log-likelihood penalized by the IS parameter count), not a conventional in-sample BIC; use it only to compare models fit on the same split, and do not interpret it as a textbook BIC. The backends (rsdc_estimate) hardcode a 70/30 split, so the OOS log-likelihood is only fully consistent with rsdc_forecast when the default threshold = 0.7 is used.

y

$T^\ast \times K$ residual matrix aligned to the forecast horizon.

See Also

rsdc_hamilton, rsdc_estimate, rsdc_minvar, rsdc_maxdiv

Examples

set.seed(123)
T <- 60; K <- 3; N <- 2
y <- scale(matrix(rnorm(T*K), T, K))
vols <- matrix(0.2 + 0.05*abs(sin(seq_len(T)/7)), T, K)
rho <- rbind(c(0.10, 0.05, 0.00), c(0.60, 0.40, 0.30))
Pfix <- matrix(c(0.9, 0.1, 0.2, 0.8), 2, 2, byrow = TRUE)
rsdc_forecast("noX", N, y, NULL,
              list(correlations = rho, transition_matrix = Pfix, log_likelihood = -200),
              vols, 1:K)

---

Multi-step-ahead regime and correlation forecasts

Description

Produces genuine $h$-step-ahead forecasts of the regime distribution and the implied correlation matrix from the end of the estimation sample. Starting from the terminal filtered regime probabilities $\xi_T$, the regime distribution is propagated through the Markov chain, $\pi_{T+k} = \pi_{T+k-1} P_{T+k}$, and the forecast correlations are the regime-probability-weighted averages of the regime correlation matrices, $\hat\rho_{T+k} = \sum_s \pi_{T+k}(s)\,\rho^{(s)}$.

Usage

rsdc_forecast_ahead(
  object,
  horizon = 10L,
  residuals = NULL,
  X = NULL,
  X_future = NULL
)

Arguments

object	An "rsdc_fit" object from rsdc_estimate.
horizon	Integer forecast horizon $h \ge 1$ (default 10).
residuals	Optional $T \times K$ residuals used to obtain the terminal filtered state; defaults to the matrix stored on the fit.
X	Optional in-sample covariate matrix ("tvtp"); defaults to the fit's stored matrix.
X_future	Optional $h \times p$ matrix of future covariates ("tvtp" only). If NULL, the last observed covariate row is reused at every horizon.

Details

For "noX"/"const" the transition matrix is constant, so $\pi_{T+k} = \xi_T P^k$. For "tvtp" each step uses a covariate row: supply future covariates in X_future (an $h \times p$ matrix); if it is NULL the last in-sample covariate row is held constant and a message is emitted.

Value

A list with:

horizon

Integer vector 1:h.

regime_probs

$h \times N$ matrix of forecast regime probabilities.

predicted_correlations

$h \times C$ matrix ($C = K(K-1)/2$) of forecast correlations, one column per asset pair.

See Also

rsdc_forecast for the in-sample / 70-30 path.

Examples

y <- scale(matrix(rnorm(200 * 2), 200, 2))
fit <- rsdc_estimate("noX", residuals = y, N = 2)
rsdc_forecast_ahead(fit, horizon = 5)

---

Hamilton Filter (Fixed P or TVTP)

Description

Runs the Hamilton (1989) filter for a multivariate regime-switching correlation model. Supports either a fixed (time-invariant) transition matrix $P$ or time-varying transition probabilities (TVTP) built from exogenous covariates X: a logistic link for $N=2$ or a softmax link for $N \ge 3$. Returns filtered/smoothed regime probabilities and the log-likelihood.

Usage

rsdc_hamilton(
  y,
  X = NULL,
  beta = NULL,
  rho_matrix,
  K,
  N,
  P = NULL,
  xi_init = NULL,
  engine = c("cpp", "r")
)

Arguments

y	Numeric matrix $T \times K$ of observations (e.g., standardized residuals/returns). Columns are treated as mean-zero with unit variance; only the correlation structure is modeled.
X	Optional numeric matrix $T \times p$ of covariates for TVTP. Required if beta is supplied.
beta	Optional numeric matrix of TVTP coefficients. For $N=2$: an $N \times p$ matrix; row $i$ gives the logistic coefficient vector so that $p_{ii,t} = \mathrm{logit}^{-1}(X_t^\top \beta_i)$. For $N \ge 3$: an $N \times (N-1)p$ matrix packed row-wise; row $i$ contains $N-1$ consecutive length-$p$ softmax vectors; the $N$-th logit is fixed at 0 (reference category).
rho_matrix	Numeric matrix $N \times C$ of regime correlation parameters, where $C = K(K-1)/2$. Each row is the lower-triangular part (by lower.tri) of a regime's correlation matrix.
K	Integer. Number of observed series (columns of y).
N	Integer. Number of regimes.
P	Optional $N \times N$ fixed transition matrix. Used only when X or beta is NULL.
xi_init	Optional numeric vector of length N. If supplied, the forward filter is initialized from this vector (normalized internally to sum to 1) instead of the default uniform $1/N$. Intended for two-pass out-of-sample forecasting: run the filter on the in-sample period first, extract the terminal column of filtered_probs, and pass it here to initialize the out-of-sample filter run. Must be non-negative and finite.
engine	Character; "cpp" (default) runs the filter/smoother in C++ (RcppArmadillo), "r" the pure-R reference. Both give identical results (used for equivalence testing).

Details

• Correlation rebuild: For regime $m$, a correlation matrix $R_m$ is reconstructed from rho_matrix[m, ] (lower-triangular fill + symmetrization). Non-PD proposals are penalized.

• Transition dynamics:

  • Fixed P: If X or beta is missing, a constant $P$ is used (user-provided via P; otherwise uniform $1/N$ rows).

  • TVTP: With X and beta: for $N=2$, $p_{ii,t} = \mathrm{logit}^{-1}(X_t^\top \beta_i)$ and $p_{ij,t} = 1 - p_{ii,t}$ ($j \ne i$); for $N \ge 3$, a softmax over $N-1$ free logit vectors per row with the $N$-th logit fixed at 0 (reference category); log-sum-exp stabilized.

• Numerical safeguards: A small ridge is added before inversion; if filtering degenerates at a time step, log_likelihood = -Inf is returned.

Value

A list with:

filtered_probs

$N \times T$ matrix of filtered probabilities $\Pr(S_t = j \mid \Omega_t)$.

smoothed_probs

$N \times T$ matrix of smoothed probabilities $\Pr(S_t = j \mid \Omega_T)$.

log_likelihood

Scalar log-likelihood of the model given y.

Note

For $N=2$, TVTP uses a logistic link on the diagonal; off-diagonals follow by complement. For $N \ge 3$, a full softmax is used with $(N-1)$ free logit vectors per row; all entries are independently determined.

X-timing: uses X[t, ] (contemporaneous) to form $P_t$, consistent with rsdc_simulate.

References

Hamilton JD (1989). “A New Approach to the Economic Analysis of Nonstationary Time Series and the Business Cycle.” Econometrica, 57(2), 357–384. doi:10.2307/1912559.

See Also

rsdc_likelihood and rsdc_estimate.

Examples

set.seed(1)
T <- 50; K <- 3; N <- 2
y <- scale(matrix(rnorm(T * K), T, K), center = TRUE, scale = TRUE)

# Example rho: two regimes with different average correlations
rho <- rbind(c(0.10, 0.05, 0.00),
             c(0.60, 0.40, 0.30))  # lower-tri order for K=3

# Fixed-P filtering
Pfix <- matrix(c(0.9, 0.1,
                 0.2, 0.8), nrow = 2, byrow = TRUE)
out_fix <- rsdc_hamilton(y = y, X = NULL, beta = NULL,
                         rho_matrix = rho, K = K, N = N, P = Pfix)
str(out_fix$filtered_probs)

# TVTP filtering (include an intercept yourself)
X <- cbind(1, scale(seq_len(T)))
beta <- rbind(c(1.2, 0.0),
              c(0.8, -0.1))
out_tvtp <- rsdc_hamilton(y = y, X = X, beta = beta,
                          rho_matrix = rho, K = K, N = N)
out_tvtp$log_likelihood

---

Negative Log-Likelihood for Regime-Switching Correlation Models

Description

Computes the negative log-likelihood for a multivariate correlation-only regime-switching model, with either a fixed (time-invariant) transition matrix or time-varying transition probabilities (TVTP) driven by exogenous covariates. Likelihood evaluation uses the Hamilton (1989) filter.

Usage

rsdc_likelihood(params, y, exog = NULL, K, N)

Arguments

params	Numeric vector of model parameters packed as: No exogenous covariates (exog = NULL): first $N(N-1)$ transition parameters (for the fixed transition matrix), followed by $N \times K(K-1)/2$ correlation parameters, stacked row-wise by regime in lower.tri order. With exogenous covariates (exog provided): first $N \times p$ TVTP coefficients (beta) for $N=2$, or $N \times (N-1)p$ coefficients for $N \ge 3$ (row $i$ holds $N-1$ length-$p$ softmax vectors; last logit = 0 reference), followed by $N \times K(K-1)/2$ correlation parameters, stacked row-wise by regime in lower.tri order.
y	Numeric matrix $T \times K$ of observations (e.g., standardized residuals). Columns are treated as mean-zero, unit-variance; only correlation is modeled.
exog	Optional numeric matrix $T \times p$ of exogenous covariates. If supplied, a TVTP specification is used.
K	Integer. Number of observed series (columns of y).
N	Integer. Number of regimes.

Details

• Transition dynamics:

  • Fixed P (no exog): params begins with transition parameters. For $N=2$, the implementation maps them to $P=\begin{pmatrix} p_{11} & 1-p_{11}\\ 1-p_{22} & p_{22}\end{pmatrix}$.

  • TVTP: with exog, for $N=2$, $p_{ii,t} = \mathrm{logit}^{-1}(X_t^\top \beta_i)$ and $p_{ij,t} = 1 - p_{ii,t}$ ($j \ne i$); for $N \ge 3$, a softmax over $N-1$ free logit vectors per row with the $N$-th logit fixed at 0 (reference category).

• Correlation build: per regime, the lower-triangular vector is filled into a symmetric correlation matrix. Non-positive-definite proposals or $|\rho|\ge 1$ are penalized via a large objective value.

• Evaluation: delegates to rsdc_hamilton; if the filter returns log_likelihood = -Inf, a large penalty is returned.

Value

Numeric scalar: the negative log-likelihood to be minimized by an optimizer.

Note

The function is written for use inside optimizers; it performs inexpensive validation and returns large penalties for invalid parameterizations instead of stopping with errors.

See Also

rsdc_hamilton (filter), optim, and DEoptim

Examples

# Small toy example (N = 2, K = 3), fixed P (no exog)
set.seed(1)
T <- 40; K <- 3; N <- 2
y <- scale(matrix(rnorm(T * K), T, K), center = TRUE, scale = TRUE)

# Pack parameters: trans (p11, p22), then rho by regime (lower-tri order)
p11 <- 0.9; p22 <- 0.8
rho1 <- c(0.10, 0.05, 0.00)  # (2,1), (3,1), (3,2)
rho2 <- c(0.60, 0.40, 0.30)
params <- c(p11, p22, rho1, rho2)

nll <- rsdc_likelihood(params, y = y, exog = NULL, K = K, N = N)
nll

# TVTP example: add X and beta (pack beta row-wise, then rho)
X <- cbind(1, scale(seq_len(T)))
beta <- rbind(c(1.2, 0.0),
              c(0.8, -0.1))
params_tvtp <- c(as.vector(t(beta)), rho1, rho2)
nll_tvtp <- rsdc_likelihood(params_tvtp, y = y, exog = X, K = K, N = N)
nll_tvtp

---

Maximum-Diversification Portfolio (Rolling Weights)

Description

Computes rolling maximum-diversification (MaxDiv) portfolio weights from a sequence of per-period covariance matrices implied by forecasted volatilities and correlations. Falls back to equal weights if the nonlinear solver fails.

Usage

rsdc_maxdiv(
  sigma_matrix,
  value_cols,
  predicted_corr,
  y,
  long_only = TRUE,
  lag = FALSE
)

Arguments

sigma_matrix	Numeric matrix $T \times K$ of forecasted standard deviations.
value_cols	Character/integer vector naming columns in sigma_matrix (asset order).
predicted_corr	Numeric matrix $T \times \binom{K}{2}$ of pairwise correlations in combn(K, 2) column order.
y	Numeric matrix $T \times K$ of asset returns (for realized stats).
long_only	Logical. If TRUE, impose $w \ge 0$ and $\sum_i w_i = 1$; otherwise bounds are $-1 \le w_i \le 1$ with $\sum_i w_i = 1$.
lag	Logical. If FALSE (default), the period-$t$ weights (built from $\Sigma_t$) are applied to the same period's returns y[t, ] (in-sample evaluation). If TRUE, weights chosen at $t-1$ earn the period-$t$ return (sum(y[t, ] * weights[t-1, ])) and the first period's return is NA; use lag = TRUE for a look-ahead-free out-of-sample backtest.

Details

• Covariance build: For each $t$, reconstruct $R_t$ from the pairwise vector; set $D_t=\mathrm{diag}(\sigma_{t,1},\dots,\sigma_{t,K})$ and $\Sigma_t = D_t R_t D_t$.

• Objective (MaxDiv): maximize $\mathrm{DR}(w) = \frac{\sum_i w_i \sigma_{t,i}}{\sqrt{w^\top \Sigma_t w}}$ subject to $\sum_i w_i = 1$ and bounds on $w$. Implemented by minimizing the negative ratio.

• Solver: Rsolnp::solnp with equality $\sum_i w_i = 1$ and bounds by long_only; on error, weights default to $1/K$.

Value

weights

$T \times K$ matrix of weights.

returns

Vector of realized portfolio returns sum(y[t,] * weights[t,]).

diversification_ratios

Vector of realized diversification ratios.

mean_diversification

Average diversification ratio.

K

Number of assets.

assets

Asset names.

volatility

Standard deviation of realized portfolio returns.

See Also

rsdc_minvar, solnp

Examples

# Toy example with K = 3
if (requireNamespace("Rsolnp", quietly = TRUE)) {
  T <- 50; K <- 3
  set.seed(42)
  vols <- matrix(0.2 + 0.05*abs(sin(seq_len(T)/7)), T, K)
  colnames(vols) <- paste0("A", 1:K)
  # simple, stationary correlations (order: (2,1), (3,1), (3,2))
  pred_corr <- cbind(rep(0.20, T), rep(0.10, T), rep(0.05, T))
  rets <- matrix(rnorm(T*K, sd = 0.01), T, K); colnames(rets) <- colnames(vols)

  mx <- rsdc_maxdiv(sigma_matrix   = vols,
                    value_cols     = colnames(vols),
                    predicted_corr = pred_corr,
                    y              = rets,
                    long_only      = TRUE)
  head(mx$weights)
  mx$mean_diversification
}

---

Minimum-Variance Portfolio (Rolling Weights)

Description

Computes rolling minimum-variance (MV) portfolio weights from a sequence of per-period covariance matrices implied by forecasted volatilities and pairwise correlations. Supports long-only or unconstrained MV. If the QP solver fails at a time step, the routine falls back to equal weights.

Usage

rsdc_minvar(
  sigma_matrix,
  value_cols,
  predicted_corr,
  y,
  long_only = TRUE,
  lag = FALSE
)

Arguments

sigma_matrix	Numeric matrix $T \times K$ of forecasted volatilities (standard deviations), one column per asset.
value_cols	Character or integer vector giving the columns in sigma_matrix to use as assets (order defines the asset order).
predicted_corr	Numeric matrix $T \times P$ of pairwise correlations, where $P = \binom{K}{2}$ and the columns correspond to combn(K, 2) order.
y	Numeric matrix $T \times K$ of asset returns aligned with sigma_matrix. Used only to compute the realized portfolio volatility.
long_only	Logical. If TRUE (default), imposes long-only MV with the full-investment constraint $\sum_i w_i = 1$ and $w_i \ge 0$. If FALSE, solves unconstrained MV with only $\sum_i w_i = 1$.
lag	Logical. If FALSE (default), the period-$t$ weights (built from $\Sigma_t$) are applied to the same period's returns y[t, ] (in-sample evaluation). If TRUE, weights chosen at $t-1$ earn the period-$t$ return (sum(y[t, ] * weights[t-1, ])) and the first period's return is NA; use lag = TRUE for a look-ahead-free out-of-sample backtest.

Details

• Covariance build: For each $t$, a correlation matrix $R_t$ is reconstructed from predicted\_corr[t, ] (columns in combn(K, 2) order) by placing each pairwise correlation in the corresponding off-diagonal entries of a $K \times K$ identity matrix. Let $D_t = \mathrm{diag}(\sigma_{t,1},\dots,\sigma_{t,K})$ and $\Sigma_t = D_t R_t D_t$.

• Optimization: Minimize $\tfrac{1}{2} w^\top \Sigma_t w$ subject to $\mathbf{1}^\top w = 1$ and, if long_only, $w \ge 0$ (solved with quadprog::solve.QP).

• Failure handling: If the QP fails at time \(t\), weights default to equal allocation \(w_i = 1/K\).

Value

An object of class "minvar_portfolio":

weights

$T \times K$ matrix of MV weights (one row per time).

cov_matrices

List of length $T$ with the per-period $K \times K$ covariance matrices.

volatility

Realized standard deviation of portfolio returns (same units as y).

y

The input y matrix (coerced to $T \times K$).

K

Number of assets.

See Also

rsdc_maxdiv (maximum diversification), solve.QP

Examples

# Toy example with K = 3 (requires the suggested 'quadprog' package)
if (requireNamespace("quadprog", quietly = TRUE)) {
  T <- 50; K <- 3
  set.seed(42)
  vols <- matrix(0.2 + 0.05*abs(sin(seq_len(T)/7)), T, K)
  colnames(vols) <- paste0("A", 1:K)
  # simple, stationary correlations
  pred_corr <- cbind(rep(0.20, T), rep(0.10, T), rep(0.05, T))  # order: (2,1), (3,1), (3,2)
  rets <- matrix(rnorm(T*K, sd = 0.01), T, K); colnames(rets) <- colnames(vols)

  mv <- rsdc_minvar(sigma_matrix  = vols,
                    value_cols    = colnames(vols),
                    predicted_corr= pred_corr,
                    y             = rets,
                    long_only     = TRUE)
  head(mv$weights)
  mv$volatility
}

---

Simulate Multivariate Regime-Switching Data (TVTP)

Description

Simulates a multivariate time series from a regime-switching model with time-varying transition probabilities (TVTP) driven by covariates X. Transition probabilities are generated via a multinomial logistic (softmax) link; observations are drawn from regime-specific Gaussian distributions.

Usage

rsdc_simulate(n, X, beta, mu, sigma, N, seed = NULL)

Arguments

n	Integer. Number of time steps to simulate; must be at least 2.
X	Numeric matrix $n \times p$ of covariates used to form the transition probabilities. Row X[1, ] forms $P_1$, used only to draw the initial state from $t(P_1)\,\pi_0$; rows X[t, ] for $t \ge 2$ form $P_t$.
beta	Numeric matrix $N \times (N-1) \cdot p$. Transition coefficients packed row-wise by regime, matching the parameterization of rsdc_hamilton and rsdc_estimate. For $N=2$, each row is a length-$p$ logistic vector: $p_{ii,t} = \mathrm{logit}^{-1}(X_t^\top \beta_i)$. For $N=3$, each row packs $N-1=2$ softmax vectors of length $p$; the $N$-th logit is the reference (zero).
mu	Numeric matrix $N \times K$. Regime-specific mean vectors.
sigma	Numeric array $K \times K \times N$. Regime-specific covariance (here, correlation/variance) matrices; each $K \times K$ slice must be symmetric positive definite.
N	Integer. Number of regimes.
seed	Optional integer. If supplied, sets the RNG seed for reproducibility.

Details

• Initial state and first draw: The initial regime $S_1$ is drawn from $t(P_1)\,\pi_0$ with $\pi_0 = (1/N,\dots,1/N)$ and $P_1$ built from $X_1$ via the same logistic/softmax link, matching the Hamilton-filter t=1 prior in rsdc_hamilton; the first observation $y_1$ is drawn from $\mathcal{N}(\mu_{S_1}, \Sigma_{S_1})$.

• TVTP transition: For $t \ge 2$ and $N=2$, $p_{ii,t} = \mathrm{logit}^{-1}(X_t^\top \beta_i)$; for $N \ge 3$, a softmax is used with $N-1$ free logit vectors per row and the $N$-th logit fixed at zero (reference category). Both use log-sum-exp stabilization.

• Sampling: Given $S_{t-1}$, draw $S_t$ from the categorical distribution with probabilities $P_t(S_{t-1}, \cdot)$ and $y_t \sim \mathcal{N}(\mu_{S_t}, \Sigma_{S_t})$.

• First-slice NA: transition_matrices[,,1] is always NA. At $t=1$ there is no predecessor state, so no transition matrix is computed. The assignment loop begins at $t=2$. Users iterating over the returned array should skip the first slice or use transition_matrices[,,2:n].

Value

A list with:

states

Integer vector of length n; the simulated regime index at each time.

observations

Numeric matrix $n \times K$; the simulated observations.

transition_matrices

Array $N \times N \times n$; the transition matrix $P_t$ used at each time step (with $P_1$ undefined by construction; see Details).

Note

Requires mvtnorm for multivariate normal sampling (called as mvtnorm::rmvnorm).

X-timing: uses X[t, ] (contemporaneous) to form $P_t$, consistent with rsdc_hamilton.

See Also

rsdc_hamilton (filter/evaluation), rsdc_estimate (estimators), rsdc_forecast (forecasting)

Examples

set.seed(123)
n <- 200; K <- 3; N <- 2; p <- 2
X <- cbind(1, scale(seq_len(n)))

beta <- matrix(0, nrow = N, ncol = (N - 1) * p)
beta[1, ] <- c(1.2,  0.0)
beta[2, ] <- c(1.0, -0.1)

mu <- rbind(c(0, 0, 0),
            c(0, 0, 0))
rho <- rbind(c(0.10, 0.05, 0.00),
             c(0.60, 0.40, 0.30))
Sig <- array(0, dim = c(K, K, N))
for (m in 1:N) {
  R <- diag(K); R[lower.tri(R)] <- rho[m, ]; R[upper.tri(R)] <- t(R)[upper.tri(R)]
  Sig[, , m] <- R
}
sim <- rsdc_simulate(n = n, X = X, beta = beta, mu = mu, sigma = Sig, N = N, seed = 99)

---

Most likely regime path (Viterbi decoding)

Description

Returns the single most likely sequence of regimes $(S_1, \dots, S_T)$ under a fitted "rsdc_fit" model, via the Viterbi algorithm. This is the maximum a posteriori joint state path and complements the marginal filtered/smoothed regime probabilities from rsdc_forecast: the smoothed probabilities are decoded pointwise, whereas Viterbi enforces a globally coherent transition path.

Usage

rsdc_viterbi(object, residuals = NULL, X = NULL)

Arguments

object	An "rsdc_fit" object from rsdc_estimate.
residuals	Optional $T \times K$ residual matrix; defaults to the matrix stored on the fit.
X	Optional $T \times p$ covariate matrix ("tvtp" only); defaults to the matrix stored on the fit.

Value

An integer vector of length $T$ giving the most likely regime at each time, in the model's ascending-correlation regime labeling (1 = lowest-correlation state).

See Also

rsdc_forecast for marginal regime probabilities.

Examples

y <- scale(matrix(rnorm(200 * 2), 200, 2))
fit <- rsdc_estimate("noX", residuals = y, N = 2)
table(rsdc_viterbi(fit))

---

Tidy a fitted RSDC model

Description

broom-style tidier returning a one-row-per-parameter coefficient table.

Usage

## S3 method for class 'rsdc_fit'
tidy(x, ...)

Arguments

x	An "rsdc_fit" object.
...	Unused; for generic compatibility.

Value

A data frame with columns term, estimate, std.error, statistic and p.value.

Examples

y <- scale(matrix(rnorm(200 * 2), 200, 2))
fit <- rsdc_estimate("noX", residuals = y, N = 2)
generics::tidy(fit)

---