The code below shows the computation of new coefficient as well as the simulations outlined in the article (i.e. vis-a-vis the ordinary correlation):
Code: Select all
'Reference: Erdem, O., Ceyhan, E. and Varli, Y. (2014). A new correlation coefficient for bivariate time-series data, Physica A, 414, pp. 274–284
mode quiet
!nobs = 100 '# observations
!nsim = 10000 '# simulations
wfcreate u !nobs 'create an unstructured workfile
!case= 1 '(1 -> Stationary, 2 -> Random walk, 3- > Cointegrated)
matrix(!nsim,2) compare 'store the results
'generate null series
series x
series y
'Run the simulations
smpl @first+1 @last
!rho = 0
for !i = 1 to !nsim
if !case = 1 then
x(1) = @rnorm 'random initialization
y(1) = @rnorm
x = @rnorm
y = @rnorm
else if !case = 2 then
x(1) = @rnorm
y(1) = @rnorm
x = x(-1) + @rnorm
y = y(-1) + @rnorm
else
x(1) = @rnorm
y(1) = @rnorm
!a=0.8
x = x(-1) + nrnd
y = !a*x + nrnd
endif
endif
'Compute the ordinary correlation
compare(!i,1) = @cor(x,y)
'Compute the proposed new correlation
call newcor(x,y,!rho)
compare(!i,2) = !rho
'Watch the progress
statusline !i of !nsim
next
'Histograms of simulated correlations
freeze(mode=overwrite,histograms) compare.distplot
show histograms
smpl @all
'Local subroutine that computes the proposed correlation coefficent
subroutine local newcor(series x,series y,scalar rho)
series dx = d(x)
series dy = d(y)
series dxy = dx*dy
group g.add dx dy
local smpl
smpl @all if @rnas(g)=0 'data may include missing values
!Ax = @sqrt(@sumsq(dx)/@obssmpl)
!Ay = @sqrt(@sumsq(dy)/@obssmpl)
!Axy = @sum(dxy)/@obssmpl
rho = !Axy/(!Ax*!Ay)
endsub