Title:	Fit the Vector Autoregressive Model for Multiple Individuals
Version:	1.0.2
Description:	Fit the vector autoregressive model for multiple individuals using the 'OpenMx' package (Hunter, 2017 <doi:10.1080/10705511.2017.1369354>).
URL:	https://github.com/jeksterslab/fitVARMxID, https://jeksterslab.github.io/fitVARMxID/
BugReports:	https://github.com/jeksterslab/fitVARMxID/issues
License:	MIT + file LICENSE
Encoding:	UTF-8
Depends:	R (≥ 4.1.0), OpenMx (≥ 2.22.10)
Imports:	stats, simStateSpace (≥ 1.2.14)
Suggests:	knitr, rmarkdown, testthat
RoxygenNote:	7.3.3.9000
NeedsCompilation:	no
Packaged:	2026-02-22 18:31:45 UTC; root
Author:	Ivan Jacob Agaloos Pesigan [aut, cre, cph]
Maintainer:	Ivan Jacob Agaloos Pesigan <r.jeksterslab@gmail.com>
Repository:	CRAN
Date/Publication:	2026-02-27 19:50:09 UTC

fitVARMxID: Fit the Vector Autoregressive Model for Multiple Individuals

Description

Fit the vector autoregressive model for multiple individuals using the 'OpenMx' package (Hunter, 2017 doi:10.1080/10705511.2017.1369354).

Author(s)

Maintainer: Ivan Jacob Agaloos Pesigan r.jeksterslab@gmail.com (ORCID) [copyright holder]

See Also

Useful links:

• https://github.com/jeksterslab/fitVARMxID

• https://jeksterslab.github.io/fitVARMxID/

• Report bugs at https://github.com/jeksterslab/fitVARMxID/issues

Fit the First-Order Vector Autoregressive Model by ID

Description

The function fits the first-order vector autoregressive model for each unit ID.

Usage

FitVARMxID(
  data,
  observed,
  id,
  time = NULL,
  ct = FALSE,
  center = TRUE,
  mu_fixed = FALSE,
  mu_free = NULL,
  mu_values = NULL,
  mu_lbound = NULL,
  mu_ubound = NULL,
  alpha_fixed = FALSE,
  alpha_free = NULL,
  alpha_values = NULL,
  alpha_lbound = NULL,
  alpha_ubound = NULL,
  beta_fixed = FALSE,
  beta_free = NULL,
  beta_values = NULL,
  beta_lbound = NULL,
  beta_ubound = NULL,
  psi_diag = FALSE,
  psi_fixed = FALSE,
  psi_d_free = NULL,
  psi_d_values = NULL,
  psi_d_lbound = NULL,
  psi_d_ubound = NULL,
  psi_d_equal = FALSE,
  psi_l_free = NULL,
  psi_l_values = NULL,
  psi_l_lbound = NULL,
  psi_l_ubound = NULL,
  nu_fixed = TRUE,
  nu_free = NULL,
  nu_values = NULL,
  nu_lbound = NULL,
  nu_ubound = NULL,
  theta_diag = TRUE,
  theta_fixed = TRUE,
  theta_d_free = NULL,
  theta_d_values = NULL,
  theta_d_lbound = NULL,
  theta_d_ubound = NULL,
  theta_d_equal = FALSE,
  theta_l_free = NULL,
  theta_l_values = NULL,
  theta_l_lbound = NULL,
  theta_l_ubound = NULL,
  mu0_fixed = TRUE,
  mu0_func = TRUE,
  mu0_free = NULL,
  mu0_values = NULL,
  mu0_lbound = NULL,
  mu0_ubound = NULL,
  sigma0_fixed = TRUE,
  sigma0_func = TRUE,
  sigma0_diag = FALSE,
  sigma0_d_free = NULL,
  sigma0_d_values = NULL,
  sigma0_d_lbound = NULL,
  sigma0_d_ubound = NULL,
  sigma0_d_equal = FALSE,
  sigma0_l_free = NULL,
  sigma0_l_values = NULL,
  sigma0_l_lbound = NULL,
  sigma0_l_ubound = NULL,
  robust = FALSE,
  seed = NULL,
  tries_explore = 100,
  tries_local = 100,
  max_attempts = 10,
  silent = FALSE,
  ncores = NULL
)

Arguments

