Covariance Prediction Using Kernel Regression

Introduction

This is the home page of the article "Estimation of Conditional Covariance Using Kernel Regression", written by Jussi Klemelä.

Supplementary Material

Supplementary material about Dow Jones 30 stocks is in covapred-supple-dj30.pdf.

Supplementary material about evaluation of the predictors with Dow Jones 30 stocks data using sum of squared prediction errors is in covapred-supple-dj30-eval.pdf.

Supplementary material about evaluation of the predictors with Dow Jones 30 stocks data using Markowith portfolio selection is in covapred-supple-dj30-wealth.pdf.

Reproducing the computations

Install software

source("http://jklm.fi/denpro/denpro.R")
source("http://jklm.fi/regpro/regpro.R")
source("http://jklm.fi/finatool/finatool.R")
source("http://jklm.fi/art/covapred/covapred.R")

Read Daily SP 500 and Nasdaq-100 Data

file<-"http://jussiklemela.com/art/covapred/sp500-ndx100.csv"
data<-read.csv(file=file)

n<-dim(data)[1]-1
d<-dim(data)[2]
r<-matrix(0,n,d)
for (j in 1:d) for (i in 1:n) r[i,j]<-(data[i+1,j]-data[i,j])/data[i,j]

times0<-matrix(0,dim(data)[1],1)
delta<-1/251.63
alku<-1985+9/12
for (i in 1:length(times0)) times0[i]<-alku+(i-1)*delta
times<-times0[2:length(times0)]

plot(r,xlab="S&P 500",ylab="Nasdaq-100")

plot(times,r[,1])
	
plot(times,r[,2])	

Read Daily SP 500 and 10 Year Bond Data

file<-"http://jussiklemela.com/art/covapred/sp500-dgs10.csv"
data<-read.csv(file=file)

r<-data
n<-dim(r)[1]  

times<-matrix(0,dim(data)[1],1)
delta<-1/251.63
alku<-1962+1/12
for (i in 1:length(times)) times[i]<-alku+(i-1)*delta

plot(r,xlab="S&P 500",ylab="10 year bond")

plot(times,r[,1])
	
plot(times,r[,2])	

Estimation in the DCC Model

SP 500 and Nasdaq-100

file<-"http://jussiklemela.com/art/covapred/sp500-ndx100.csv"
data<-read.csv(file=file)
n<-dim(data)[1]-1
d<-dim(data)[2]
r<-matrix(0,n,d)
for (j in 1:d) for (i in 1:n) r[i,j]<-(data[i+1,j]-data[i,j])/data[i,j]

library("rmgarch")

t0<-250
parat<-matrix(0,n,8)
sigmat<-matrix(0,n,4)
cort<-matrix(0,n,1)
cort.eps<-matrix(0,n,1)
sigmat1<-matrix(0,n,1)
sigmat2<-matrix(0,n,1)
       
step<-50
ota<-seq(t0+1,n,step)

