* ARCH model
* Data = excess

cal 1964 8 12             ;* set calendar to begin 1964 August
all 0 1991:11              ;* data runs until 1991 November
open data excess.txt
data(format=prn,org=obs) / excess

set y = excess
gra 1 ; # y

* ACF plot of the raw data
cor(partial=pacf,qstats,number=36,span=12,dfc=4) y
source(noecho) bjident.src
@bjident y
* Since it appears stationary, we do not need to take a first difference.

com sbcmin = 100000000., sbcp = 0, sbcq = 0
do i=0,5
do j=0,5
 boxjenk(ar=i,ma=j,constant,noprint) y
 compute aic = log(%rss/%nobs) + 2*(%nreg)*(%nobs)
 compute sbc = log(%rss/%nobs) + %nreg*log(%nobs)/(%nobs)
 display 'p ='  i  ' q =' j  '  AIC = ' aic     '  SBC = ' sbc
 if sbc < sbcmin
 {
  compute sbcp = i
  compute sbcq = j
  compute sbcmin = sbc
  }
endo j
endo i
disp
disp 'Selected values of p and q using SBC = ' sbcp  sbcq
disp

* Now, Re-Estimating the ARMA(1,2) model
boxjenk(define=eq1,ar=1,ma=2, constant) y / u
compute sbc = log(%rss/%nobs) + %nreg*log(%nobs)/(%nobs)
display  '  SBC = ' sbc

cor(partial=pacf,qstats,number=36,span=12,dfc=4) u
source(noecho) bjident.src
@bjident u

* Testing for arch effects
set u2 = u*u
lin u2
# constant u2{1 to 5}

* LM stat =  (do it yourself!)

* Estimating GARCH models

source(noecho) garch.src
@GARCH y / res1 Ht1   ; * Choose ar=1, ma=2, p=1, q=1, det=const