Details

Measurement Model

By default, the measurement model is given by

\mathbf{y}_{i, t} = \boldsymbol{\eta}_{i, t} .

However, the full measurement model can be parameterized as follows

\mathbf{y}_{i, t} = \boldsymbol{\nu}_{i} + \boldsymbol{\Lambda} \boldsymbol{\eta}_{i, t} + \boldsymbol{\varepsilon}_{i, t}, \quad \mathrm{with} \quad \boldsymbol{\varepsilon}_{i, t} \sim \mathcal{N} \left( \mathbf{0}, \boldsymbol{\Theta}_{i} \right)

where \mathbf{y}_{i, t}, \boldsymbol{\eta}_{i, t}, and \boldsymbol{\varepsilon}_{i, t} are random variables and \boldsymbol{\nu}_{i}, \boldsymbol{\Lambda}, and \boldsymbol{\Theta}_{i} are model parameters. \mathbf{y}_{i, t} represents a vector of observed random variables, \boldsymbol{\eta}_{i, t} a vector of latent random variables, and \boldsymbol{\varepsilon}_{i, t} a vector of random measurement errors, at time t and individual i. \boldsymbol{\nu}_{i}, denotes a vector of intercepts (fixed to a null vector by default), \boldsymbol{\Lambda} a matrix of factor loadings, and \boldsymbol{\Theta}_{i} the covariance matrix of \boldsymbol{\varepsilon}. In this model, \boldsymbol{\Lambda} is an identity matrix and \boldsymbol{\Theta}_{i} is a diagonal matrix.

Discrete-Time Dynamic Structure

The dynamic structure is given by

\boldsymbol{\eta}_{i, t} = \boldsymbol{\alpha}_{i} + \boldsymbol{\beta}_{i} \boldsymbol{\eta}_{i, t - 1} + \boldsymbol{\zeta}_{i, t}, \quad \mathrm{with} \quad \boldsymbol{\zeta}_{i, t} \sim \mathcal{N} \left( \mathbf{0}, \boldsymbol{\Psi}_{i} \right)

where \boldsymbol{\eta}_{i, t}, \boldsymbol{\eta}_{i, t - 1}, and \boldsymbol{\zeta}_{i, t} are random variables, and \boldsymbol{\alpha}_{i}, \boldsymbol{\beta}_{i}, and \boldsymbol{\Psi}_{i} are model parameters. Here, \boldsymbol{\eta}_{i, t} is a vector of latent variables at time t and individual i, \boldsymbol{\eta}_{i, t - 1} represents a vector of latent variables at time t - 1 and individual i, and \boldsymbol{\zeta}_{i, t} represents a vector of dynamic noise at time t and individual i. \boldsymbol{\alpha}_{i} denotes a vector of intercepts, \boldsymbol{\beta}_{i} a matrix of autoregression and cross regression coefficients, and \boldsymbol{\Psi}_{i} the covariance matrix of \boldsymbol{\zeta}_{i, t}.

If center = TRUE, the dynamic structure is parameterized as follows

\boldsymbol{\eta}_{i, t} = \boldsymbol{\mu}_{i} + \boldsymbol{\beta}_{i} \left( \boldsymbol{\eta}_{i, t - 1} - \boldsymbol{\mu}_{i} \right) + \boldsymbol{\zeta}_{i, t}

where \boldsymbol{\mu}_{i} is equilibrium level of the latent state toward which the system is pulled over time.

Continuous-Time Dynamic Structure

The continuous-time parameterization, when ct = TRUE, for the dynamic structure is given by

\mathrm{d} \boldsymbol{\eta}_{i, t} = \left( \boldsymbol{\alpha}_{i} + \boldsymbol{\beta}_{i} \boldsymbol{\eta}_{i, t - 1} \right) \mathrm{d} t + \boldsymbol{\Psi}_{i}^{\frac{1}{2}} \mathrm{d} \mathbf{W}_{i, t}

note that \mathrm{d}\boldsymbol{W} is a Wiener process or Brownian motion, which represents random fluctuations.

If center = TRUE, the dynamic structure is parameterized as follows

\mathrm{d} \boldsymbol{\eta}_{i, t} = \boldsymbol{\beta}_{i} \left( \boldsymbol{\eta}_{i, t - 1} - \boldsymbol{\mu}_{i} \right) \mathrm{d} t + \boldsymbol{\Psi}_{i}^{\frac{1}{2}} \mathrm{d} \mathbf{W}_{i, t}

