CUHK- Shenzhen STAT 401 - Seasonal and Related Models

Unformatted text preview:

What is this all about?Data/Theory goals of the topicSeasonal modelsSeasonal unit roots and related discussionPeriodic modelsSome final notesReferencesSeasonal and Related ModelsVladas Pipiras, STOR @ UNC-CHFebruary, 2022What is this all about?Consider the following time series and its analysis.ts0 <- scan("airpass.txt")ts <- log(ts0)par(mfrow = c(1, 2))plot.ts(ts0)plot.ts(ts)Timets00 20 40 60 80 100 120 140100 300 500Timets0 20 40 60 80 100 120 1405.0 5.5 6.0 6.5n <- length(ts)tt <- seq(1,n)tt2 <- ttˆ2cos12 <- cos(2*pi*tt/12)sin12 <- sin(2*pi*tt/12)cos6 <- cos(2*pi*tt/6)sin6 <- sin(2*pi*tt/6)cos4 <- cos(2*pi*tt/4)sin4 <- sin(2*pi*tt/4)rg <- lm(ts ~ tt + tt2 + cos12 + sin12 + cos6 + sin6 + cos4 + sin4)summary(rg)#### Call:## lm(formula = ts ~ tt + tt2 + cos12 + sin12 + cos6 + sin6 + cos4 +## sin4)#### Residuals:1## Min 1Q Median 3Q Max## -0.173540 -0.035102 0.002494 0.038911 0.128827#### Coefficients:## Estimate Std. Error t value Pr(>|t|)## (Intercept) 4.736e+00 1.417e-02 334.223 < 2e-16 ***## tt 1.319e-02 4.508e-04 29.247 < 2e-16 ***## tt2 -2.148e-05 3.011e-06 -7.132 5.45e-11 ***## cos12 -1.415e-01 6.582e-03 -21.499 < 2e-16 ***## sin12 -4.927e-02 6.594e-03 -7.471 8.93e-12 ***## cos6 -2.275e-02 6.582e-03 -3.457 0.000732 ***## sin6 7.872e-02 6.584e-03 11.957 < 2e-16 ***## cos4 2.731e-02 6.582e-03 4.149 5.87e-05 ***## sin4 -8.706e-03 6.582e-03 -1.323 0.188144## ---## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1#### Residual standard error: 0.05584 on 135 degrees of freedom## Multiple R-squared: 0.9849, Adjusted R-squared: 0.984## F-statistic: 1100 on 8 and 135 DF, p-value: < 2.2e-16ts2 <- residuals(rg)par(mfrow = c(1, 2))plot.ts(ts)lines(fitted.values(rg),col="red")plot.ts(ts2)Timets0 20 40 60 80 100 120 1405.0 5.5 6.0 6.5Timets20 20 40 60 80 100 120 140−0.15 −0.05 0.05par(mfrow = c(1, 2))acf(ts2,lag.max = 30)pacf(ts2,lag.max = 30)20 5 10 15 20 25 30−0.4 0.0 0.4 0.8LagACFSeries ts20 5 10 15 20 25 30−0.3 −0.1 0.1 0.3LagPartial ACFSeries ts2ts20 <- diff(ts,lag = 12)ts2a <- diff(ts20)par(mfrow = c(1, 2))plot.ts(ts20)plot.ts(ts2a)Timets200 20 40 60 80 100 1200.0 0.1 0.2 0.3Timets2a0 20 40 60 80 100 120−0.15 −0.05 0.05 0.15par(mfrow = c(1, 2))acf(ts2a,lag.max = 30)pacf(ts2a,lag.max = 30)0 5 10 15 20 25 30−0.4 0.0 0.4 0.8LagACFSeries ts2a0 5 10 15 20 25 30−0.3 −0.1 0.1LagPartial ACFSeries ts2a3Note that the resulting series exhibit strong (partial) correlations at some lags which are multiples of theperiod s. That is, bρ(s), bρ(2s), . . ., bα(s), bα(2s), . . .,are large.Data/Theory goals of the topicIntroduce (seasonal) models that can capture the above behavior. Discuss a few related issues, and alsoanother model with seasonal (periodic) features.Seasonal modelsMotivating seasonal models: Consider time series data for r years and 12 months:Year\Month 1 2 . . . 121 Y1Y2. . . Y122 Y13Y14. . . Y24...............r Y12(r−1)+1Y12(r−1)+2. . . Y12rThere arej= 1, . . . ,12 column seriesYj+12t,t= 0, . . . , r −1. Assume the same AR model for each:j = 1, . . . , 12,Yj+12t− Φ1Yj+12(t−1)− . . . − ΦPYj+12(t−P )= Uj+12twith {Uj+12t} ∼ WN(0, σ2). This is the same asΦ(B12)Yt= Utwith Φ(z) = 1 − Φ1z − . . . − ΦPzP. Terminology: between-year model.B Example: P = 1, Φ1= 0.7.To incorporate dependence between months, take another AR model for {Ut}:φ(B)Ut= Ztwith φ(z) = 1 − φ1z − . . . − φpzpand {Zt} ∼ WN(0, σ2Z). Combining the two models:φ(B)Φ(B12)Yt= Zt.What one obtains is known as a seasonal ARIMA model.Definition d, Dnon-negative integers:{Xt}is a seasonal ARIMA(p, d,0)×(P, D,0)stime series with periods if the differenced time series Yt= (I − B)d(I − Bs)DXtis a causal ARMA time series defined byφ(B)Φ(Bs)Yt= Ztwith φ(z) = 1 − φ1z − . . . − φpzp, Φ(z) = 1 − Φ1z − . . . − ΦPzPand {Zt} ∼ WN(0, σ2Z).MA parts can be added similarly.B Note the multiplicative nature of SARIMA models.Example: SARIMA(1, 0, 0) × (1, 0, 0)12with φ1= 0.7, Φ1= 0.6.library(astsa)set.seed(123)ts.e1 <- sarima.sim(ar=0.7, sar=0.6, S=12,n=400)ts.e1 <- as.numeric(ts.e1)4par(mfrow = c(1, 3))plot.ts(ts.e1)acf(ts.e1,lag.max = 30)pacf(ts.e1,lag.max = 30)Timets.e10 100 200 300 400−4 −2 0 2 40 5 10 15 20 25 300.0 0.2 0.4 0.6 0.8 1.0LagACFSeries ts.e10 5 10 15 20 25 30−0.2 0.0 0.2 0.4 0.6LagPartial ACFSeries ts.e1set.seed(1)ts.e2 <- arima.sim(list(order=c(13,0,0),ar=c(.7,rep(0,10),.6,-.7*.6)),n=400)par(mfrow = c(1, 1))plot.ts(ts.e2)Timets.e20 100 200 300 400−4 −2 0 2par(mfrow = c(1, 2))acf(ts.e2,lag.max = 30)pacf(ts.e2,lag.max = 30)library(forecast)50 5 10 15 20 25 300.0 0.4 0.8LagACFSeries ts.e20 5 10 15 20 25 30−0.2 0.2 0.4 0.6LagPartial ACFSeries ts.e2ar_acf<-ARMAacf(ar=c(.7,rep(0,10),.6,-.7*.6),lag.max=30,pacf=F)ar_pacf<-ARMAacf(ar=c(.7,rep(0,10),.6,-.7*.6),lag.max=30,pacf=T)par(mfrow = c(1, 2))plot(ar_acf, type='h')plot(ar_pacf, type='h')0 5 10 15 20 25 300.2 0.4 0.6 0.8 1.0Indexar_acf0 5 10 15 20 25 30−0.4 0.0 0.4Indexar_pacfH Homework 5 will ask you to explore similarly another SARIMA model.A number of R functions used earlier (model selection, fitting, forecasting, etc.) allow for seasonal models.library(forecast)set.seed(123)ts.e1 <- sarima.sim(ar=0.7, sar=0.6, S=12,n=400)auto.mod <- auto.arima(ts.e1, max.p = 2, max.q = 2, max.d = 0,start.p = 0, start.q = 0,max.P = 2, max.Q = 2, max.D = 0,start.P = 0, start.Q = 0,allowdrift = FALSE, ic = "aic")auto.mod## Series: ts.e1## ARIMA(1,0,0)(1,0,0)[12] with zero mean#### Coefficients:## ar1 sar16## 0.6061 0.5298## s.e. 0.0396 0.0426#### sigma^2 estimated as 0.9612: log likelihood=-560.87## AIC=1127.73 AICc=1127.79 BIC=1139.71sarima.model <- Arima(ts.e1,order=c(1,0,0),seasonal=c(1,0,0),method = "ML",include.mean = TRUE)sarima.model## Series: ts.e1## ARIMA(1,0,0)(1,0,0)[12] with non-zero mean#### Coefficients:## ar1 sar1 mean## 0.6057 0.5295 0.0851## s.e. 0.0396 0.0426 0.2541#### sigma^2 estimated as 0.9634: log likelihood=-560.81## AIC=1129.62 AICc=1129.72 BIC=1145.59h = 12ts.e1.forecast <- forecast(auto.mod, h)plot(ts.e1.forecast)Forecasts from ARIMA(1,0,0)(1,0,0)[12] with zero mean0 5 10 15 20 25 30 35−4 −2 0 2 4Back to the time series motivating this topic.ts0 <- scan("airpass.txt")ts1 <- log(ts0)ts <- ts(data=ts1,frequency=12)auto.ts <- auto.arima(ts, max.p = 2, max.q = 2, max.d = 1,start.p = 0, start.q = 0,max.P = 2, max.Q = 2, max.D = 1,start.P = 0, start.Q = 0,allowdrift = FALSE, ic = "aic")auto.ts## Series: ts7## ARIMA(0,1,1)(0,1,1)[12]####


View Full Document

CUHK- Shenzhen STAT 401 - Seasonal and Related Models

Download Seasonal and Related Models
Our administrator received your request to download this document. We will send you the file to your email shortly.
Loading Unlocking...
Login

Join to view Seasonal and Related Models and access 3M+ class-specific study document.

or
We will never post anything without your permission.
Don't have an account?
Sign Up

Join to view Seasonal and Related Models 2 2 and access 3M+ class-specific study document.

or

By creating an account you agree to our Privacy Policy and Terms Of Use

Already a member?