starma: Space-time series estimation procedure

In starma: Modelling Space Time AutoRegressive Moving Average (STARMA) Processes

Description

starma fits a STARMA model to a space-time series. It is the central function for the estimation part of the three-stage iterative model building procedure.

Usage

starma(data, wlist, ar, ma, iterate=1)

## S3 method for class 'starma'
print(x, ...)

Arguments

data	a matrix or data frame containing the space-time series: row-wise should be the temporal observations, with each column corresponding to a site.
wlist	a list of the weight matrices for each k-th order neighbours, first one being the identity.
ar	either an integer specifying the maximum time lag of the AR part, or a matrix filled with 0 or 1 indicating whether 'row'-th tlag, 'col'-th slag AR parameter should be estimated.
ma	either an integer specifying the maximum time lag of the MA part, or a matrix filled with 0 or 1 indicating whether 'row'-th tlag, 'col'-th slag MA parameter should be estimated.
iterate	an integer specifying how many times the Kalman filter should be re-run on itself (see Details).
x	a starma class object.
...	unused

Details

The definition here used for STARMA models is the following:

z(t) = phi[10]z(t-1) + phi[11]W[1]z(t-1) + ... + phi[pl]W[l]z(t-p) + theta[10]e(t-1) + theta[11]W[1]e(t-1) + ... + theta[qm]W[m]e(t-q) + e(t)

starma uses a Kalman filter algorithm (Cipra and Motykova, 1987): the parameters are set as the state vector of the state space system, making the iterations of the algorithm estimate directly the parameters. Thus, no optimization routine is required, making the algorithm extremely efficient time wise and computationally wise. Furthermore, the code has been written in C++ using Rcpp and RcppArmadillo (Eddelbuettel and Sanderson, 2014).

Note that, as the residuals must be iteratively estimated when running the Kalman filter, a single run might lead to poor results when estimating an MA parameter. Re-running the Kalman filter at least once, using the previously estimated parameters to add prior knowledge on the residuals leads to better estimates. For STAR model (when no MA parameter needs be estimated), the function ignores the iterate argument.

One of the strength of this estimation function is that the user can to estimate as few parameters as needed, even at high time and/or space lags, since the possibility to input a 1/0 matrix as AR and MA orders is given.

Value

A list of class starma containing the following elements:

phi	The estimated AR parameters
phi_sd	The corresponding standard errors
theta	The estimated MA parameters
theta_sd	The corresponding standard errors
sigma2	The white noise variance matrix estimated by the Kalman filter. Note that, to achieve parcimony, only the mean of the diagonal elements should be kept (since the noise is supposed to be Gaussian anyway)
residuals	The estimated residuals of the model
loglik	The conditional log likelihood of the model
bic	The corresponding BIC
call	The function call
df	Degrees of freedom of the model: (nb of obs) - (nb of parameters)

Author(s)

Felix Cheysson

References

- Cipra, T., & Motykova, I. (1987). Study on Kalman filter in time series analysis. Commentationes Mathematicae Universitatis Carolinae, 28(3).

- Dirk Eddelbuettel, Conrad Sanderson (2014). RcppArmadillo: Accelerating R with high-performance C++ linear algebra. Computational Statistics and Data Analysis, Volume 71, March 2014, pages 1054-1063. URL http://dx.doi.org/10.1016/j.csda.2013.02.005

Examples

data(nb_mat)	# Get neighbourhood matrices

# Simulate a STARMA model
eps <- matrix(rnorm(94*200), 200, 94)
sim <- eps
for (t in 3:200) {
	sim[t,] <- (.4*diag(94) + .25*blist[[2]]) %*% sim[t-1,] +
		(.25*diag(94)                ) %*% sim[t-2,] +
		(            - .3*blist[[2]]) %*% eps[t-1,] +
		eps[t, ]
}

sim <- sim[101:200,]
sim <- stcenter(sim)	# Center and scale the dataset

# Autocorrelation functions
stacf(sim, blist)
stpacf(sim, blist)

# Select parameters to estimate
ar <- matrix(0, 2, 2)
ar[ ,1] <- 1	# phi10 and phi20
ar[1,2] <- 1	# phi11
ma <- matrix(c(0,1), 1, 2)	# theta11

# Run the Kalman filter algorithm
model <- starma(sim, blist, ar, ma)
summary(model)

Example output

Call:
starma.default(data = sim, wlist = blist, ar = ar, ma = ma)

         Estimate Std..Error t.value   p.value    
phi10    0.400318   0.010056  39.810 < 2.2e-16 ***
phi11    0.261516   0.022803  11.468 < 2.2e-16 ***
phi20    0.257946   0.010004  25.785 < 2.2e-16 ***
theta11 -0.333273   0.031684 -10.519 < 2.2e-16 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

starma documentation built on May 2, 2019, 1:26 p.m.