Value

Returns an object of class varmxid which is a list with the following elements:

call

Function call.

args

List of function arguments.

fun

Function used ("FitVARMxID").

model

A list of generated OpenMx models.

output

A list of fitted OpenMx models.

converged

A logical vector indicating converged cases.

robust

A list of output from OpenMx::imxRobustSE() with argument details = TRUE for each id if robust = TRUE.

Author(s)

Ivan Jacob Agaloos Pesigan

References

Hunter, M. D. (2017). State space modeling in an open source, modular, structural equation modeling environment. Structural Equation Modeling: A Multidisciplinary Journal, 25(2), 307–324. doi:10.1080/10705511.2017.1369354

Neale, M. C., Hunter, M. D., Pritikin, J. N., Zahery, M., Brick, T. R., Kirkpatrick, R. M., Estabrook, R., Bates, T. C., Maes, H. H., & Boker, S. M. (2015). OpenMx 2.0: Extended structural equation and statistical modeling. Psychometrika, 81(2), 535–549. doi:10.1007/s11336-014-9435-8

See Also

Other VAR Functions: LDL(), Softplus()

Examples

# Generate data using the simStateSpace package-------------------------
library(simStateSpace)
set.seed(42)
n <- 5
time <- 100
p <- 2
alpha <- rep(x = 0, times = p)
beta <- 0.50 * diag(p)
psi <- 0.001 * diag(p)
psi_l <- t(chol(psi))
mu0 <- simStateSpace::SSMMeanEta(
  beta = beta,
  alpha = alpha
)
sigma0 <- simStateSpace::SSMCovEta(
  beta = beta,
  psi = psi
)
sigma0_l <- t(chol(sigma0))
sim <- SimSSMVARFixed(
  n = n,
  time = time,
  mu0 = mu0,
  sigma0_l = sigma0_l,
  alpha = alpha,
  beta = beta,
  psi_l = psi_l
)
data <- as.data.frame(sim)

# Fit the model---------------------------------------------------------
# center = TRUE
library(fitVARMxID)
fit <- FitVARMxID(
  data = data,
  observed = paste0("y", seq_len(p)),
  id = "id",
  center = TRUE
)
print(fit)
summary(fit)
coef(fit)
vcov(fit)
converged(fit)

# Fit the model---------------------------------------------------------
# center = FALSE
library(fitVARMxID)
fit <- FitVARMxID(
  data = data,
  observed = paste0("y", seq_len(p)),
  id = "id",
  center = FALSE
)
print(fit)
summary(fit)
coef(fit)
vcov(fit)
converged(fit)

LDL' Decomposition of a Symmetric Positive-Definite Matrix

Description

Performs an LDL' factorization of a symmetric positive-definite matrix X, such that

X = L D L^\prime,

where L is unit lower-triangular (ones on the diagonal) and D is diagonal.

Usage

LDL(x, epsilon = 1e-10)

InvLDL(s_l, uc_d)

Arguments

Details

LDL() returns both the unit lower-triangular factor L and the diagonal factor D. The strictly lower-triangular part of L is also provided for convenience. The function additionally computes an unconstrained vector uc_d such that Softplus(uc_d) = d. This uses a numerically stable inverse softplus implementation based on log(expm1(d)) (and a large-d rewrite), rather than the unstable expression \log(\exp(d) - 1).

InvLDL() returns a symmetric positive definite matrix from the strictly lower-triangular part of L and the unconstrained vector uc_d. The reconstructed matrix is symmetrized as (\Sigma + \Sigma^\prime)/2 to reduce numerical asymmetry.

Value

• LDL(): a list with components:

  • l: a unit lower-triangular matrix L

  • s_l: a strictly lower-triangular part of L

  • d: a vector of diagonal entries of D

  • uc_d: unconstrained vector with \mathrm{softplus}(uc\_d) = d

  • x: input matrix (with diagonal zeros possibly replaced by epsilon)

  • epsilon: the epsilon value used

• InvLDL(): a symmetric positive definite matrix

See Also

Other VAR Functions: FitVARMxID(), Softplus()

Examples

