################################################################# # Breitung rank cointegration test implementation Codes in R ## # R-code by Aviral Kumar Tiwari: aviral.eco@gmail.com ## # 18/03/2014 ## # translated from Gretl codes of Alexander Ludwig ## ################################################################# library(agricolae) #library(dynlm) # SPECIFICATIONS: function(x, y, trend) # time series: x and y # trend: raw data (0), de-meaned data (1), de-trended data (2) breitung.rank.coint <- function(x, y, trends) { trend<-seq(x) if(trends==0) {x <- x} else if(trends==1) {x <- x - mean(x)} else if(trends==2) {x <- lm(x~trend)$resid} if(trends==0) {y <- y} else if(trends==1) {y <- y - mean(y)} else if(trends==2) {y <- lm(y~trend)$resid} x_detr <- x y_detr <- y ranks_x<- rank(x_detr) ranks_y<- rank(y_detr) ranks_diff<- (ranks_x - ranks_y) nrobs<- length(ranks_diff) kappa_T<- (1/nrobs*max(abs(ranks_diff))) ranks_diff_sq<- (ranks_diff^2) zeta_T<- (1/nrobs^3 * sum(ranks_diff_sq)) ranks_diff<-ts(ranks_diff) ranks_diff.l<-lag(ranks_diff,(-1),na.pad=T) d_ranks_diff<- (ranks_diff-ranks_diff.l) sigma2_d<- (1/nrobs^2 * sum(d_ranks_diff^2)) kappa_T_1<-(kappa_T/sqrt(sigma2_d)) zeta_T_1<-(zeta_T/sigma2_d) d_ranks_x<-diff(ranks_x);d_ranks_y<-diff(ranks_y) #method: "pearson", "kendall", "spearman" rho_T<-correl(d_ranks_x, d_ranks_y, method = "pearson")$rho #$stat kappa_T_2<-(kappa_T_1 / (1-0.174*rho_T^2)) zeta_T_2<-(zeta_T_1 / (1-0.462*rho_T)) # Attention: here without squared (according to Liew et al., 2012) uhat <- lm(ranks_y~ranks_x + 0)$resid #Linear Regression Without an Intercept ZETA_T<-(1/nrobs^3 * sum(uhat^2)) uhat<-ts(uhat) uhat.l<-lag(uhat,(-1),na.pad=T) d_uhat<-(uhat - uhat.l) sigma2_u<-(1/nrobs^2 * sum(d_uhat^2)) ZETA_T_1<-(zeta_T/sigma2_u) #kappa_T;zeta_T;kappa_T_1;zeta_T_1;kappa_T_2;zeta_T_2;ZETA_T_1;rho_T results <- c(kappa_T,zeta_T,kappa_T_1,zeta_T_1,kappa_T_2,zeta_T_2,ZETA_T_1,rho_T) results } ## Example x<-rnorm(100,1,2) y<-rnorm(100,2,5) breitung.rank.coint(x,y,2)