## Chapter 9 ## Sec 9.1: Time Series -- Some Basic Ideas ## ss 9.1.1: Preliminary graphical explorations ## Plot depth measurements: ts object LakeHuron (datasets) plot(LakeHuron, ylab="depth (in feet)", xlab = "Time (in years)") lag.plot(LakeHuron, lags=3, do.lines=FALSE) ## ss 9.1.2: The autocorrelation function acf(LakeHuron) ## ss 9.1.3: Autoregressive (AR) models ## The {\rm AR(1)} model ## Footnote code ## series.diff <- diff(series) ## Footnote code LH.yw <- ar(x = LakeHuron, order.max = 1, method = "yw") # autocorrelation estimate # order.max = 1 for the AR(1) model LH.yw$ar # autocorrelation estimate of alpha LH.mle <- ar(x = LakeHuron, order.max = 1, method = "mle") # maximum likelihood estimate LH.mle$ar # maximum likelihood estimate of alpha ## Estimates of the series mean and of sigma^2 can be obtained by typing LH.mle$x.mean # estimated series mean LH.mle$var.pred # estimated innovation variance ## The general {\rm AR(}p{\rm )} model ## Footnote code LH.mle2 <- ar(LakeHuron, order.max=2, method="mle") ## Footnote code ar(LakeHuron, method="mle") # AIC is used by default if # order.max is not specified ## Footnote code acf(LakeHuron, type="partial") ## ss 9.1.4: {\kern-4pt}$^*$~Autoregressive moving average (ARMA) models -- theory ## Sec 9.2: {\kern-4pt}$^*$~Regression modeling with moving average errors ## Plot time series avrain and SOI: ts object bomsoi (DAAG) plot(ts(bomsoi[, 15:14], start=1900), panel=function(y,...)panel.smooth(bomsoi$Year, y,...)) xbomsoi <- with(bomsoi, data.frame(SOI=SOI, cuberootRain=avrain^0.33)) xbomsoi$trendSOI <- lowess(xbomsoi$SOI)$y xbomsoi$trendRain <- lowess(xbomsoi$cuberootRain)$y rainpos <- pretty(bomsoi$avrain, 5) with(xbomsoi, {plot(cuberootRain ~ SOI, xlab = "SOI", ylab = "Rainfall (cube root scale)", yaxt="n") axis(2, at = rainpos^0.33, labels=paste(rainpos)) ## Relative changes in the two trend curves lines(lowess(cuberootRain ~ SOI)) lines(lowess(trendRain ~ trendSOI), lwd=2, col="gray40") }) xbomsoi$detrendRain <- with(xbomsoi, cuberootRain - trendRain + mean(trendRain)) xbomsoi$detrendSOI <- with(xbomsoi, SOI - trendSOI + mean(trendSOI)) ## Footnote code oldpar <- par(mfrow=c(1,2), pty="s") plot(cuberootRain ~ SOI, data = xbomsoi, ylab = "Rainfall (cube root scale)", yaxt="n") axis(2, at = rainpos^0.33, labels=paste(rainpos)) with(xbomsoi, lines(lowess(cuberootRain ~ SOI))) plot(detrendRain ~ detrendSOI, data = xbomsoi, xlab="Detrended SOI", ylab = "Detrended rainfall", yaxt="n") axis(2, at = rainpos^0.33, labels=paste(rainpos)) with(xbomsoi, lines(lowess(detrendRain ~ detrendSOI))) par(oldpar) attach(xbomsoi) xbomsoi.ma0 <- arima(detrendRain, xreg=detrendSOI, order=c(0,0,0)) xbomsoi.ma12 <- arima(detrendRain, xreg=detrendSOI, order=c(0,0,12)) xbomsoi.ma12s <- arima(detrendRain, xreg=detrendSOI, seasonal=list(order=c(0,0,1), period=12)) xbomsoi.ma12 xbomsoi.maSel <- arima(x = detrendRain, order = c(0, 0, 12), xreg = detrendSOI, fixed = c(0, 0, 0, NA, rep(0, 4), NA, 0, NA, NA, NA, NA)) xbomsoi.maSel xbomsoi0.maSel <- arima(x = detrendRain, order = c(0, 0, 12), fixed = c(0, 0, 0, NA, rep(0, 4), NA, 0, NA, NA, NA)) detach(xbomsoi) ## Other checks and computations Box.test(resid(lm(detrendRain ~ detrendSOI, data = xbomsoi)), type="Ljung-Box", lag=20) attach(xbomsoi) xbomsoi2.maSel <- arima(x = detrendRain, order = c(0, 0, 12), xreg = poly(detrendSOI,2), fixed = c(0, 0, 0, NA, rep(0, 4), NA, 0, rep(NA,5))) xbomsoi2.maSel ## Examine normality of estimates of "residuals" ("innovations") qqnorm(resid(xbomsoi.maSel, type="normalized")) detach(xbomsoi) ## Sec 9.3: *Nonlinear time series x <- numeric(999) # x will contain the ARCH(1) errors x0 <- rnorm(1) for (i in 1:999){ x0 <- rnorm(1, sd=sqrt(.25 + .95*x0^2)) x[i] <- x0 } library(tseries) garch(x, order = c(0, 1), trace=FALSE) ## Sec 9.4: Other time series packages ## Sec 9.5: Further Reading ## Spatial statistics Box.test(resid(lm(detrendRain ~ detrendSOI, data = xbomsoi)), type="Ljung-Box", lag=20) xx <- matrix(x, ncol=1000} library(tseries) data(ice.river) # Needed with version 0.9-30 of tseries river1 <- diff(log(ice.river[, 1]))