This is an R Markdown document. Markdown is a simple formatting syntax for authoring HTML, PDF, and MS Word documents. For more details on using R Markdown see http://rmarkdown.rstudio.com.

When you click the Knit button a document will be generated that includes both content as well as the output of any embedded R code chunks within the document. You can embed an R code chunk like this:

summary(cars)
##      speed           dist       
##  Min.   : 4.0   Min.   :  2.00  
##  1st Qu.:12.0   1st Qu.: 26.00  
##  Median :15.0   Median : 36.00  
##  Mean   :15.4   Mean   : 42.98  
##  3rd Qu.:19.0   3rd Qu.: 56.00  
##  Max.   :25.0   Max.   :120.00

You can also embed plots, for example:

Note that the echo = FALSE parameter was added to the code chunk to prevent printing of the R code that generated the plot.

download data

install.packages(“Quandl”)

library(“Quandl”) #(1)model identification

rcpg <- Quandl(“FRED/WCOILBRENTEU”, type=“zoo”)

summary(rcpg)

plot(rcpg)

str(rcpg)

head(rcpg)

tail(rcpg)

transformed data

g=log(rcpg)

Y=diff(g) Y plot(Y)

plot ACF for y up to lag 24

acf(as.data.frame(Y),type=‘correlation’,lag=24)

plot PACF for y up to lag 24

acf(as.data.frame(Y),type=‘partial’,lag=24, main=“”)

(2) Model Estimation

find the order of AR by two ways 1- PACF here we find the order is 2 from PACF plot, 2-information cariteria AIC here an AR(3)is selected based on AIC

estimate an AR(2) model - there is 2 significant coefficients in the PACF plot for y

m1 <- arima(y, order=c(1,0,0))

show the structure of object m1

str(m1)

diagnostics for the AR(2) model -

tsdiag(m1,gof.lag=24) # we prefer MA(2) to MA(1) B1=arima(Y,order=c(0,0,2),method=“ML”) B1 str(B1)

B2=arima(Y,order=c(0,0,1),method=“ML”) B2 str(B2)

tsdiag(A,gof.lag=24) # estimate an AR(2) model m2 <- arima(Y, order=c(2,0,0)) m2 # diagnostics for the AR(2) model - tests for corelation 1- L-B, 2-B-p ,from L-B test we find that transformed data are not correlated(not significant) tsdiag(m2,gof.lag=24)

estimate an AR(3) model since PACF for lag 2 and 3 are comparable in size

m3 <- arima(Y, order=c(3,0,0)) m3 # diagnostics for the AR(3) model,tests for corelation 1- L-B, 2-B-p ,from L-B test we find that transformed data are correlated( significant) tsdiag(m3,gof.lag=24)

z-statistics for coeffici

m1=ar(Y,method=“mle”)

m1

m1$order

A=arima(Y,order=c(1,0,0),method=“ML”)

m2=arima(Y, order=c(2,0,0))

m2

m3=arima(Y, order=c(3,0,0))

m3

z-statistics for coefficients of AR(3) model - phi2 is signifficant at 5% level, phi3 is marginally insignifficant

m3\(coef/sqrt(diag(m3\)var.coef)) # p values (1-pnorm(abs(m3\(coef)/sqrt(diag(m3\)var.coef))))*2

use AIC to choose order p of the AR model

m <- ar(Y,method=“mle”) str(m) m m\(order # note that AIC prefers AR(3) to AR(2) m\)aic

note that BIC prefers AR(2) to AR(3), in general BIC puts a larger penalty on additional coefficients than AIC

BIC(m2) BIC(m3)

(3)checking for adequacy of the model

Ljung-Box test - for residuals of a model adjust the degrees of freedom m by subtracting the number of parameters g

m2.LB.lag12 <- Box.test(m2\(residuals, lag=12, type="Ljung") str(m2.LB.lag12) m2.LB.lag12 1-pchisq(m2.LB.lag12\)statistic,10)

m2.LB.lag16 <- Box.test(m2\(residuals, lag=16, type="Ljung") m2.LB.lag16 1-pchisq(m2.LB.lag16\)statistic,14)

the constant term is obtained

(1-0.1866+0.0215-0.0907)*0.0004

so the fitted model is Yt=0.00029768+0.1866Yt-1-0.0215Yt-2+0.0907Yt-3+at

residual standard error

sqrt(m2$sigma2)

caracteristic equation

p1=c(1, -m2$coef[1:3])

find solution

roots=polyroot(p1)

roots

compute the absolote values of the solutions

Mod(roots)

to compute average length of business cycles(quarters) it is about 3 years

k=2*pi/acos(1.989104/2.354474)

k

AIC

ord=ar(Y,method=“mle”)

ord$aic

ord$order

m3=arima(Y,order=c(3,0,0))

m3

compute the intercept pi(0)

(1-0.1866+0.0215-0.0907)*mean(Y)

compute stander error of residuals

sqrt(m3$sigma2)

Box.test(m3$residuals,lag=12,type=“Ljung”)

compute p-value using 9 degrees of freedom

pv=1-pchisq(9.9132,9)

pv