#

#
STATISTICAL TIME SERIES ANALYSIS (STSA) module for O-Matrix.

#

# Written by Dr.
Dimitrios D. Thomakos

#

# e-mail:
thomakos@uop.gr

#

# STSA Version 2.0

#

# Example file #11

#

# In
this example we illustrate the estimation and forecasting of ARFIMA-type
models.

# This
class of models extends traditional, short-memory ARMA models, to incorporate
the effects

# of
long-memory.

#

# In
the example, we illustrate how to use the arfima_estimate and arfima_forecast
STSA functions

# and
we compare the forecasting performance of ARFIMA forecasts vs. forecasts
generated by a linear

#
AR model.

#

# Clear workspace
and command window, initialize graphics

clear

clc

ginit

ranseed

# Input simulated
data

fname = "data\fract.prn"

N = 400.

y = read(fname,"double",N,1)

# We
can estimate the fractional order and the autoregressive coefficient of the y
series

# using an ARFIMA(1,d,0) model...Note that the fractional order of the
y series is the

# same
as the fractional order of the x series.

#

sv = { 0.1, 0.1 } # Set starting values

# Function call, get
the estimates

[theta,S,cov_theta,e,iter]
= arfima_estimate(y,1,0,20,sv,5,1)

print

print "True parameter values were =
",[0.5,0.4]

print

# Now,
forecast the last few values of the series using the ARFIMA model and an AR
model

# and
compare their performance...

#

nf = 10. #
Number of points to forecast; this can be changed...

yin = y(1::end-nf)

yot = y(end-nf+1::end)

# Estimate and
forecast using the AR model

maxp = 10. #
max lag to consider; set it high enough as the ACF
decays slowly

sm = 0. #
flag; do not substract the sample mean

method = 1. #
flag; use LS for estimation

[aic,bic]
= ar_order(yin,maxp,sm,method)

paic = aic(1,1) #
get the selected lag order

sv = ones(paic,1)*0.1

[phi,e]
= arma_estimate(yin,paic,0,sv)

s2 = (e'*e)/(N-nf-paic) # get the
residual variance

[yf_ar,msef_ar] = arma_forecast(yin,paic,0,phi,s2,nf)

#

# Now
estimate and forecast using ARFIMA model

#

sv = { 0.1, 0.1 }

[theta,S,cov_theta,e,iter]
= arfima_estimate(yin,1,0,20,sv,5,1)

yf_arf = arfima_forecast(yin,1,0,theta,20,nf)

#

# Compare
forecasting performance

#

print

fe_ar = yot - yf_ar

fe_arf= yot - yf_arf

[mae_ar,mse_ar] = foreval(fe_ar)

[mae_arf,mse_arf] = foreval(fe_arf)

print "MAE of AR model = ", mae_ar

print "MAE of ARFIMA model = ", mae_arf

print

print "MSE of AR model = ", mse_ar

print "MSE of ARFIMA model = ", mse_arf

# Now,
plot the actual values and the forecasts from the two models

gcolor({"black","red","green"})

gstyle("star")

gtitle("Actuals (stars) and forecasts (red =
AR, green = ARFIMA)")

gplot(yot)

gstyle("solid")

gplot(yf_ar)

gplot(yf_arf)