#

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

#

# Written by Dr. Dimitrios D. Thomakos

# University of Peloponesse

# e-mail: thomakos@uop.gr

#

# STSA Version 2.0

#

# Example file #24

#

# In this file we illustrate the use of the Holt-Winters method for forecasting; the data are

# simulated and the model used to create them is commented out...

#

# Initialize

clear

ginit

clc

ranseed

# Input simulated data

fname = "data\holt.prn"

N = 200.

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

#begin##

# Simulate a structural model of deterministic plus stochastic trend plus seasonal plus noise

N = 200.

sigma = 0.4

u = snormal(N+1,1)*sigma

# Add the deterministic component here

A = 0.7

B = -1.

omega = 12.

t = seq(int(N+1))

r = 0.01*t + A%sin(omega*t+B)

# Now, generate the stochastic trend and add the noise

theta = 0.2

e = snormal(N+1,1)*theta

x = zeros(N+1,1)

y = zeros(N+1,1)

for i = 2 to N+1 begin

x(i) = x(i-1) + e(i)

y(i) = x(i) + r(i) + u(i)

end

# Crop and plot the data

y = y(2::end)

#end##

pcapt = "Simulated Series"

ptitle = "Deterministic+Stochastic Trend plus Seasonal plus Noise"

pcolor = "red"

pstyle = "solid"

pgrid  = "on"

xtitle = "Observation"

ytitle = ""

fast_plot(y,pcapt,ptitle,pcolor,pstyle,pgrid,xtitle,ytitle)

# Now, optimize the Holt-Winters parameters and fit the series

s = 12.

h = 1.

sv = { 0.1, 0.1, 0.1 }

ht = 1e-3

tol = 1e-4

miter = 150.

level = 10.

out   = 1.

[lambda,rss,cov] = holt_winters_optimize(y,s,h,sv,ht,tol,miter,level,out)

[m,b,c,f,r] = holt_winters_filter(y,s,h,lambda)

# Finally, plot the series along with the estimated level, slope and seasonal components

#

# Note that as the series levels out after about observation 100 the estimated slope

# component has a negative slope...

#

gaddwin("Actual values and underlying level component")

gplot([y,m])

gaddwin("Slope component")

gplot(b(s+1::end))

gaddwin("Cyclical component")

gplot(c(s+1::end))