I prepared the following notes for a consulting client, and I thought they might be of interest to some other people too.

Let denote the value of the time series at time , and suppose we wish to fit a trend with correlated errors of the form

where represents the possibly nonlinear trend and is an autocorrelated error process.

For example, if is a linear function, then we can simply set and define

In matrix form we can write

where , , and . Note that I have left the intercept out of the vector so that the matrix matches the required xreg argument in auto.arima.

This model can be estimated by setting the xreg argument to be a matrix with one column:

x1 <- 1:length(y)
fit <- auto.arima(y, xreg=x1)

The associated coefficient is the slope of the trend line.

Here is a simple example of a linear trend fitted to the Asian sheep data from the fpp package :

library(fpp)
T <- length(livestock)
x1 <- seq(T)
fit <- auto.arima(livestock, xreg=x1)
fc <- forecast(fit, xreg=T+seq(10))
b0 <- coef(fit)["intercept"]
b1 <- coef(fit)["x1"]
t <- seq(T+10)
trend <- ts(b0 + b1*t, start=start(livestock))
 
plot(fc, main="Linear trend with AR(1) errors")
lines(trend, col='red')

A linear trend fitted to the Asian sheep data. The automatically selected error term is an AR(1) process.

A more flexible approach is to use a piecewise linear trend which bends at some time. If the trend bends at time , then it can be specified by including the following predictors in the model.

(1)  

In auto.arima, set xreg to be a matrix with two columns:

fit <- auto.arima(y, xreg=cbind(x1, pmax(0,x1-tau))

If the associated coefficients of and are and , then gives the slope of the trend before time , while the slope of the line after time is given by .

This can be extended to allow any number of “bend points” known as knots. Just add additional columns with 0s before each knot, and values 1, 2, … after the knot.

Here is a piecewise linear trend fitted to the Asian sheep data with knots at years 1990 and 1992:

x2 <- pmax(0, x1-30)
x3 <- pmax(0, x1-32)
fit <- auto.arima(livestock, xreg=cbind(x1,x2,x3))
fc <- forecast(fit, 
   xreg=cbind(max(x1)+seq(10), max(x2)+seq(10), max(x3)+seq(10)))
b0 <- coef(fit)["intercept"]
b1 <- coef(fit)["x1"]
b2 <- coef(fit)["x2"]
b3 <- coef(fit)["x3"]
trend <- ts(b0 + b1*t + b2*pmax(0,t-30) + b3*pmax(0,t-32),
    start=start(livestock))
 
plot(fc, main="Piecewise linear trend with AR(1) errors")
lines(trend, col='red')

A piecewise-linear trend fitted to the Asian sheep data

If there is to be no trend before the first knot, but a piecewise linear trend thereafter, leave out the first column of the above matrix .

If there is to be a piecewise linear trend up to the last knot, but no trend thereafter, a slightly modified set up can be used. For one knot at time , we can set

xreg <- pmin(0, x1-tau)

where the first 0 in the column is in row . Additional knots can be handled in the same way. For example, if there are two knots, then will be the slope of the trend up to the first knot, and will be the slope between the first and second knots.