for (t in 2:t0){
   now<-r[1:t,]
   sigmat[t,1]<-var(now[,1])
   sigmat[t,4]<-var(now[,2])		       
   sigmat[t,2]<-cov(now[,1],now[,2])		       
   sigmat[t,3]<-cov(now[,1],now[,2])		       
   sigmat1[t]<-sd(now[,1])
   sigmat2[t]<-sd(now[,2])
   cort[t]<-cor(now[,1],now[,2])
   cort.eps[t]<-cort[t]  #cor(now[,1]/sigmat1,now[2:t,2]/sigmat2)
}
sigmat[1,]<-sigmat[2,]
cort[1]<-cort[2]
cort.eps[1]<-cort.eps[2]
sigmat1[1]<-sigmat1[2]
sigmat2[1]<-sigmat2[2]			  
# 
q<-cov(now[,1]/sigmat1[1:t],now[,2]/sigmat2[1:t])   #sigmat[t,2]
q11<-var(now[,1]/sigmat1[1:t])   
q22<-var(now[,2]/sigmat2[1:t])   
for (i in 1:length(ota)){   # for (t in (t0+1):n){
    t<-ota[i]
    if (i==length(ota)) up<-n else up<-ota[i+1]-1		  
    #
    now<-ts(r[1:t,])
    dc<-dccfit(dspec,now)
    #
    parat[t:up,1]<-as.numeric(coef(dc)[1])  #alpha0  
    parat[t:up,2]<-as.numeric(coef(dc)[2])  #alpha1 
    parat[t:up,3]<-as.numeric(coef(dc)[3])  #beta
    parat[t:up,4]<-as.numeric(coef(dc)[4])  #alpha0
    parat[t:up,5]<-as.numeric(coef(dc)[5])  #alpha1
    parat[t:up,6]<-as.numeric(coef(dc)[6])  #beta
    parat[t:up,7]<-as.numeric(coef(dc)[7])  #e
    parat[t:up,8]<-as.numeric(coef(dc)[8])  #f
}
for (t in (t0+1):n){
    alpha0<-parat[t,1]
    alpha1<-parat[t,2]
    beta<-parat[t,3]
    sigmat1[t]<-sqrt(alpha0+alpha1*r[t,1]^2+beta*sigmat1[t-1]^2)
    alpha0<-parat[t,4]
    alpha1<-parat[t,5]
    beta<-parat[t,6]
    sigmat2[t]<-sqrt(alpha0+alpha1*r[t,2]^2+beta*sigmat2[t-1]^2)
    #
    e1<-parat[t,7]
    f<-parat[t,8]
    e0<-(1-e1-f)*cov(r[1:t,1]/sigmat1[1:t],r[1:t,2]/sigmat2[1:t])
    q<-e0+e1*(r[t,1]/sigmat1[t])*(r[t,2]/sigmat2[t])+f*q
    q11<-e0+e1*(r[t,1]/sigmat1[t])^2+f*q11
    q22<-e0+e1*(r[t,2]/sigmat2[t])^2+f*q22
    #
    cort[t]<-q/sqrt(q11*q22)
    cort.eps[t]<-q
    sigmat[t,1]<-sigmat1[t]^2    
    sigmat[t,2]<-cort[t]*sigmat1[t]*sigmat2[t]			      
    sigmat[t,3]<-cort[t]*sigmat1[t]*sigmat2[t]
    sigmat[t,4]<-sigmat2[t]^2		      
}

Estimation in the BEKK Model

SP 500 and Nasdaq-100

file<-"http://jussiklemela.com/art/covapred/sp500-ndx100.csv"
data<-read.csv(file=file)
n<-dim(data)[1]-1
d<-dim(data)[2]
r<-matrix(0,n,d)
for (j in 1:d) for (i in 1:n) r[i,j]<-(data[i+1,j]-data[i,j])/data[i,j]

library("mgarchBEKK")

t0<-250
parat<-matrix(0,n,11)
sigmat<-matrix(0,n,4)

step<-20
ota<-seq(t0+1,n,step)

for (t in 2:t0){
   now<-r[1:t,]
   sigmat[t,1]<-var(now[,1])
   sigmat[t,4]<-var(now[,2])		       
   sigmat[t,2]<-cov(now[,1],now[,2])		       
   sigmat[t,3]<-cov(now[,1],now[,2])		       
}
S<-matrix(0,2,2)
S[1,1]<-sigmat[t,1]
S[1,2]<-sigmat[t,2]
S[2,1]<-sigmat[t,3]
S[2,2]<-sigmat[t,4]
for (i in 1:length(ota)){   # for (t in (t0+1):n){
    t<-ota[i]
    if (i==length(ota)) up<-n else up<-ota[i+1]-1		  
    #
    now<-ts(r[1:t,])
    bek<-BEKK(now)   #-mean(y))  
    C<-bek$est.params$"1"
    A<-bek$est.params$"2"
    B<-bek$est.params$"3"
    #
    parat[t:up,1]<-C[1,1]
    parat[t:up,2]<-C[1,2]
    parat[t:up,3]<-C[2,2]
    parat[t:up,4]<-A[1,1]
    parat[t:up,5]<-A[1,2]
    parat[t:up,6]<-A[2,1]
    parat[t:up,7]<-A[2,2]
    parat[t:up,8]<-B[1,1]
    parat[t:up,9]<-B[1,2]
    parat[t:up,10]<-B[2,1]
    parat[t:up,11]<-B[2,2]
}
for (t in (t0+1):n){
    C<-matrix(0,2,2)
    A<-matrix(0,2,2)
    B<-matrix(0,2,2)
    C[1,1]<-parat[t,1]
    C[1,2]<-parat[t,2]
    C[2,1]<-0
    C[2,2]<-parat[t,3]
    A[1,1]<-parat[t,4]
    A[1,2]<-parat[t,5]
    A[2,1]<-parat[t,6]
    A[2,2]<-parat[t,7]
    B[1,1]<-parat[t,8]
    B[1,2]<-parat[t,9]
    B[2,1]<-parat[t,10]
    B[2,2]<-parat[t,11]
    #
    e<-matrix(r[t,],2,1)
    E<-e%*%t(e)		      
    S<-t(C)%*%C+t(A)%*%E%*%A+t(B)%*%S%*%B
    sigmat[t,1]<-S[1,1]			      
    sigmat[t,2]<-S[1,2]			      
    sigmat[t,3]<-S[2,1]		
    sigmat[t,4]<-S[2,2]			      
    #sigmat[t]<-sqrt(alpha0+alpha1*r[t]^2+beta*sigmat[t-1]^2)
}
sigmat[1,]<-sigmat[2,]
				       

