ARIMAXff: VGLTSMs Family functions: The Order-(p, d, q) Autoregressive...

View source: R/TSfamilyARIMAX.R

ARIMAXffR Documentation

VGLTSMs Family functions: The Order–(p, d, q) Autoregressive Integrated Moving Average Model (ARIMA(p, d, q)) with covariates


Maximum likelihood estimation of the drift, standard deviation or variance of the random noise, and coefficients of an autoregressive integrated moving average process of order-(p, d, q) with covariates by MLE using Fisher scoring. No seasonal terms handled yet. No seasonal components handled yet.


      ARIMAXff(order   = c(1, 1, 0),
               zero     = c("ARcoeff", "MAcoeff"),
               diffCovs = TRUE,
               xLag     = 0,
               include.current = FALSE,
               type.EIM = c("exact", "approximate")[1], 
               var.arg  = TRUE,
               nodrift  = FALSE,
               noChecks = FALSE,
               ldrift   = "identitylink", 
               lsd      = "loglink",
               lvar     = "loglink",
               lARcoeff = "identitylink",
               lMAcoeff = "identitylink", 
               idrift   = NULL,
               isd      = NULL,
               ivar     = NULL,
               iARcoeff = NULL, 
               iMAcoeff = NULL) 



Integer vector with three components, (p, d, q): The AR order (p), the degree of differencing (d), and the MA order (q).


Integer or character–strings vector. Name(s) or position(s) of the parameters/linear predictors to be modeled as intercept-only. Details at zero.


Logical. The default is diffCovs = TRUE, which means that the order–d differences of the covariates (if entered) are internally computed and then incorporated in the conditional–mean model. Otherwise, only the current actual values of the covariates are included.


Integer, non–negative. If xLag > 0, the covariates at current time, \boldsymbol{x}_t, plus the lagged covariates up to order 'xLag' are embedded into the design matrix (as covariates too). Leave xLag = 0 and only the current value, \boldsymbol{x}_t, will be considered. See more details below.


Logical. Same as ARIMAX.errors.ff.


The type of expected information matrix (EIM) of the ARMA process to be utilized in Fisher scoring. type.EIM = "exact" (default) enables the exact IM (Porat, 1986), otherwise the approximate version is utilized.


Logical. If FALSE (default), then the standard deviation of the random noise is estimated. Else, the variance estimate is returned.


Logical. nodrift = TRUE supresses estimation of the intercept (the drift in the ARMA case), which is set to zero internally.


Logical. If FALSE (default), this family function internally checks stationarity (AR case) and invertibility (MA case) of the the estimated model. A warning is correspondingly displayed.

ldrift, lsd, lvar, lARcoeff, lMAcoeff

Link functions applied to the intercept, the random noise standard deviation (or optionally, the variance), and the coefficients in the ARMA–type conditional–mean model.

idrift, isd, ivar, iARcoeff, iMAcoeff

Optional initial values for the intercept (drift), noise SD (or variance), and ARMA coeffcients (a vector of length p + q). If failure to converge occurs then try different values and monitor convergence by using trace = TRUE in the vglm() call.


Let \boldsymbol{x}_t be a (probably time–varying) vector of suitable covariates. The ARIMAX model handled by ARIMAXff is

\nabla^d Y_t = \mu^{\star} + \boldsymbol{\beta}^T \nabla^d \boldsymbol{x}_t + \theta_1 \nabla^d Y_{t - 1} + \cdots + \theta_p \nabla^d Y_{t - p} + \phi_1 \varepsilon_{t - 1} + \cdots + \phi_q \varepsilon_{t - q} + \varepsilon,

with \nabla^d (\cdot) the operator differencing a series d times. If diffCovs = TRUE, this function differencing the covariates d times too.

Similarly, ARMAXff manages

\nabla^d Y_t = \mu^{\star} + \boldsymbol{\beta}^T \boldsymbol{x}_t + \theta_1 Y_{t - 1} + \cdots + \theta_p Y_{t - p} + \phi_1 \varepsilon_{t - 1} + \cdots + \phi_q \varepsilon_{t - q} + \varepsilon,


\varepsilon_{t | \Phi_{t - 1}} \sim N(0, \sigma_{\varepsilon_t | \Phi_{t - 1}}^2).

Note, \sigma_{\varepsilon | \Phi_{t - 1}}^2 is conditional on \Phi_{t - 1}, the information of the joint process \left(Y_{t - 1}, \boldsymbol{x}_t \right), at time t, and hence may be modelled in terms of \boldsymbol{x}_t, if required.

ARIMAXff() and ARMAXff() handle multiple responses, thus a matrix can be used as the response. Note, no seasonal terms handled. This feature is to be incorporated shortly.

The default linear predictor is

\boldsymbol{\eta} = \left( \mu, \log \sigma^{2}_{\varepsilon_{t | \Phi_{t - 1}}}, \theta_1, \ldots, \theta_p, \phi_1, \ldots, \phi_q \right)^T.

Other links are also handled. See Links.

Further choices for the random noise, besides Gaussian, will be implemented over time.

As with ARXff and MAXff, choices for the EIMs are "exact" and "approximate". Covariates may be incorporated in the fit for any linear predictor above. Hence, ARIMAXff supports non–stationary processes (\sigma_{\varepsilon_t | \Phi_{t - 1}}^2) may depend on \boldsymbol{X}_t. Also, constraint matrices on the linear predictors may be entered through cm.ARMA or using the argument constraints, from vglm.

