PC session 2

[Liquor] Clear the workspace

rm(list=ls())

[Liquor] Read the data file (decimal point)

mydata <- read.table(file=file.choose(),header=TRUE)

[Liquor] Use variable names directly

attach(mydata)

[Liquor] Declare monthly time series starting July 1985

liquor_ts <- ts(liquor,frequency=12,start=c(1985,7))

[Liquor] Plot the series

ts.plot(liquor_ts)

[Liquor] Log-transform the series

logliquor_ts <- log(liquor_ts)

[Liquor] Plot the log series

ts.plot(logliquor_ts)

[Liquor] Log-differences (growth)

dlogliquor_ts <- diff(log(liquor_ts))

[Liquor] Plot the log-differenced series

ts.plot(dlogliquor_ts)

[Liquor] Seasonal plot of the log-differenced series

monthplot(dlogliquor_ts)

[Liquor] Log, first difference and seasonal difference (lag 12)

dlogsliquor_ts <- diff(diff(log(liquor_ts), lag = 12))

[Liquor] Plot the log series in first and seasonal differences

ts.plot(dlogsliquor_ts)

[Liquor] Load the unit root test package

library(CADFtest)

[Liquor] Maximum lag = square root of sample size of dlogliquor

max.lag <- round(sqrt(length(dlogliquor_ts)))

[Liquor] Unit root test with intercept, type with drift, lag selection by BIC, max.lag.y

CADFtest(dlogsliquor_ts, type="drift", criterion="BIC", max.lag.y=max.lag)

[Liquor] Put two plots above each other

par(mfrow=c(2,1))

[Liquor] Correlogram

acf(dlogsliquor_ts)

[Liquor] Partial correlogram of dlogsliquor

pacf(dlogsliquor_ts)

[Liquor] Ljung-Box (Q) test for white noise

Box.test(dlogsliquor_ts, lag=max.lag, type="Ljung-Box")

[Liquor] Fit ARIMA(3,1,0)(0,1,0) on the log series

fit_ar <- arima(logliquor_ts, order=c(3,1,0), seasonal=list(order=c(0,1,0)))

[Liquor] Look at the model output of fit_ar

fit_ar

[Liquor] Absolute t-values of all coefficients (significant if > 1.96)

abs(fit_ar$coef/sqrt(diag(fit_ar$var.coef)))

[Liquor] Plot the residuals of ar

ts.plot(fit_ar$residuals)

[Liquor] Correlogram of the residuals

acf(fit_ar$residuals)

[Liquor] Q-test on the residuals

Box.test(fit_ar$residuals, lag=max.lag, type="Ljung-Box")

[Liquor] Fit alternative ARIMA(0,1,4)(0,1,0)

fit_ma <- arima(logliquor_ts, order=c(0,1,4), seasonal=list(order=c(0,1,0)))

[PPI] Clear the workspace

rm(list=ls())

[PPI] Read the data file (decimal comma)

mydata <- read.table(file=file.choose(),header=TRUE,dec=",")

[PPI] Use variable names directly

attach(mydata)

[PPI] Declare monthly time series starting January 1981

ppi_ts <- ts(PPI,frequency=12,start=c(1981,1))

[PPI] Log-transform the series

logppi_ts <- log(ppi_ts)

[PPI] Plot the log series

ts.plot(logppi_ts)

[PPI] Check the number of rows and columns

dim(mydata)

[PPI] Create a deterministic trend variable

TREND <- 1:372

[PPI] Regress log(PPI) on a constant and trend

fit <- lm(logppi_ts ~ TREND)

[PPI] Look at the regression output

summary(fit)

[PPI] Plot the regression residuals

ts.plot(fit$residuals)

[PPI] Load the unit root test package

library(CADFtest)

[PPI] Maximum lag = square root of sample size

max.lag <- round(sqrt(372))

[PPI] Unit root test with intercept and trend

CADFtest(logppi_ts, type="trend", criterion="BIC", max.lag.y=max.lag)

[PPI] Log-differences

dlogppi_ts <- diff(log(ppi_ts))

[PPI] Plot the log-differenced series

ts.plot(dlogppi_ts)

[PPI] Take subsample from 1997 to 2011