Prediction Using Moving Averages

SP 500 and Nasdaq-100

file<-"http://jussiklemela.com/art/covapred/sp500-ndx100.csv"
data<-read.csv(file=file)
n<-dim(data)[1]-1
d<-dim(data)[2]
r<-matrix(0,n,d)
for (j in 1:d) for (i in 1:n) r[i,j]<-(data[i+1,j]-data[i,j])/data[i,j]

# Exponentially weighted moving average				       

t0<-250
hoot<-c(8,10,12,14,16,18,20)
ce<-covapred.ewma(r,hoot,t0=t0)

# Exponentially weighted moving average	correlation predictor

file<-"http://jussiklemela.com/art/covapred/garch-sp500-ndx100.csv"
sigmat<-read.csv(file=file)
sigmat<-sigmat[,c(2,3)]

t0<-250	 
hoot<-c(10,50,200,400,800,1500,3000,6000)
cec<-covapred.ewma.cor(r,hoot,sigmat,t0=t0)

# Comparison of cumulative sums of squared prediction errors      

t00<-1000      
cspe.ma<-matrix(0,n,1)
cspe<-matrix(0,n,1)
for (i in t00:n){
    ine<-which.min(cspe.ma[i-1,])
    cspe.ma[i]<-cspe.ma[i-1]+ce$spe[i,ine]
    #
    ine<-which.min(cspe[i-1,])
    cspe[i]<-cspe[i-1]+cec$spe[i,ine]
}
ero<-cspe-cspe.ma
jako<-cspe/cspe.ma
       
times<-matrix(0,dim(data)[1],1)
delta<-1/251.63
alku<-1985+9/12
for (i in 1:n) times[i]<-alku+(i-1)*delta

ota<-c(1000:(n-10))
plot(times[ota],ero[ota],type="l",lty=1)
segments(0,0,20000,0,col="green")

ota<-c(1000:(n-10))
plot(times[ota],jako[ota],type="l",lty=1,ylim=c(0.98,1.02))
segments(0,1,20000,1,col="green")
     

Prediction Using Kernel Regression

SP 500 and Nasdaq-100

file<-"http://jussiklemela.com/art/covapred/sp500-ndx100.csv"
data<-read.csv(file=file)
n<-dim(data)[1]-1
d<-dim(data)[2]
r<-matrix(0,n,d)
for (j in 1:d) for (i in 1:n) r[i,j]<-(data[i+1,j]-data[i,j])/data[i,j]

# Kernel prediction

t0<-1000				       
hoot<-c(0.6,0.7,0.8,0.9,1,1.1,1.2)    
ck<-covapred.kernel(r,hoot,t0=t0,g1=20,g2=20,g3=10)

# Kernel correlation prediction

t0<-1000     
hoot<-c(0.6,0.7,0.8,0.9,1,1.1,1.2)    
ckc<-covapred.kernel.cor(r,hoot,sigmat,t0=1000,g1=20,g3=10)

# Comparison of cumulative sums of squared prediction errors      

t00<-1000      
cspe.ck<-matrix(0,n,1)
cspe.ckc<-matrix(0,n,1)
for (i in t00:n){
    ine<-which.min(cspe.ck[i-1,])
    cspe.ck[i]<-cspe.ck[i-1]+ck$spe[i,ine]
    #
    ine<-which.min(cspe.ckc[i-1,])
    cspe.ckc[i]<-cspe.ckc[i-1]+ckc$spe[i,ine]
}
ero<-cspe.ckc-cspe.ck
jako<-cspe.ckc/cspe.ck
       
times<-matrix(0,dim(data)[1],1)
delta<-1/251.63
alku<-1985+9/12
for (i in 1:n) times[i]<-alku+(i-1)*delta

ota<-c(1000:(n-10))
plot(times[ota],ero[ota],type="l",lty=1)
segments(0,0,20000,0,col="green")

ota<-c(1000:(n-10))
plot(times[ota],jako[ota],type="l",lty=1) #,ylim=c(0.9,1.2))
segments(0,1,20000,1,col="green")
        				       

July 2017