set.seed(123)
x <- crossprod(matrix(rnorm(16), 4, 4)) + diag(1e-6, 4)
ldl <- LDL(x = x)
ldl
inv_ldl <- InvLDL(s_l = ldl$s_l, uc_d = ldl$uc_d)
inv_ldl
max(abs(x - inv_ldl))

Softplus and Inverse Softplus Transformations

Description

The softplus transformation maps unconstrained real values to the positive real line. This is useful when parameters (e.g., variances) must be positive. The inverse softplus transformation recovers the unconstrained value from a strictly positive input.

Usage

Softplus(x)

InvSoftplus(x)

Arguments

Details

Mathematical definitions:

• Softplus(x) = log(1 + exp(x))

• InvSoftplus(x) = log(exp(x) - 1)

Numerical implementation (stable forms):

• Softplus(x) = max(x, 0) + log1p(exp(-abs(x)))

• InvSoftplus(y) uses log(expm1(y)), and for large y uses the rewrite y + log1p(-exp(-y)) for stability.

For numerical stability, these functions use log1p() and expm1() internally. InvSoftplus() requires strictly positive input and will error if any values are ⁠<= 0⁠.

Value

• Softplus(): numeric vector or matrix of nonnegative values (mathematically strictly positive for finite inputs, but can underflow to 0 for very negative values).

• InvSoftplus(): numeric vector or matrix of unconstrained values.

Author(s)

Ivan Jacob Agaloos Pesigan

See Also

Other VAR Functions: FitVARMxID(), LDL()

Examples

# Apply softplus to unconstrained values
x <- c(-5, 0, 5)
y <- Softplus(x)

# Recover unconstrained values
x_recovered <- InvSoftplus(y)

y
x_recovered

Parameter Estimates

Description

Parameter Estimates

Usage

## S3 method for class 'varmxid'
coef(
  object,
  mu = TRUE,
  alpha = TRUE,
  beta = TRUE,
  nu = TRUE,
  psi = TRUE,
  theta = TRUE,
  ncores = NULL,
  ...
)

Arguments

Value

Returns a list of vectors of parameter estimates.

Author(s)

Ivan Jacob Agaloos Pesigan

Check Model Convergence

Description

Determines whether each fitted OpenMx model in a varmxid object meets convergence criteria based on (a) acceptable optimizer status and gradient size, (b) a positive-definite Hessian, and (c) parameters not being at their bounds.

Usage

converged(object, ...)

## S3 method for class 'varmxid'
converged(object, prop = FALSE, ...)

Arguments

Value

For the varmxid method: If prop = FALSE, a named logical vector indicating, for each individual fit, whether the convergence criteria are met. If prop = TRUE, the proportion of cases that converged.

Author(s)

Ivan Jacob Agaloos Pesigan

Print Method for Object of Class varmxid

Description

Print Method for Object of Class varmxid

Usage

## S3 method for class 'varmxid'
print(
  x,
  means = FALSE,
  mu = TRUE,
  alpha = TRUE,
  beta = TRUE,
  nu = TRUE,
  psi = TRUE,
  theta = TRUE,
  digits = 4,
  ...
)

Arguments

Value

Prints means or raw estimates depending the the value of the argument means.

Author(s)

Ivan Jacob Agaloos Pesigan

Summary Method for Object of Class varmxid

Description

Summary Method for Object of Class varmxid

Usage

## S3 method for class 'varmxid'
summary(
  object,
  means = FALSE,
  mu = TRUE,
  alpha = TRUE,
  beta = TRUE,
  nu = TRUE,
  psi = TRUE,
  theta = TRUE,
  digits = 4,
  ncores = NULL,
  ...
)

Arguments

Value

Returns means or raw estimates depending the the value of the argument means.

Author(s)

Ivan Jacob Agaloos Pesigan

Sampling Covariance Matrix of the Parameter Estimates

Description

Sampling Covariance Matrix of the Parameter Estimates

Usage

## S3 method for class 'varmxid'
vcov(
  object,
  mu = TRUE,
  alpha = TRUE,
  beta = TRUE,
  nu = TRUE,
  psi = TRUE,
  theta = TRUE,
  robust = FALSE,
  ncores = NULL,
  ...
)

Arguments

Value

Returns a list of sampling variance-covariance matrices.

Author(s)

Ivan Jacob Agaloos Pesigan