Checks on stationarity and invertibility on the esitmated process are performed by default. Set noChecks = TRUE to dismiss this step.


An object of class "vglmff" (see vglmff-class) to be used by VGLM/VGAM modelling functions, e.g., vglm or vgam.


zero can be a numeric or a character–strings vector specifying the position(s) or the name(s) of the parameter(s) modeled as intercept–only. Numeric values can be set as usual (See CommonVGAMffArguments). If names are entered, the parameter names in this family function are:

c("drift.mean", "noiseVar" || "noiseSD", "ARcoeff", "MAcoeff").

Manually modify this if required. For simplicity, the second choice is recommended.


No seasonal components handled yet.

If no covariates, \boldsymbol{x}_t, are incorporated in the analysis, then ARIMAXff fits an ordinary ARIMA model. Ditto with ARMAXff.

If nodrift = TRUE, then the 'drift' is removed from the vector of parameters and is not estimated.

By default, an ARMA model of order–c(1, 0) with order–1 differences is fitted. When initial values are entered (isd, iARcoeff, etc.), they are recycled according to the number of responses.

Also, the ARMA coefficients are intercept–only (note, zero = c("ARcoeff", "MAcoeff")) This may altered via zero, or by constraint matrices (See constraints) using cm.ARMA.

Checks on stationarity and/or invertibility can be manually via checkTS.VGAMextra.


Victor Miranda and T. W. Yee


Miranda, V. and Yee, T.W. (2018) Vector Generalized Linear Time Series Models. In preparation.

Porat, B., and Friedlander, B. (1986) Computation of the Exact Information Matrix of Gaussian Time Series with Stationary Random Components. IEEE Transactions on Acoustics, Speech and Signal Processing. ASSp-34(1), 118–130.

See Also

ARXff, MAXff, checkTS.VGAMextra, cm.ARMA, CommonVGAMffArguments, constraints, vglm.


nn      <- 90
theta   <- c(0.12, 0.17)  # AR coefficients
phi     <- c(-0.15, 0.20)  # MA coefficients.
sdWNN <- exp(1.0)          # SDs
mu    <- c(1.25, 0.85)     # Mean (not drift) of the process.
covX  <- runif(nn + 1)     # A single covariate.
mux3  <- mu[1] + covX
## Simulate ARMA processes. Here, the drift for 'tsd3' depends on covX.
tsdata <- data.frame(TS1 = mu[1] + arima.sim(model = list(ar = theta, ma = phi,
                             order = c(2, 1, 2)), n = nn, sd = sdWNN ),
                      TS2 = mu[2] + arima.sim(model = list(ar = theta, ma = phi,
                            order = c(2, 1, 2)), n = nn, sd = exp(2 + covX)),
                      TS3 =  mux3 + arima.sim(model = list(ar = theta, ma = phi,
                            order = c(2, 1, 2)), n = nn, sd = exp(2 + covX) ),
                      x2 = covX)

### EXAMPLE 1. Fitting a simple ARIMA(2, 1, 2) using vglm(). 
# Note that no covariates involved. 
fit.ARIMA1 <-  vglm(TS1 ~ 1, ARIMAXff(order = c(2, 1, 2), var.arg = FALSE,
                       # OPTIONAL INITIAL VALUES
                       # idrift = c(1.5)*(1 - sum(theta)), 
                       # ivar = exp(4), isd = exp(2),
                       # iARcoeff = c(0.20, -0.3, 0.1),
                       # iMAcoeff = c(0.25, 0.35, 0.1),
                        type.EIM = "exact"), 
                 data = tsdata, trace = TRUE, crit = "log")
vcov(fit.ARIMA1, untransform = TRUE)
# Fitting same model using arima().
( fitArima  <- arima(tsdata$TS1, order = c(2, 1, 2)) ) 

### EXAMPLE 2. Here only the ARMA coefficients and drift are intercept-only.
# The random noise variance is not constant.
fit.ARIMA2 <-  vglm(TS2 ~ x2,  ARIMAXff(order = c(2, 1, 2), var.arg = TRUE,
                     lARcoeff = "rhobitlink", lMAcoeff = "identitylink",
                     type.EIM = c("exact", "approximate")[1], 
                     # NOTE THE ZERO ARGUMENT. 
                     zero = c("drift.mean", "ARcoeff", "MAcoeff")),
              data = tsdata, trace = TRUE)

coef(fit.ARIMA2, matrix = TRUE)

### EXAMPLE 3. Here only ARMA coefficients are intercept-only.
# The random noise variance is not constant.
# Note that the "drift" and the "variance" are "generated" in 
# terms of 'x2' above for TS3.

fit.ARIMA3 <- vglm(TS3 ~ x2,  ARIMAXff(order = c(1, 1, 2), var.arg = TRUE,
                     lARcoeff = "identitylink", lMAcoeff = "identitylink",
                     type.EIM = c("exact", "approximate")[1], nodrift = FALSE,
                     zero = c( "ARcoeff", "MAcoeff")), # NOTE THE ZERO ARGUMENT. 
              data = tsdata, trace = TRUE)
coef(fit.ARIMA3, matrix = TRUE)

VGAMextra documentation built on Nov. 2, 2023, 5:59 p.m.