#requires tseries   for the function arma()
## Of course you need to install this library first
library(tseries)

##start function code

eacf.tables=function(x,max.p=15,max.q=15) #A function to produce eacf tables in the Tsay book
# maximum orders to search is default 15 for both
# x is a vector not a time series object
{
pq=max.p+max.q
x1=x-mean(x)
n=length(x1)
coef=matrix(rep(0,pq*pq),pq,pq)
for (i in 1:pq)
{
 md=arma(x1,order=c(i,0), include.intercept=F)
 for (j in 1:i){
 coef[i,j]=md$coef[j]
 }
};
eacf=matrix(rep(0,max.p*max.q),max.p,max.q);
m1=acf(x1,lag.max=max.q,plot=F);
for (j in 1:max.q){eacf[1,j]=m1$acf[j+1]};


for (j in 1:max.q)
{
tmp=coef
 for (i in 1:(pq-j)){
  for (ii in 1:i){
   if (ii == 1) coef[i,ii]= tmp[i+1,ii]+tmp[i+1,i+1]/tmp[i,i]
   if (ii >=2) coef[i,ii]=tmp[i+1,ii]-tmp[i+1,i+1]*tmp[i,ii-1]/tmp[i,i]
                  }
                    };
for (i in 2:max.p){
  y=matrix(rep(0,n),n,1)
  for (it in i:n) {
    w=x1[it]
    for (ii in 1:(i-1)){
      w=w-coef[(i-1),ii]*x1[it-ii]
      }
    y[it,1]=w
    }
  m1=acf(y,lag.max=max.q,plot=F)
  eacf[i,j]=m1$acf[j+1]
  }
};
seacf=matrix(rep(0,max.p*max.q),max.p,max.q);
sd = 2/sqrt(n);
for (i in 1:max.p){
 for (j in 1:max.q){
  if (abs(eacf[i,j]) >= sd) seacf[i,j]=2
 }
}
print('EACF table')
print(eacf)
print(' ')
print('Simplified EACF: 2 denotes significance')
print(seacf)
}
## end function code

## call the function with:
eacf.tables( sp500overlapping.returns)