ppi_cut_ts <- window(ppi_ts,start=1997,end=2011)

[PPI] Log-differences of the subsample

dlogppi_cut_ts <- diff(log(ppi_cut_ts))

[PPI] Plot the subsample

ts.plot(dlogppi_cut_ts)

[PPI] Seasonal plot of the subsample

monthplot(dlogppi_cut_ts)

[PPI] Subsample in first and seasonal differences

dlogsppi_cut_ts <- diff(diff(log(ppi_cut_ts),lag=12))

[PPI] Maximum lag for the subsample

max.lag <- round(sqrt(length(dlogsppi_cut_ts)))

[PPI] Unit root test with intercept on the subsample, type is drift

CADFtest(dlogsppi_cut_ts, type="drift", criterion="BIC", max.lag.y=max.lag)

[PPI] Two plots above each other

par(mfrow=c(2,1))

[PPI] Correlogram of the subsample

acf(dlogsppi_cut_ts)

[PPI] Partial correlogram of the subsample

pacf(dlogsppi_cut_ts)

[PUNEW] Clear the workspace

rm(list=ls())

[PUNEW] Read the data file (decimal point)

mydata <- read.table(file=file.choose(),header=TRUE,dec=".")

[PUNEW] Use variable names directly

attach(mydata)

[PUNEW] Declare quarterly time series starting 1989 Q1

PUNEW_ts <- ts(PUNEW,frequency=4,start=c(1989,1))

[PUNEW] Log-transform the series

logPU_ts <- log(PUNEW_ts)

[PUNEW] Plot the log series

ts.plot(logPU_ts)

[PUNEW] Load the unit root test package

library(CADFtest)

[PUNEW] Check the number of rows

dim(mydata)

[PUNEW] Maximum lag = square root of sample size

max.lag <- round(sqrt(65))

[PUNEW] Unit root test with intercept and trend

CADFtest(logPU_ts, type="trend", criterion="BIC", max.lag.y=max.lag)

[PUNEW] Log-differences

dlogPU_ts <- diff(log(PUNEW_ts))

[PUNEW] Plot the log-differenced series

plot.ts(dlogPU_ts)

[PUNEW] Seasonal plot

monthplot(dlogPU_ts)

[PUNEW] Unit root test with intercept on the differenced series

CADFtest(dlogPU_ts, type="drift", criterion="BIC", max.lag.y=max.lag)

[PUNEW] Two plots above each other

par(mfrow=c(2,1))

[PUNEW] Correlogram (suggests MA order)

acf(dlogPU_ts)

[PUNEW] Partial correlogram (suggests AR order)

pacf(dlogPU_ts)

[PUNEW] Fit ARIMA(3,1,0) without seasonal part

fit_ar <- arima(logPU_ts, order=c(3,1,0), seasonal=c(0,0,0))

[PUNEW] Look at the model output

fit_ar

[PUNEW] Absolute t-values of all coefficients

abs(fit_ar$coef/sqrt(diag(fit_ar$var.coef)))

[PUNEW] Plot the residuals

plot.ts(fit_ar$residuals)

[PUNEW] Correlogram of the residuals

acf(fit_ar$residuals)

[PUNEW] Q-test on the residuals

Box.test(fit_ar$residuals, lag=max.lag, type="Ljung-Box")

[PUNEW] Fit ARIMA(3,1,1)

fit_arma <- arima(logPU_ts, order=c(3,1,1), seasonal=c(0,0,0))

[PUNEW] Correlogram of the ARMA residuals

acf(fit_arma$residuals)

[PUNEW] Partial correlogram of the ARMA residuals

pacf(fit_arma$residuals)

[PUNEW] Q-test on the ARMA residuals

Box.test(fit_arma$residuals, lag=max.lag, type="Ljung-Box")

[PUNEW] AIC of the AR(3) model

AIC(fit_ar)

[PUNEW] AIC of the ARMA(3,1) model

AIC(fit_arma)

[PUNEW] SIC of the AR(3) model (n = 64)

AIC(fit_ar,k=log(64))

[PUNEW] SIC of the ARMA(3,1) model (n = 64)

AIC(fit_arma,k=log(64))