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

Description

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.

Usage

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)

Value

An object of class "vglmff"

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

Arguments

order

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

zero

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

diffCovs

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.

xLag

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.

include.current

Logical. Same as ARIMAX.errors.ff.

type.EIM

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, et.al. 1986), otherwise the approximate version is utilized.

var.arg

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

nodrift

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

noChecks

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.

Author

Victor Miranda and T. W. Yee

Warning

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.

Details

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, $$

where $$ \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.

References

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.

Examples

Run this code

set.seed(3)
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")
Coef(fit.ARIMA1)
summary(fit.ARIMA1)
vcov(fit.ARIMA1, untransform = TRUE)
#------------------------------------------------------------------------#
# Fitting same model using arima().
#------------------------------------------------------------------------#
# COMPARE to EXAMPLE1 
( 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)
summary(fit.ARIMA2)
constraints(fit.ARIMA2)



### 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)
summary(fit.ARIMA3)
constraints(fit.ARIMA3)

Run the code above in your browser using DataLab

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