VAR FPCA analysis paper replication

I've been trying to replicate the following publication: CONSISTENT FUNCTIONAL PCA FOR FINANCIAL TIME-SERIES, Sebastian Jaimungal, Eddie K. H. Ng, 2007

but I havent been able to get the same results estimating the $\beta_{k}$ parameters.

First, I got the data from Bloomberg (CL Comdty contracts) from the period specified in the paper but I don't get how is the $\tau$ parameter specified. I using it as number of years to maturity but with the 70 contracts available in bbg seems like I won't fit the $\tau$ that appears in the publication's charts.

Also, assuming I'm correct, when i fit the $\beta_{k}$ using MLS with the common formula ($\beta_{k} = (X^{T}X)^{-1}Xy)$ I get, for the first beta, a similar chart but with a different scale:

I guess that the parameter calculation is correct, because when i plot the curve against the data it seems very close and the error is quite small, but cant get the same $\beta_{k}$ scale:

could it be maybe that the betas are calibrated with the specifict tenors of $\phi$ function?

code:

# Libraries:

install.packages("matrixcalc",  repos="http://R-Forge.R-project.org")  #If not already installed

library(vars)
library(astsa)
library("XLConnect")

rm(list=ls())
"%^%" <- function(S, power)
{
with(eigen(S), vectors %*% (values^power * t(vectors)))
}

#Parameters

# A_k parametros
ak <- c(4,2,1,0.2)
#Beta Matrix
beta_hat <- matrix(data =NA, nrow = length(data[,1]), ncol = 5)
# Matrix Price
price_hat <- matrix(data = NA, nrow = length(data[,1]), ncol = length(data[1,]))

#Beta

for(i in 1:length(data[,1]))
{
#Clean from NA
count = 0
for(j in seq(1,length(data[1,]),by=2))
{
if( (is.na(data[i,j]) == FALSE) & (as.numeric(tau[i,j]) > 0))
{
count = count + 1
}
}

#Vectors Length
tau_tmp <- c(1:count)*0
price_tmp <- c(1:count)*0

count = 0
for(j in seq(1,length(data[1,]),by=2))
{
if( (is.na(data[i,j]) == FALSE) & (as.numeric(tau[i,j]) > 0))
{
count = count +1
tau_tmp[count] <- as.numeric(tau[i,j])/360
price_tmp[count] <- log(as.numeric(data[i,j]),base =exp(1))
#price_tmp[count] <- as.numeric(data[i,j])
}
}

#Phi Length
phi_tmp <- matrix(nrow =5, ncol = length(tau_tmp))
#price_tmp <- (price_tmp)/sd(price_tmp)

#Phi Values
phi_tmp[1,] <- 1
for(j in 1:length(tau_tmp))
{
for (k  in 1: length(ak))
{
phi_tmp[k+1,j] <- (1- exp(-ak[k]*tau_tmp[j]))/ak[k]
}

}
#Beta (beta_hat = (x'*x)^-1*x'*y)
A <- ginv(t(phi_tmp) %*% phi_tmp) %*% t(phi_tmp)
beta_hat[i,] <- t(price_tmp%*%A)

for(j in 1:length(phi_tmp[1,]))
{
price_hat[i,j] <- sum(beta_hat[i,]*phi_tmp[,j])
}

}

#Plot beta
colnames(beta_hat) <- c("Beta1","Beta2","Beta3","Beta4","Beta5")
plot.ts(beta_hat,base = exp(1) , main = "Betas", xlab = "")

#VAR(1)
fitvar2 = VAR(beta_hat, p=1, type="both")

#Plot E
e_hat <- residuals(fitvar2)
A <- e_hat[,1]
#Lag E
B <- c(1:length(e_hat[,1]))
for (i in 1:length(e_hat[,1]))
{
if (i == 1)
{
B[i] <- NA
}
if (i >1)
{
B[i] <- e_hat[i-1,1]
}
}
#Model against data
b <- 857
ph_plot <- na.omit(price_hat[b,])
pr_plot <- na.omit(log(as.numeric(data[b,]),base=exp(1)))
tau_plot <- na.omit(as.numeric(tau[b,])/360)
plot(ph_plot,col='red')
lines(pr_plot,col='blue')

#FPCA
R <- 0.05
space <- 0.05
#Parameters
phi_tmp <- matrix(data = 0, nrow =5, ncol = length(seq(R,7, by=space)))
N= length(e_hat[,1])
mat = seq(R,7, by=space)
phi_tmp[1,] <- 1
for(j in 1:length(seq(R,7, by=space)))
{
for (k  in 1: length(ak))
{
phi_tmp[k+1,j] <- (1- exp(-ak[k]*mat[j]))/ak[k]
#phi_tmp[k+1,j] <- (1- exp(-ak[k]*j))/ak[k]

}
}

W = phi_tmp %*% t(phi_tmp)
uk = eigen((1/N)*(W %^% 0.5) %*% t(e_hat) %*% e_hat %*% (W %^% 0.5))
z <- (W %^% -0.5) %*% uk$vectors epsilon <- matrix( data = 0, nrow = 5, ncol = length(phi_tmp[1,])) for (i in 1:length(phi_tmp[1,])) { epsilon[,i] = t(z) %*% phi_tmp[,i] } plot(mat,epsilon[1,], ylim = c(-1,1)) lines(mat,epsilon[2,]) lines(mat,epsilon[3,]) print(uk$values/sum(uk\$values))
yt = ginv(z) %*% t(beta_hat)
plot.ts(t(yt))

• Should you post an answer and accept it to indicate the problem has been solved? – Richard Hardy Aug 3 '16 at 18:34
• Mmm is just that something in that publication doesn't seem to fit. I continued with the log transformation and tried to replicate the next steps and got completely different results, so i'm still trying to guess if there is a mistake in my procedure or maybe the data that the author is different than mine. I'll update the post with the code.... – Jose Pedro Melo Aug 3 '16 at 23:21
• any ideas? or should i ask in other fórum? – Jose Pedro Melo Sep 4 '16 at 19:35
• @JosePedroMelo, could you provide a work link on the paper? – Nick Oct 2 '16 at 1:48
• @Nick utstat.utoronto.ca/sjaimung/papers/VAR-FPCA.pdf thats the link. – Jose Pedro Melo Oct 6 '16 at 12:55