Source-linked AI summary
A Tutorial on Estimating Time-Varying Vector Autoregressive Models
Jonas M B Haslbeck, Laura F Bringmann, Lourens J Waldorp
TL;DR
Psychological time series support within-person modeling, but standard VAR models assume constant parameters even when psychological dynamics may change over time. The paper compares spline-based GAM and kernel-smoothing approaches, with and without regularization, and finds that method choice depends mainly on sample size and sparsity. It also provides practical recommendations and tutorials for applying these methods.
Problem
Standard VAR models assume constant parameters, although time-varying parameters are plausible and potentially important in psychological research.
Method
The paper evaluates GAM spline and kernel-smoothing methods, with and without regularization or significance testing, using simulations and applied tutorials.
Results
Method performance depends on sample size and sparsity: KS methods generally perform better for sparse models and moderate samples, whereas GAM methods perform better when many parameters are nonzero and samples are large.
Takeaways & Limitations
Researchers can select among the evaluated methods using their expected VAR sparsity and available sample size.
Takeaways & Limitations
Kernel-smoothing estimates have lower sensitivity at the beginning and end of the time series because the kernel is truncated.
Abstract
from arXiv · showhide
Time series of individual subjects have become a common data type in psychological research. These data allow one to estimate models of within-subject dynamics, and thereby avoid the notorious problem of making within-subjects inferences from between-subjects data, and naturally address heterogeneity between subjects. A popular model for these data is the Vector Autoregressive (VAR) model, in which each variable is predicted as a linear function of all variables at previous time points. A key assumption of this model is that its parameters are constant (or stationary) across time. However, in many areas of psychological research time-varying parameters are plausible or even the subject of study. In this tutorial paper, we introduce methods to estimate time-varying VAR models based on splines and kernel-smoothing with/without regularization. We use simulations to evaluate the relative performance of all methods in scenarios typical in applied research, and discuss their strengths and weaknesses. Finally, we provide a step-by-step tutorial showing how to apply the discussed methods to an openly available time series of mood-related measurements.
1 Introduction
Intensive longitudinal data support within-person modeling, but standard VAR models assume parameters remain constant across time. This tutorial introduces spline- and kernel-based time-varying VAR methods, evaluates them in simulations, and provides practical tutorials and recommendations.
- Intensive longitudinal data enable models of within-subject dynamics while addressing heterogeneity between subjects.
- Stationary VAR models can be inappropriate when psychological parameters change over time.The paper illustrates this issue with time-varying effects among depressed mood, anxiety, and worrying.
- Time-varying models can represent changing symptom interactions and support the study of transitions between healthy and unhealthy states.
- The tutorial presents GAM spline methods and penalized kernel-smoothing methods for estimating time-varying VAR parameters.
- The paper evaluates method performance in ESM-like situations, discusses strengths and weaknesses, compares stationary counterparts, and provides R-package tutorials.
2 Estimating Time-Varying VAR Models
The paper defines time-varying VAR models and estimates their parameters using either spline-based GAMs or locally weighted kernel-smoothing. Both approaches rely on smooth parameter functions, while kernel-smoothing additionally trades sensitivity against reliability through bandwidth choice.
- 2.1 Vector Autoregressive (VAR) Model: A VAR(1) predicts each variable at time t from all variables at time t − 1, including autoregressive and cross-lagged effects.
- 2.1 Vector Autoregressive (VAR) Model: Time-varying VAR models assign parameter matrices that change across time while assuming the parameter functions are smooth deterministic functions.Smoothness means the first derivative exists everywhere.
- 2.2 The GAM Method: GAM estimates each time-varying parameter as a function of time using weighted spline basis functions combined through regression.The paper uses thin plate splines and penalized likelihood to control flexibility and reduce overfitting.
- 2.3 The Kernel-smoothing Method: Both methods assume local stationarity, meaning nearby parameter values are similar over time.Kernel-smoothing additionally requires a bandwidth small enough to capture time variation.
- 2.3 The Kernel-smoothing Method: Kernel-smoothing estimates a sequence of local models by weighting observations according to their temporal distance from each target time point.Smaller bandwidths increase sensitivity to parameter changes but use less data and produce less reliable estimates.
- 2.3 The Kernel-smoothing Method: Kernel estimates are less sensitive to effects at the beginning and end of a time series because the kernel is truncated there.
3 Evaluating Performance via Simulation
Two simulations evaluate time-varying VAR estimators across parameter shapes, sample sizes, and sparsity-relevant structure recovery. Performance depends on whether parameters are constant or time-varying, sample size, regularization, and estimator family.
- Simulation design: The first simulation generates time-varying VAR models from random graphs with fixed sparsity to compare methods across parameter functions and sample sizes.It evaluates the methods’ strengths and weaknesses in scenarios typical of applied research.
- Absolute estimation error: For constant nonzero parameters, stationary unregularized GLM has the lowest error, while all methods improve as n increases.Stationary methods are appropriate because the true parameter is constant; regularization reduces performance for nonzero parameters.
- Absolute estimation error: For linearly varying parameters, regularized methods perform best at very small n, unregularized methods dominate at intermediate n, and time-varying methods outperform stationary methods after n exceeds about 155.KS performs best among time-varying methods up to roughly n = 530, after which time-varying methods perform similarly.
- Absolute estimation error: For step functions, time-varying methods outperform stationary methods between n = 46 and n = 69, while regularized KS and thresholded GAM benefit from recovering zero-valued intervals.The stronger time variation makes the best constant approximation less accurate, and exact zeros favor thresholding or regularization.
- Absolute estimation error: For constant zero parameters, regularized GLM and KS have the lowest error, whereas unregularized GAM performs poorly because of high variance.Regularization and thresholding help when the true parameter is zero, while unregularized methods are favored for nonzero parameters.
- Structure recovery and computation: In structure recovery, unthresholded methods have sensitivity 1 and precision 0 because they always return nonzero estimates, making them unsuitable for identifying nonzero structure.For the remaining methods, sensitivity tends toward 1 as n increases; computation time is not a major concern in typical psychological applications.
3.2 Simulation B: Varying Sparsity
Simulation B varied local sparsity through indegree in an upper-triangular time-varying VAR structure. Regularization helped at low indegree, whereas GAM methods were more robust than kernel-smoothing methods as indegree increased.
- Data generation: The simulation varied local sparsity by changing each node’s indegree from 1 to 20 in an upper-triangular VAR structure.The row sums correspond to node indegrees, and the 20-node graph preserved the illustrated six-variable structure.
- Results: All methods became less accurate as indegree increased because more parameters were nonzero and predictors were more correlated.Figure 9 reports mean absolute error across indegrees 1, 10, and 20, both averaged across parameter types and separately by function.
- Results: Regularized methods performed best at low indegree but worst at high indegree because shrinkage toward zero is beneficial when most parameter functions are zero.When most functions are nonzero, the same zero bias increases estimation error.
- Results: GAM methods showed a smaller performance decline than kernel-smoothing methods as indegree increased.The paper’s overall discussion likewise reports that GAM methods outperform kernel-smoothing methods when indegree is high.
- Overall discussion: For practical choice, kernel-smoothing methods are recommended for sparse models with moderate samples, whereas GAM methods are likely better for largely nonzero models with large samples.The recommendation specifies moderate samples as n < 200–300 and favors regularization when few parameter functions are nonzero.
- Limitations: The simulations varied signal-to-noise conditions, temporal variation, and model size, so performance depends on these settings and on effective sample size.The studies used p = 10 variables and one lag; planned and unplanned missingness can make effective sample size lower than the number of measurements.
4 Estimating time-varying VAR model on Mood Time Series
The tutorial applies the KS(L1) method to estimate a time-varying VAR model for mood measurements, while addressing missingness, preprocessing, tuning, reliability, visualization, and model selection. The example uses 12 mood variables and yields 20 estimation points from 876 usable VAR-design rows.
- Dataset and preprocessing: The tutorial fits a time-varying VAR model to 12 mood variables measured over 1476 time points from an individual diagnosed with major depression.Measurements were collected at 10 pseudo-randomized intervals per day across 238 consecutive days during a medication dose-reduction period.
- Model specification: The analysis uses 20 equally spaced estimation points, bandwidth b = 0.34 selected by out-of-bag cross-validation, and a first-order ℓ1-regularized VAR with cross-validated penalties.Variables are standardized before fitting because ℓ1 regularization otherwise depends on predictor variance.
- Model specification: The beep and day indicators preserve measurement succession so estimated effects represent the intended previous measurement rather than overnight or other longer lags.This is necessary because the dataset contains both random missingness and a designed day-night break.
- Dataset and preprocessing: Missing measurements and day-night breaks reduce the VAR design matrix to 876 usable rows from 1475 possible lagged observations.The first observation lacks a previous time point, and additional observations are excluded when successive measurements are unavailable.
- Reliability and interpretation: The tutorial assesses estimate reliability with bootstrapped confidence intervals for KS and credible intervals for GAM, while noting that ℓ1 shrinkage biases estimates toward zero.The resulting model is a time-varying mixed VAR model with 12 nodes and 20 estimation points.
- Interpretation and model selection: Time-varying parameter visualizations can reveal changes after medication reduction and support analyses of intercepts, individual interactions, or network density.The authors suggest that the largest changes after the second dose reduction could be related to that reduction, while treating model selection separately.
5 Discussion
The discussion compares GAM and kernel-smoothing methods with and without regularization, links their performance to sample size and sparsity, and provides reproducible software tutorials. It also emphasizes that small samples limit recovery of strong time dependence and identifies Bayesian and broader model extensions as future work.
- Discussion: The simulations compare GAM and kernel-smoothing methods with and without regularization across sample sizes and assumptions about VAR sparsity.The comparison is intended to help researchers select among the considered methods for situations typical of psychological applications.
- Discussion: The paper provides step-by-step tutorials for kernel-smoothing and GAM implementations using the mgm and tvvarGAM R packages.The tutorials are intended to support reproducible application of both method families.
- Discussion: For the studied time-varying functions, the best time-varying method outperformed stationary methods when n > 46.The authors suggest this threshold indicates applicability to typical experience-sampling-method data.
- Discussion: With small samples, time-varying methods produce estimates similar to stationary methods and cannot recover most strong dependence of parameters on time.Thus, apparent stationarity in small samples need not reflect the underlying parameter function.
- Future work: Future work could compare Bayesian approaches with the presented frequentist methods and extend the methods beyond standard VAR models.Examples include mixed VAR models and unified structural equation models.
A Sampling Variation around Aggregated Absolute Errors
This section examines sampling variation around mean absolute estimation errors across methods and sample sizes. It shows that method differences can be difficult to predict for an individual sample and that regularization changes the shape of sampling quantiles.
- Aggregated absolute errors: Mean absolute errors summarize each method’s expected error across time points and iterations, while population quantiles describe sampling variation around those means.The figure compares six stationary and time-varying estimation methods as a function of observations n.
- Aggregated absolute errors: For n = 103, GAM has a larger population mean error than KS in the example, but their sampling distributions largely overlap.The overlap indicates that the difference between methods has substantial variance for a specific sample.
- Aggregated absolute errors: For n = 103, it is difficult to predict for a specific sample whether GAM or KS will have the larger error.This follows from the large variance of the difference between their errors despite the population-level mean ordering.
- Aggregated absolute errors: Unregularized methods have wide confidence intervals at small n that narrow as n increases, whereas regularized quantiles first increase and then decrease.At small n, regularization sets most estimates to zero, producing initially narrow quantiles.
B Sampling Variation around Absolute Errors over Time
This section examines sampling variation in parameter estimates for constant and linearly increasing functions across sample sizes and estimation methods. Sampling variance is small at n = 103 and approaches zero as sample size grows.
- Sampling variation over time: Figure 12 reports mean estimates and 10% and 90% quantiles for constant and linearly increasing parameters across five methods and three sample sizes.The rows correspond to n = 103, 530, and 1803, while columns distinguish constant from linearly increasing parameters.
- Sampling variation over time: Sampling variance is small for n = 103 and approaches zero as n becomes large.The figure summarizes this pattern separately for constant and linearly increasing parameter functions.
C Computational Cost
KS(L1) has nearly constant fitting cost as observations increase, whereas GAM(st) and KS(L1) bandwidth selection grow roughly linearly. In the current 10-variable setting, bandwidth selection makes KS(L1) more expensive overall.
- Computational Cost: KS(L1) fitting has computational complexity O(p|E|p log p|L|q), where p is the number of variables, |E| the estimation points, and |L| the included lags.The supplied complexity expression is reproduced in plain-text notation.
- Computational Cost: Bandwidth selection has complexity O(p|F||Fs|p log p|L|q), while standard GAM fitting has complexity O(pnq2q), with q increasing as basis functions increase.Here, |F| is the number of folds, |Fs| the leave-out points per fold, n the modeled time points, and q the total coefficients.
- Computational Cost: KS(L1) computational cost hardly increases with the number of observations n, while GAM(st) cost increases roughly linearly.The comparison concerns fitting time-varying VAR models as n grows.
- Computational Cost: KS(L1) bandwidth selection also increases roughly linearly with n, making total KS(L1) cost larger than GAM(st) for p = 10 variables.This comparison includes the data-driven bandwidth-selection step required by KS(L1).
D Code to select Appropriate Bandwidth in KS(L1) Method
The KS(L1) bandwidth is selected by fitting models across candidate values and choosing the value with the smallest average out-of-sample prediction error. In this example, the selected bandwidth is 0.34.
- D Code to select Appropriate Bandwidth in KS(L1) Method: bwSelect() fits KS(L1) models across training sets and selects the bandwidth minimizing prediction error across variables and hold-out sets.The procedure evaluates candidate bandwidths using out-of-sample prediction error.
- D Code to select Appropriate Bandwidth in KS(L1) Method: The example evaluates ten equally spaced candidate bandwidths between 0.01 and 1 using specified training/test splits.The split count, test-set size, and candidate sequence are passed through bwFolds, bwFoldsize, and bwSeq.
- D Code to select Appropriate Bandwidth in KS(L1) Method: The bandwidth 0.34 minimized the average out-of-sample prediction error in the fitted bandwidth path.The full bandwidth path is shown in Figure 14.
- D Code to select Appropriate Bandwidth in KS(L1) Method: Prediction error was large at bandwidth 0.01, decreased toward 0.34, and increased afterward.If an endpoint bandwidth minimizes error, the search should be repeated beyond the current range.
E Estimating time-varying VAR model via GAM(st)
This section demonstrates estimating a time-varying VAR model with the GAM(st) method using the reproducible tvvarGAM R package and executable code.
- E Estimating time-varying VAR model via GAM(st): The tutorial estimates a time-varying VAR model with the GAM(st) method using the R package tvvarGAM.The shown code is described as fully reproducible and is also available in an R file on GitHub.
- E Estimating time-varying VAR model via GAM(st): Readers can execute the reproducible code while reading, and the code is available in an R file on GitHub.The passage identifies reproducibility as a feature of the tutorial demonstration.
E.1 Load R-packages and dataset
The tutorial loads the mood dataset, selects its 12 mood-related variables, and prepares the required R packages for the GAM(st) analysis.
- E.1 Load R-packages and dataset: The dataset is loaded from the mgm package and subset to 12 mood-related variables.The selected variables and their labels are stored for analysis.
- E.1 Load R-packages and dataset: The analysis uses mgm version 1.2-8 and tvvarGAM version 0.1.1.The tvvarGAM package is installed from GitHub in the shown setup.
- E.1 Load R-packages and dataset: The mood matrix receives the selected mood labels as column names, while time information is loaded from the dataset's time component.The code assigns the variable labels and extracts time_data.
E.2 Estimating time-varying VAR model
The tvvarGAM() function estimates a time-varying VAR model from intensive longitudinal data using measurement-order information and a specified number of basis functions. Its output includes parameter estimates and 5% and 95% confidence intervals at each time point.
- tvvarGAM() estimates the time-varying VAR model using the data, notification numbers, day numbers, and desired number of basis functions.The notification and day variables are used to construct the VAR design matrix.
- The function returns parameter estimates in a p × (p + 1) × timepoints array, including a first row of estimated intercepts.
- The output also contains 5% and 95% confidence intervals with the same dimensions as the estimate array.The tvvarGAM package therefore provides interval estimates without a separate resampling scheme.
E.3 Visualize time-varying VAR model
Figure 15 visualizes thresholded time-varying VAR estimates at three estimation points and plots selected parameters over time. The GAM(ks) results show fewer significant edges than KS(L1), while some parameter trajectories differ across methods.
- Figure 15 shows thresholded VAR models at estimation points 2, 10, and 18, displaying only significant arrows with p-value < 0.05.Blue arrows indicate positive relationships, red arrows negative relationships, and arrow width represents the absolute parameter value.
- The figure’s bottom row plots three parameters over time using unthresholded point estimates and shaded 5% and 95% credible intervals.Self-loops in the top row represent autocorrelations.
- The GAM(ks) method produces fewer edges than KS(L1), indicating more conservative estimates.
- The effect of Relaxed on itself tends to decrease over the measured time interval, consistent with the KS(L1) results.
- GAM(ks) estimates an increase in Guilty’s effect on Satisfied in the second half, whereas KS(L1) estimates a decrease.The paper states that the GAM(st) second-half estimates appear incorrect because stationary VAR estimates show a clearly negative effect there.