#

# 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 #09

#

# In this example we illustrate the differences between the two available kernels for nonparametric

# estimation and forecasting. Then, we illustrate the differences in the fit between a linear and a

# nonparametric regression when the true model is nonlinear but unknown.

#

# Clear workspace and command window, initialize graphics, seed the random number generator

clear

clc

ginit

ranseed

# Plot the optimal and Gaussian kernels

x = -3.5:0.1:3.5         # grid of x values

x0 = 1.             # evaluation point; you can change this

h  = 0.1            # bandwidth; you can change this

z  = (x-x0)/h       # standardized variate

y2 = gaussian_kernel(z)  # the Gaussian kernel

y1 = optimal_kernel(z)   # the optimal kernel

# Now plot them

gaddwin("The optimal (green) and Gaussian (red) kernel functions")

format double "f4.2"

title = ["Evaluation point is x0 = ",ntoa(x0)," and bandwidth is h = ",ntoa(h)]

gtitle(title)

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

gxaxis("linear",0.,1.5)

gplot(x,[y1,y2])

# Now, repeat for a different bandwidth value

h  = 1.             # bandwidth; you can change this

z  = (x-x0)/h       # standardized variate

y2 = gaussian_kernel(z)  # the Gaussian kernel

y1 = optimal_kernel(z)   # the optimal kernel

gaddwin("The kernels with a 10-times larger bandwidth")

format double "f4.2"

title = ["Evaluation point is x0 = ",ntoa(x0)," and bandwidth is h = ",ntoa(h)]

gtitle(title)

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

gxaxis("linear",-2.5,4.)

gplot(x,[y1,y2])

# Simulate a nonlinear regression model and estimate it using

# nonparametric regression

N = 200.                 # length of realization

s = 0.2                  # standard deviation of error variates

e = snormal(N,1)*s       # the regression error variates

x = snormal(N,1)              # the explanatory variates

a = 1.                   # parameter; you can change this

b = 5.                   # parameter; you can change this

y = a*tanh(-b*x) + e          # the dependent variable

gaddwin("Simulated y = a*tanh(-b*x) + e")

gplot(y)

h = 0.1                  # set the bandwidth; you can change this

# Draw the scatterplot to show that there is nonlinearity

scatterplot(x,y,0,h,"q","x","y")

# Estimate the nonparametric regression and a linear one for comparison

z = [ones(N,1),x]

[m,u,stat] = np_cmean(y,z,0,h)          # nonparametric

[b,seb,covb,rss,s2,prd,u] = lsq(z,y)    # linear

# Now plot the observations and the two fits

gaddwin("Observations, linear fit and nonparametric fit")

gtitle("Circle = observations, red = linear fit, green = nonparametric fit")

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

gstyle("circle")

gplot(y)

gstyle("solid")

gplot([prd,m])