marino89
|
|
| Total Posts: 16 |
| Joined: Jul 2012 |
| |
|
Hi,
I'm currently working on the Crude Oil WTI forward curves and searching on Google, I found this very interesting article : http://www.utstat.utoronto.ca/sjaimung/papers/VAR-FPCA.pdf.
I've tried to implement the FCPA prestented in the paper with R (fda package) but I didn't found the same results (my eigenvalues functions don't have the same shape and my first beta coefficient (Figure 1 in the paper) doesn't reflect a deterministic trend).
Has anyone of you already worked on that paper?
Thanks! |
|
|
|
 |
|
Yes,
I think I understand you problem,
If you want I can explain you by e-mail ?
It's very long to explain here...
Phynance 471 |
|
|
 |
marino89
|
|
| Total Posts: 16 |
| Joined: Jul 2012 |
| |
|
Hi,
Of course you can explain it to me by email : marino1089gmailcom, it will be very of you, I really don't understand why although using the same data as used in the paper I don't get the same results.
I'm using R fad package which is wonderful and very practical and perhaps my problem comes from the fact that I cannot specify an exponential basis like the one in the paper. Otherwise, I follow exactly the article, I've got more than 3000 forward curves with for each 18 tenors, I fit my data using an exponential basis (specifying the same coefficients), I then fit a VAR(1) model to the coefficients I get from the decomposition of my curves, I then keep the residuals in order to filter and detrend my data and finally carry out a functional PCA.
Here are my problems :
- I don't have a beta 1 time-series that has a deterministic trend as in the paper, my time-series is rather non-stationnary with many gaps.
- The trend regressors in the VAR model is not significant and my residuals are not normally distributed
- and the MAIN problem : my first eigen function explains more than 99.8% of the total variance of my functional sample whereas in the paper the first eigen function accounts for only 92.5% of the total variance because of the previous VAR fitting and detrending.
So to sum up, I'm a little bit lost... It's quite a pity because this article is really really interesting.
thanks a lot for your help!!
See you by email |
|
|
|
|
marino89
|
|
| Total Posts: 16 |
| Joined: Jul 2012 |
| |
|
| I think the problem has something to do with the functional basis (exponential basis) I chose on R . |
|
|
|
| |
|
pj
|
|
| Total Posts: 2720 |
| Joined: Jun 2004 |
| |
|
yet another couple of phrogs  |
вакансия "Программист Психологической службы"
-але! у нас ошибко! не работает бля-бля-бля
-вы хотите об этом поговорить? |
|
|
marino89
|
|
| Total Posts: 16 |
| Joined: Jul 2012 |
| |
|
?
that's a very serious question actually, I've been working on that project for 10 days. I really need some help. |
|
|
|
|
jslade
|
|
| Total Posts: 745 |
| Joined: Feb 2007 |
| |
|
Hey Marino:
This is an interesting paper, but a very complex R package. I had a hard time making it do anything sensible. When I did get it doing sensible things, I got the same result as you did. It seems likely we are forgetting something. Or else, perhaps the exponential basis functions don't work the way we think. Lord knows I can't figure it out from reading the FDA source code.
Maybe starting from the package YieldCurve, you can get something better. |
"Learning, n. The kind of ignorance distinguishing the studious." |
|
|
FatChoi
|
|
| Total Posts: 106 |
| Joined: Feb 2008 |
| |
|
| Possibly stupid question but have you tried with only 5 years data? If you are using 3000 curves you might have the gulf war or the financial crisis in your data set which J&N don't. |
|
|
|
|
marino89
|
|
| Total Posts: 16 |
| Joined: Jul 2012 |
| |
|
Thanks for yout answers. It's that the R package fda is complex, and I'm quite sure the problems comes from the way exponential basis are defined in it as I don't get the same  time-series as the one Jaimungal and Ng get in their paper . Basis exponential are defined as  in the R package whereas J&N's basis exponential functions look like :  . The fact is that we can't specify our own functional basis in the package...
As far data are concerned, I used exactly the same data as in the article (my dataset is made of forward curves from 01/01/2000 to 06/01/2012).
Last point, I've tried to reduce the number of exponential basis function, and do the same work with only 3 functions and I got some results pretty close as the one in the article. Using 3 R exponential basis function, my first PC accounts for 88.8% of the variability of the data and my 3 top eigen functions have better shape. I found out that when adding more exp basis function, the top 3 eigen functions turns out to be periodic functions. Finally, My loading coefficients are not as smooth as in the paper but look like rather as a white noise.
Anyway, plenty of mysteries.
|
|
|
|
FatChoi
|
|
| Total Posts: 106 |
| Joined: Feb 2008 |
| |
|
More back seat driving.
The different bases shouldn't be an obstacle as they span the same space so the fit should be nearly independent and you can transform between them.
Untransformed the PCs will not be the same because the weightings are different. In particular the J&N basis is zero at  for  whereas the fda basis is all 1 there. That will make the  in the fda basis be the sum of the loads in the J&N basis with possibly quite different results. Dividing by  will also produce significantly different results.
If you can't get at the factor loadings in fda to transform them it shouldn't be too hard to fit these by hand then you are just looking at VAR(1) which is easily dealt with in R.
I'm not sure how you get a periodic eigenfunction because there are none bar the constant function in the space of functions defined and the eigenvalues of the residual covariance matrix are real and non-negative. It is possible that the A matrix has complex eigenvalues but likely not exactly periodic behaviour and in any case that governs the weightings of the eigenfunctions not the eigenfunctions themselves.
I would suggest if doing this myself that it would be worth thinking through both the basis and the weightings used in the fit. For my money the dimension reduction achieved through projection onto the function space is as important as the two dimensions shed in the PCA step.
|
|
|
|
|
jslade
|
|
| Total Posts: 745 |
| Joined: Feb 2007 |
| |
|
Marino:
It turns out, it is pretty easy to alter the basis functions if you want to. It doesn't produce satisfying results, but it begins to make a little more sense.
expon2 <- function (x, ratevec=1, nderiv=0) {
x <- as.vector(x)
n <- length(x)
nrate <- length(as.vector(ratevec))
expval <- matrix(0,n,nrate)
for (irate in 1:nrate) {
rate <- ratevec[irate]
if(rate==0) {
if(nderiv==0) {
expval[,irate] <- rep(1,n)
} else {
expval[,irate] <- rep(0,n)
}
} else {
if(nderiv==0) {
expval[,irate] <- (1 - exp(rate*x))/abs(rate)
} else {
expval[,irate] <- -(rate)^(nderiv) * exp(rate*x)/abs(rate)
}
}
}
return (expval)
}
unlockBinding("expon", as.environment("package:fda"))
assignInNamespace("expon", expon2, ns="fda", envir=as.environment("package:fda"))
assign("expon", expon2, envir=as.environment("package:fda"))
lockBinding("expon", as.environment("package:fda"))
|
"Learning, n. The kind of ignorance distinguishing the studious." |
|
|
marino89
|
|
| Total Posts: 16 |
| Joined: Jul 2012 |
| |
|
@ jslade : thanks a lot for your answer, I've tried your code and it works perfectly (I didn't know packages could be modified). Nevertheless, I don't get the same results as J&N when going on with the FPCA (the first eigen function accounts for about 100% of the total variability).
Here is my code : forward_matrix is a 2d-array containing my forward curves (one per column).
library(fda)
library(vars)
nbasis<-5
exp_basis<-create.exponential.basis(c(0,1),nbasis,c(0,-4,-2,-1,-0.2))
# Study period :
study<- which(annees>=2008 & annees<=2011)
# Décomposition des courbes forward dans la base fonctionnelle choisie :
functional_form<-Data2fd(seq(0,1,len=18),as.array(forward_matrix[,study]),exp_basis)
ts.plot(projection[1,])
decomposition<-functional_form$coefs
par(mfrow=c(5,1))
ts.plot(forward_matrix[1,study],main="Spot Crude Oil 2009-2011",ylab="Spot")
ts.plot(Data2fd(seq(0,1,len=18),forward_matrix[,study],exp_base)$coef[1,],main="1st coefficient",ylab="Coef 1", col="blue")
ts.plot(Data2fd(seq(0,1,len=18),forward_matrix[,study],exp_base)$coef[2,],main="2nd coefficient",ylab="Coef 2", col="dark red")
ts.plot(Data2fd(seq(0,1,len=18),forward_matrix[,study],exp_base)$coef[3,],main="3rd coefficient",ylab="Coef 3", col="dark green")
ts.plot(Data2fd(seq(0,1,len=18),forward_matrix[,study],exp_base)$coef[4,],main="4th coefficient",ylab="Coef 4", col="dark blue ")
# VAR fitting : estimation and residuals
model_VAR<-VAR(ts(t(decomposition)),p=1,type="both")
residuals<-resid(model_VAR)
plot(ts(residals))
# FPCA carried on the residuals ts :
FPCA<-pca.fd(Data2fd(seq(0,1,len=nbasis),t(residus),exp_basis),nharm=3)
FPCA$values
FPCA$varprop
par(mfrow=c(1,1))
plot(FPCA$harmonics)
loadings<-FPCA$scores.
|
|
|
|
|
FatChoi
|
|
| Total Posts: 106 |
| Joined: Feb 2008 |
| |
|
How does it look with
study<- which(annees>=2002 & annees<=2006)
?
|
|
|
|
marino89
|
|
| Total Posts: 16 |
| Joined: Jul 2012 |
| |
|
Here are plotted the data curve fitting coefficient (the basis used is the one jslade coded). The  time-series is very closed to the one in J&N paper but the problem is that I don't have loadings coefficients as the one plotted in the paper (figure 5) nor the same top 3 eigen function (figure 6).
Here is the R output of the FPCA :
> FPCA$values
[1] 1.293868e+10 1.561545e+00 4.855210e-02 1.268641e-03 4.979731e-05
and the variability explained by the top 3 eigen functions :
> FPCA2$varprop
[1] 1.000000e+00 4.278497e-10 2.502134e-11
I'm supposed to find that the 1st eigen function accounts for 92.5% of the total variability!
|
|
|
|
|
FatChoi
|
|
| Total Posts: 106 |
| Joined: Feb 2008 |
| |
|
| It looks like the coefficients 2, 3 & 4 are quite precisely related to each other so I would look at the projection onto the function space. Should be possible to compare decomposition directly to verify. |
|
|
|
marino89
|
|
| Total Posts: 16 |
| Joined: Jul 2012 |
| |
|
| And then, I don't see the point... |
|
|
|
|
jslade
|
|
| Total Posts: 745 |
| Joined: Feb 2007 |
| |
|
Hey Marino:
How come you only got 18 contracts? I have 60. I also noticed you didn't normalize the exponents to the daily decay time; perhaps you are using monthly data? Or maybe I am misunderstanding something.
I also got encouraging VAR models, which is mostly what I am looking for, but I can't say the output of the rest makes much sense to me. This is perhaps because the method of fitting the coeffecients is different in the FDA package than it is in the paper. Possibly, fiddling with smooth.basis() would help.
Do let us know if you figure out how to use this thing. I don't care much about oil, but such tricks might be useful elsewhere. |
"Learning, n. The kind of ignorance distinguishing the studious." |
|
|
marino89
|
|
| Total Posts: 16 |
| Joined: Jul 2012 |
| |
|
Hi jslade,
I've got daily data from 01/01/2000 to 06/01/2012 and I chose to limit my forward curves to 18 tenors (i could have got more data).
What do you mean by "normalize the exponents to the daily decay time"?
I'll try to fit the coefficients with smooth.basis(). We keep in touch. thks |
|
|
|
|
FatChoi
|
|
| Total Posts: 106 |
| Joined: Feb 2008 |
| |
|
Hi Marino, when you call
exp_basis<-create.exponential.basis(c(0,1),nbasis,c(0,-4,-2,-1,-0.2))
and
functional_form<-Data2fd(seq(0,1,len=18),as.array(forward_matrix[,study]),exp_basis)
you are assuming that time is measured in units of 1/(time to maturity of the 18th contract) and that the contracts are evenly spaced at fixed distances in the future rather than at fixed dates. Today I think the current front contract in WTI is the September '12 contract which has last trade day on August 21 which is 27 days ahead (I'm not an oil person so that may not be quite the date of interest but the idea is the same) so you should be evaluating your basis function at 27 days for todays observation of the first contract and at 26 days for tomorrow.
|
|
|
|
