Solved – Code for detrended cross-correlation in R

correlationcross correlationrtime series

I want to code for Detrended Cross Correlation in R for time-series data but I'm still stuck. I don't know why the coefficient is not in range -1 : 1. I try to write following these equation below

Measuring correlations between non-stationary series with DCCA coefficient

Detrened cross-correlation coefficient is calculated as detrended covariance of two dataset over detrened variance of two integrated series

enter image description here (Equation 1)

For time-series {xt}, use integrated series profile

enter image description here (Equation 2)

where the data must be detrended by local trend in box of size s

enter image description here (Equation 3)

enter image description here(Equation 4)

The X_hat is linear fit value evaluated by least square method

Detrended covariance of two profiles

enter image description here (Equation 5)

Average the covariance over all boxes

enter image description here (Equation 6)

## data_1
    x= c(-1.042061,-0.669056,-0.685977,-0.067925,0.808380,1.385235,1.455245,0.540762 ,0.139570,-1.038133,0.080121,-0.102159,-0.068675,0.515445,0.600459,0.655325,0.610604,0.482337,0.079108,-0.118951,-0.050178,0.007500,-0.200622)
    ## data_2
    y= c(-2.368030,-2.607095,-1.277660,0.301499,1.346982,1.885968,1.765950,1.242890,-0.464786,0.186658,-0.036450,-0.396513,-0.157115,-0.012962,0.378752,-0.151658,0.774253,0.646541,0.311877,-0.694177,-0.412918,-0.338630,0.276635)
    ## window size = 6
    k=6
    DCCA_CC=function(x,y,k){
      ## calculate cumulative sum profile of all t
    xx<- cumsum(x - mean(x))  ## Equation 2
    yy<- cumsum(y - mean(y))  ## Equation 2

      ## Divide in to overlapping boxes of size k

  slide_win_xx = mat_sliding_window(xx,k)
  slide_win_yy = mat_sliding_window(yy,k)
  ## calculate linear fit value in each box 
  x_hat = t(apply(slide_win_xx,1,function(n) (lm(n~seq(1:length(n)))$fitted.values)))
  y_hat = t(apply(slide_win_yy,1,function(n) (lm(n~seq(1:length(n)))$fitted.values)))

##  Get detrend variance in each box with linear fit value (detrend by local trend).
  F2_dfa_x = c()
  F2_dfa_y = c()
  for(i in 1:nrow(x_hat)){
 ## Equation 4
    F2_dfa_x = c(F2_dfa_x,mean((xx[i:(i+k-1)]-x_hat[i,])^2))
  }
  for(i in 1:nrow(y_hat)){
## Equation 4
    F2_dfa_y = c(F2_dfa_y,mean((yy[i:(i+k-1)]-y_hat[i,])^2))
  }
  ## Average detrend variance over all boxes to obtain fluctuation
  F2_dfa_x = mean(F2_dfa_x) ## Equation 3
  F2_dfa_y = mean(F2_dfa_y) ## Equation 3

  ## Get detrended covariance of two profile
  F2_dcca = c()
  for(i in 1:nrow(x_hat)){
  ## Equation 5
    F2_dcca = c(F2_dcca,mean((xx[i:(i+k-1)]-x_hat[i,]) * (yy[i:(i+k-1)]-y_hat[i,]) ))
  }

## Equation 6
  F2_dcca = mean(F2_dcca)

## Calculate correlation coefficient 
  rho = F2_dcca / (F2_dfa_x * F2_dfa_y) ## Equation 1
  return(rho)
}

mat_sliding_window = function(xx,k){
## Function to generate boxes given dataset(xx) and box size (k)
  slide_mat=c()
  for (i in 1:(length(xx)-k+1)){
    slide_mat = rbind(slide_mat,xx[i:(i+k-1)] )
  }
  return(slide_mat)
}

print(DCCA_CC(x,y,k)) ##This give me 3.392302

I'm not sure if something wrong in integrated profile.

Best Answer

Alright, I finally find my mistake. It is before the final line of the DCCA_CC function.

instead of 

rho = F2_dcca / (F2_dfa_x * F2_dfa_y)

It must be

rho = F2_dcca / sqrt(F2_dfa_x * F2_dfa_y)
Related Question