CepReg: Cepstral Regression

In CepReg: A Cepstral Model for Covariate-Dependent Time Series

Cepstral Regression

Description

Performs cepstral regression to model frequency domain relationships between a functional response and scalar covariates. Supports ordinary least squares (OLS), reduced-rank regression (RRR), and envelope regression (ENV) methods. Automatically selects the number of cepstral basis functions via AIC.

Usage

CepReg(
  y,
  x,
  method = c("ols", "rrr", "env"),
  number_of_K,
  if_bootstrap = FALSE,
  level = NULL,
  nboot = NULL,
  ind = NULL,
  nrank = NULL
)

Arguments

y	Numeric matrix of dimension (time points) × (samples).
x	Numeric matrix of scalar covariates with dimensions (samples) × (covariates).
method	One of "ols", "rrr", or "env" specifying the regression method.
number_of_K	Maximum number of cepstral basis functions to consider for AIC selection.
if_bootstrap	Logical; whether to compute bootstrap confidence intervals (default FALSE).
level	Confidence level for bootstrap intervals. Required if if_bootstrap = TRUE.
nboot	Integer; the number of bootstrap samples. Required if if_bootstrap = TRUE.
ind	Integer vector; indices of covariates for which the effect functions are to be estimated and plotted. Required if if_bootstrap = TRUE.
nrank	Integer; the rank used for reduced-rank regression. Required when method = "rrr" or when bootstrapping with "rrr".

Value

A list with components:

eff

A list of estimated effect functions (e.g., alpha_effect, beta_effect).

boot

A list of bootstrap results including confidence intervals; NULL if if_bootstrap = FALSE.

fit

A list containing regression coefficients, residuals, smoothed spectral estimates, and other model outputs.

Examples

set.seed(123)
niter <- 5
len <- 10
N <- 3
p <- 2
L <- floor(len/2)-1
frq <- (1:L)/len
mu <- rep(0, p)
rho <- 0
Sigma <- generate_sig(p, rho)

X <- MASS::mvrnorm(N, mu, Sigma)
X[,1] <- runif(N, 0, 1)

spec <- matrix(0,len,N)
for(j in 1:N){
  eta1 <- rnorm(1,0,0.5)
  eta2 <- rnorm(1,0,0.5)
  eta3 <- rnorm(1,0,0.5)
  spec[,j] <- exp(
    2*cos(2*pi*(1:len)/len) +
    X[j,1]*(2*cos(4*pi*(1:len)/len)) +
    eta1 + eta2*cos(2*pi*(1:len)/len) +
    eta3*(cos(4*pi*(1:len)/len))
    )
 }

Z <- data_generater(N,len,sqrt(spec))

res_ols <- CepReg(Z, X, method = "ols", number_of_K = 2,
         if_bootstrap = TRUE, level = 0.95,
         nboot = 2, ind = 1)

eff_ols <- res_ols$eff
boot_ols <- res_ols$boot

plot(frq, eff_ols$alpha_effect, type = 'l', col = "black", xlab = "Frequency", ylab = "",
     ylim = range(c(boot_ols$alpha_ci,
     eff_ols$alpha_effect, 2*cos(2*pi*frq)+0.577)))
title(ylab = expression(alpha(omega)), line = 2, cex.lab = 1.2)
lines(frq, boot_ols$alpha_ci[, 1], col = "black")
lines(frq, boot_ols$alpha_ci[, 2], col = "black")
lines(frq, 2 * cos(2 * pi * frq) + 0.577, col = "red", lty = 1, lwd = 2)
legend("topright",legend = c("True", "Estimated", "CI"),
     col = c("red", "black", "black"), ncol = 1)

CepReg documentation built on Sept. 10, 2025, 10:38 a.m.