A structural equation formulation for general quasi-periodic Gaussian processes

Problem Statement

Quasi-periodic signals—common in fields like robotics, physiology, astronomy, and communication—are characterized by regular patterns disrupted by noise or unmodeled disturbances. While standard Periodic Gaussian Processes (PGPs) effectively model correlations within a single period, they fail to capture dependencies between different periods. Existing solutions, such as the Quasi-Periodic Gaussian Process (QPGP) introduced by [13], address this by defining a covariance function as the product of a geometrically decaying kernel and MacKay's periodic kernel. However, these approaches have significant limitations:

1. Kernel Rigidity: They are often restricted to specific kernel forms (e.g., MacKay's kernel), limiting flexibility in modeling within-period variations.
2. Computational Cost: Evaluating the likelihood and making predictions for these processes typically involves inverting large covariance matrices of order $n$, resulting in a computational complexity of $O(n^3)$ or $O(k^2p^2)$ (where $n$ is the total sample size, $p$ is the period, and $k$ is the number of periods). This creates a bottleneck for long-duration time series.
3. Uncertainty Quantification: Efficient methods for estimating standard errors and confidence intervals for parameters in these complex models are often lacking or computationally prohibitive.

Methodology

The authors propose a new structural equation formulation based on a Gaussian linear state-space model to define a broad family of QPGPs.

1. Structural Equation Model (Dynamical System):
Instead of defining the process solely through a covariance function, the authors define the QPGP $\{Y_t\}$ via a vector autoregressive recursion of order 1. The data is partitioned into $p$-dimensional periodic blocks $Y_i$. The evolution of these blocks is governed by:
$$Y_{i+1} = \omega Y_i + Z_{i+1}$$
where:

• $\omega \in (-1, 1)$ is the between-period correlation parameter.
• $Z_{i+1}$ are independent and identically distributed (i.i.d.) $p$-dimensional Gaussian vectors with a zero mean and a covariance matrix $K$ representing the within-period correlation.
• $K$ can be any valid periodic covariance kernel (e.g., MacKay, Matérn, Cosine), allowing for flexible modeling of intra-period dynamics.

This formulation decouples the between-period correlation (controlled by $\omega$) from the within-period structure (controlled by $K$).

2. Estimation Strategy:
The authors develop a two-stage fast estimation algorithm to estimate $\omega$ and the periodic covariance kernel $\kappa_p$ (represented by matrix $K$):

• Stage I (Reduced Likelihood): They utilize a reduced negative log-likelihood function that ignores the marginal contribution of the initial block (to gain flexibility for general QPGPs). By taking derivatives with respect to $\omega$ and $K$, they derive explicit update equations. An iterative alternating minimization scheme is used to solve for $\tilde{\omega}$ and $\tilde{K}$.
• Stage II (Projection): Since the unconstrained $\tilde{K}$ from Stage I may not strictly satisfy the properties of a periodic covariance kernel, the algorithm projects $\tilde{K}$ onto the space of valid periodic kernels. This is achieved by minimizing the Frobenius norm and ensuring positive definiteness via spectral truncation (Fourier transform techniques).
• Parametric Extension: If the kernel belongs to a known parametric family (e.g., Matérn with parameters $\theta, \sigma^2$), Stage II is simplified to a grid search or optimization over these hyperparameters.

3. Prediction and Uncertainty:

• Prediction: The structural equations allow for a computationally efficient formula for the best linear predictor. The complexity is reduced from $O(t^2)$ to $O(p^2)$ per step, as it only requires inverting matrices of size $p$ rather than $t$.
• Bootstrap: The rapid generation of sample paths via the structural equations enables a model-based bootstrap procedure. Residuals from the fitted model are resampled to generate new QPGP trajectories, allowing for the estimation of standard errors and confidence intervals for $\omega$ and $\kappa_p$.

Key Contributions

1. General Structural Formulation: A new dynamical equation framework that supports a broad family of QPGPs with arbitrary periodic covariance kernels, overcoming the rigidity of previous MacKay-based approaches.
2. Computational Efficiency: A significant reduction in computational complexity. Likelihood evaluations and predictions scale with $O(p^2)$ rather than $O(k^2p^2)$ or $O(n^3)$, making the method scalable for large datasets.
3. Consistent Estimation: The paper establishes the consistency of the proposed estimators ($\hat{\omega}$ and $\hat{\kappa}_p$) as the sample size increases.
4. Uncertainty Quantification: A principled bootstrap approach is provided to quantify estimation uncertainty, including standard errors and confidence intervals, which is computationally feasible due to the efficient structural formulation.

Results

The authors validate their methodology through simulation studies and three real-world case studies:

• Simulation:

  • Speed: The proposed method is orders of magnitude faster than exact likelihood evaluation (Eq. 11) and the benchmark algorithm from [13]. For $n=10,000$, likelihood evaluation using the reduced likelihood (Eq. 14) took ~0.48ms, compared to ~17.20ms for the exact likelihood (Eq. 11) and ~0.62ms for the standard reduced likelihood (Eq. 13).
  • Accuracy: The Root Mean Squared Error (RMSE) of the proposed estimator is comparable to Maximum Likelihood Estimation (MLE) via grid search, despite being significantly faster (milliseconds vs. thousands of milliseconds per run).
  • Bootstrap Performance: Bootstrap standard errors closely approximate the true standard errors, and 95% confidence intervals achieve coverage probabilities near the nominal 95% level.

• Case Studies:

  1. CO2 Emission Data: Fitted to detrended monthly data ($n=612$). The model identified a period of 12 months. The general kernel provided the best fit (lowest RMSE), with $\omega \approx 0.975$. The benchmark method [13] required 632.06ms for estimation on this dataset, whereas the proposed method required 20.24ms.
  2. Sunspot Numbers: Fitted to yearly data ($n=322$). The model identified a period of 11 years. The Periodic Matérn kernel performed best, with statistically significant parameter estimates.
  3. Water Level Tidal Data: Fitted to high-frequency data ($n=14,400$). The model identified a period of 148 (approx. 24.6 hours). The method successfully handled the large sample size, demonstrating scalability, with the Matérn kernel yielding the lowest RMSE.

In all cases, the proposed method outperformed the benchmark [13] in terms of computational time (estimation and prediction) while maintaining or improving predictive accuracy (RMSE).

Significance and Claims

The paper claims that its structural equation formulation provides a flexible and scalable framework for modeling quasi-periodic signals. By shifting from a direct covariance matrix inversion to a state-space recursion, the authors achieve a "significant reduction in complexity" that enables the analysis of large-scale and streaming time series data.

The authors emphasize that their approach offers a trade-off: it sacrifices a small, quantifiable amount of estimation accuracy (by ignoring the initial block's marginal contribution and using a reduced likelihood) to gain substantial computational efficiency. This efficiency, in turn, enables the practical application of bootstrap methods for uncertainty quantification, which is often infeasible with standard QPGP methods on large datasets. The work is presented as a generalization that allows for richer periodic dependence modeling and is suitable for applications requiring both speed and statistical rigor.