This is introduced in my cartoon video of Yukkuri Kaisetsu (Touhou Project fan-art):
Auto-Regressive Conditional Heteroscedasticity (ARCH)
ARCH was developed by an economist Robert F. Engle III having won the 2003 Nobel Memorial Prize in Economic Sciences for its achievement.
Dependent
variable: the variance error terms of the first regression:
The
explanatory variable X can be the
lagged dependent variables and/or the
other variables:
Then, it find
the coefficient γ of the lagged
squared error terms with reference to the log-likelihood:
Generalised Auto-Regressive Conditional Heteroscedasticity (GARCH)
GARCH assumes the variance of the error term symmetrically varies depending the average size of the error terms in pervious time steps. It adds the lagged variance on the explanatory variable of the second regression with reference to the log-likelihood for finding the coefficients γ and δ:
Simulation Data with Python
The following exhibits display the simulation data evaluated using ARCH to illustrate how ARCH functions.
To facilitate the visual representation of this simulation, the most basic form of Engle's ARCH, as introduced in the Wikipedia entry below, has been implemented.
This simplified simulation demonstrates motion prediction for intercepting incoming flying projectiles with erratic movements, resembling the fluctuations of a stock price.
Following is my Python codes:
# Setup Jupyter Notebook in Anaconda dark background https://www.youtube.com/watch?v=_Z6LoHPmoWQ
# pip install jupyterthemes
# !pip install --upgrade jupyterthemes
# importing necessary tools for math and plot
# For mathematics
import math
import random
import numpy as np
import statistics
import statsmodels.api as sm
# For graphing
import matplotlib.pyplot as plt
import matplotlib.gridspec as gridspec
import pandas as pd
# Setting Plot Size
plt.rcParams["figure.figsize"] = (7,3)
# SymPy for writing polynomial formula for determining the slopes
import sympy as smp
# Inline animations in Jupyter
# https://stackoverflow.com/questions/43445103/inline-animations-in-jupyter
from matplotlib.animation import FuncAnimation
import matplotlib.animation as animation
import itertools # for joining lists inside a list
from itertools import count
# Data generation
def Y_(T):
c=np.random.normal(0, 1);
m=4*random.randint(1,6); sign=random.uniform(-1,1);
Y=[c-c/len(T)*t+sign*math.sin(math.pi*(T[t])*m)+np.random.normal(0,1)*(math.cos(130*math.pi*(T[t]))) for t in range(len(T))]
return Y
# Autoregression OLS function
def AutoOLS(DataSet,L): # L: time Lag
DS_0 = DataSet[L:]; DS_lag=[]; N=len(DataSet)
# Adding lagged variables
for l in range(1,L+1):
DS_lag.append(DataSet[L-l:N-l])
ones = np.ones(len(DS_lag[0]))
X = sm.add_constant(np.column_stack((DS_lag[0], ones)))
for ele in DS_lag[1:]:
X = sm.add_constant(np.column_stack((ele, X)))
result = sm.OLS(DS_0, X).fit()
return result
' ######### Generating the reconstruction + prediction parts ######### '
def ReconPred(n):
# Setting the initial parameters
#n=260;
T=np.linspace(0.01,1,n);
P = 5; # Time lag for the first Auto-Regression
Q = 5; # Time lag for the second Auto-Regression
Y=Y_(T); #plt.plot(Y); plt.ylim(-6, 6); plt.show(block=False)
Y_p=Y[:math.floor(n/2)]; Y_f=Y[math.floor(n/2):n]
# Retriving parameters
Outputs1=AutoOLS(Y,P)
alphas=Outputs1.params; #print("alphas",alphas)
Tvalues=Outputs1.tvalues; #print("t-values",Tvalues)
epsilons=Outputs1.resid; #plt.title("Residuals"); plt.plot(epsilons)
#display(Outputs1.summary())
# Analysing the variance = squared residuals
sigma2=[epsilons[s]**2 for s in range(len(epsilons))]
Outputs2=AutoOLS(sigma2,Q)
gammas=Outputs2.params; #print("gammas",alphas)
Pvalues=Outputs2.pvalues; #print("p-values",Pvalues)
#display(Outputs2.summary())
# Reconstruction & Prediction
L=max(P,Q); # samples observation
N=len(Y_p) # Reconstruction period
tau=len(Y_f) # Future prediction period
Y_hat=[]; Y_hat=Y[:L];
Tau=N+L+tau;
for j in range(L,Tau):
Y_hat_j=0;
# Coefficients
for p in range(len(alphas)):
if Tvalues[p]>=2.86: # augmented dickey-fuller T-test
if p>=P:
if len(alphas)-p==1:
Y_hat_j+=alphas[p];
else:
if j-L>=N:
Y_hat_j+=alphas[p]*Y[j-p-N]
else:
Y_hat_j+=alphas[p]*Y[j-p]
# Variances
for q in range(len(gammas)):
if abs(Pvalues[q])<0.05: # simple T-test
if q>len(gammas):
Y_hat_j+=gammas[q]/abs(gammas[q])*np.sqrt(abs(gammas[q]))
if j-L>N:
Y_hat_j+=gammas[q]/abs(gammas[q])*np.sqrt(abs(gammas[q]))*(Y[j-q-1-N]-Y_hat[j-q-1-N])
else:
Y_hat_j+=(gammas[q])*(Y[j-q-1]-Y_hat[j-q-1])
# adding one recon/pred value
Y_hat.append(Y_hat_j)
# plt.plot(Y,"o", color = "lightskyblue")
# plt.plot(Y_hat, "-", color = "lightgreen")
# plt.plot(Y_hat, "--", color = "lightskyblue")
# plt.ylim(-5, 5);
# plt.show(block=False)
return Y, Y_hat
#Y, Y_hat=ReconPred(2*random.randint(30,130))
' ######### Defining the function for displaying animation ######### '
def TwoSparksMoving(F,F_hat):
Sphere=5; Slash=10;
N_=len(F); L=len(F_hat)-len(F)
F_hat
f=np.flip(F); f_hat=np.flip(F_hat);
plt.figure()
plt.rcParams["animation.html"] = "jshtml"
plt.rcParams['figure.dpi'] = 150
plt.ioff()
fig, ax = plt.subplots()
count_ = count(); StepSize=math.ceil(N_/60); NoFrames=math.floor(N_/StepSize)
plt.figure();
def animate(a):
counts=next(count_); Show=counts*StepSize
ax.cla()
# Sphere
ax.plot(f,"--", color='lightskyblue', linewidth=1)
ax.plot(f,"o", color='lightskyblue', linewidth=6)
ax.plot(f[0:math.floor(N_/2)-Show],"o", color='black', linewidth=6)
# Slash
ax.plot(f_hat[0:max(0,Show)], color='lightgreen', linewidth=3)
ax.plot(f_hat[0:max(0,Show-Slash)], color='black', linewidth=15)
ax.set_xlim(0,N_)
ax.set_ylim(-6,6)
ax.set_facecolor('xkcd:black')
Animation=animation.FuncAnimation(fig, animate, frames=NoFrames)
display(Animation)
# Separating figures
plt.show(block=False)
' ######### Animation show ######### '
for o in range(30):
Y, Y_hat=ReconPred(2*random.randint(30,130))
epsilon=[abs(Y[t]-Y_hat[t])/abs(Y[t]) for t in range(math.floor(len(Y)/2))]
if statistics.mean(epsilon)<1:
TwoSparksMoving(Y, Y_hat);
else:
continue
No comments:
Post a Comment