IBEX (^IBEX)

We will obtain by descomponer ( AMPLITUDE DOMAIN-FEQUENCY REGRESSION - Francisco Parra - ref:https://cran.r-project.org/web/packages/descomponer/index.html)) a smooth time series for IBEX annual, quarterly and monthly frequencies to 2012-2016 period (Yahoo finances). Since the series has 1260 data we will represent the cycle trends occurring 5 times in the time period of the data (\(256=\frac{1280}{5}\)), 20 times in the time period (\(64=\frac{1280}{20}\)), and 60 times in the time period (\(21.333=\frac{1280}{60}\))

library(descomponer)
library(quantmod)
getSymbols('^IBEX', from = '1990-01-01')
## [1] "IBEX"
t1 <- Sys.time()
ibex.des.a=descomponer(Ad(IBEX)['2012::2016'],256,1)
ibex.des.t=descomponer(Ad(IBEX)['2012::2016'],60,1)
ibex.des.m=descomponer(Ad(IBEX)['2012::2016'],21,1)
plot(ts(Ad(IBEX)['2012::2017'],frequency = 256, start = c(2012, 1)))
lines(ts(ibex.des.a$datos$TD,frequency = 256, start = c(2012, 1)),col='red')
lines(ts(ibex.des.t$datos$TD,frequency = 256, start = c(2012, 1)),col='green') 
lines(ts(ibex.des.m$datos$TD,frequency = 256, start = c(2012, 1)),col='blue') 

t2 <- Sys.time() 
print(t2-t1)
## Time difference of 19.87534 mins

John Ehlers’s Autocorrelation Periodogram

JThe Ehlers Autocorrelation Periodogram is a technical indicator proposed by John F. Ehlers; it attempts to reveal dominant cycles in market data and measure their amplitude. (Autocorrelation Periodogram - John F. Ehlers - Technical Analysis of Stocks and Commodities, Sep 2016 ref: http://technical.traders.com/content/TTlink.asp?mo=07&yr=2016).

The Autocorrelation Periodogram takes Autocorrelation of price data filtered by Ehlers Roofing Filter. The Autocorrelation function compares the data with itself a certain period back; its value of +1 would signify the perfect correlation and -1, the perfect anti-correlation. The Ehlers Autocorrelation Periodogram dynamically set a period between a minimum (21) and a maximum (256) period length.

The R function are in: https://quantstrattrader.wordpress.com/2017/02/15/ehlerss-autocorrelation-periodogram/

AGC <- function(loCutoff = 10, hiCutoff = 48, slope = 1.5) {      accSlope = -slope #acceptableSlope = 1.5 dB   ratio = 10 ^ (accSlope / 20)   if ((hiCutoff - loCutoff) > 0)
ratio = 10 ^ (accSlope / 20)
    factor <-  ratio ^ (2 / (hiCutoff - loCutoff));
  return (factor)
}

 
# Autocorrelation Periodogram - John F. Ehlers - Technical Analysis of Stocks and Commodities, Sep 2016
# ref: http://technical.traders.com/content/TTlink.asp?mo=07&yr=2016
autocorrPeriodogram <- function(x, period1 = 10, period2 = 48, avgLength = 3) {
# high pass filter
# NOTE: changed from previous version in "Cycle Analytics for Traders"
alpha1 <- (1 - sin(2 * pi / period2)) / cos(2 * pi / period2)
hp <- (1 - alpha1) / 2 * (x - lag(x))
hp <- hp[-1]
hp <- filter(hp, alpha1, method = "recursive")
hp <- c(NA, hp)
hp <- xts(hp, order.by = index(x))
# super smoother
a1 <- exp(-sqrt(2) * pi / period1)
b1 <- 2 * a1 * cos(sqrt(2) * pi / period1)
c2 <- b1
c3 <- -a1 * a1
c1 <- 1 - c2 - c3
filt <- c1 * (hp + lag(hp)) / 2
leadNAs <- sum(is.na(filt))
filt <- filt[-c(1: leadNAs)]
filt <- filter(filt, c(c2, c3), method = "recursive")
filt <- c(rep(NA, leadNAs), filt)
filt <- xts(filt, order.by = index(x))
# Pearson correlation for each value of lag
autocorr <- matrix(0, period2, length(filt))
for (lag in 2: period2) {
# Set the average length as M
if (avgLength == 0) M <- lag
else M <- avgLength
autocorr[lag, ] <- runCor(filt, lag(filt, lag), M)
}
autocorr[is.na(autocorr)] <- 0
# Discrete Fourier Transform for each autocorrelation
# The sum of the squares of each value represents relative power at each period
cosinePart <- sinePart <- sqSum <- R <- Pwr <- matrix(0, period2, length(filt))
for (period in period1: period2) {
for (N in 2: period2) {
cosinePart[period, ] = cosinePart[period, ] + autocorr[N, ] * cos(N * 2 * pi / period)
sinePart[period, ] = sinePart[period, ] + autocorr[N, ] * sin(N * 2 * pi / period)
}
sqSum[period, ] = cosinePart[period, ] ^ 2 + sinePart[period, ] ^ 2
R[period, ] <- EMA(sqSum[period, ] ^ 2, ratio = 0.2)
}
R[is.na(R)] <- 0
# Normalising Power
#K <- AGC(period1, period2, 1.5) 
maxPwr <- rep(0, length(filt))
for(period in period1: period2) {
for (i in 1: length(filt)) {
if (R[period, i] >= maxPwr[i]) maxPwr[i] <- R[period, i]
#else maxPwr[i] <- K * maxPwr[i] 
}
}
for(period in 2: period2) {
Pwr[period, ] <- R[period, ] / maxPwr
# to enhance resolution Pwr can be elevated to a power
}
# Compute the dominant cycle using the Center of Gravity of the spectrum
Spx <- Sp <- rep(0, length(filter))
for(period in period1: period2) {
Spx <- Spx + period * Pwr[period, ] * (Pwr[period, ] >= 0.5)
Sp <- Sp + Pwr[period, ] * (Pwr[period, ] >= 0.5)
}
dominantCycle <- Spx / Sp
dominantCycle[is.nan(dominantCycle)] <- 0 # 1 en el algoritmo original
dominantCycle <- xts(dominantCycle, order.by=index(x))
dominantCycle <- dominantCycle[dominantCycle > 0] # no esta en el algoritmo original
#heatmap(Pwr, Rowv = NA, Colv = NA, na.rm = TRUE, labCol = "", add.expr = lines(dominantCycle, col = 'blue'))
  return(dominantCycle) #no esta en el algoritmo original
}
# Autocorrelation Reversals (code 8.4)

The IBEX time series and The Autocorrelated Periodogram:

t1 <- Sys.time() 
par(mar=c(1,1,1,1))
out <- autocorrPeriodogram(Ad(IBEX)['2012::2016'], period1 = 30, period2 = 365, avgLength = 3) 
plot(Ad(IBEX)['2012::2016'])

plot(out['2012::2016'])

t2 <- Sys.time()
print(t2-t1)
## Time difference of 35.22101 secs
acpIndicator <- function(x, minPeriod, maxPeriod, indicatorFun = EMA, ...) {
  acpOut <- autocorrPeriodogram(x = x, period1 = minPeriod, period2 = maxPeriod)
  roundedAcpNs <- round(acpOut, 0) # round to the nearest integer
  uniqueVals <- unique(roundedAcpNs) # unique integer values
  out <- xts(rep(NA, length(roundedAcpNs)), order.by=index(roundedAcpNs))
 
  for(i in 1:length(uniqueVals)) { # loop through unique values, compute indicator
    tmp <- indicatorFun(x, n = uniqueVals[i], ...)
    out[roundedAcpNs==uniqueVals[i]] <- tmp[roundedAcpNs==uniqueVals[i]]
  }
  return(out)
}

And the result:

t1 <- Sys.time()
ehlersSMA <- acpIndicator(Ad(IBEX['2012::2016']), 21, 256, indicatorFun = SMA)
plot(Ad(IBEX)['2012::2016'])
lines(ehlersSMA['2012::2016'], col = 'red')

t2 <- Sys.time()
print (t2-t1)
## Time difference of 39.72575 secs

Comparison of results

par(mar=c(1,1,1,1))
plot(ts(Ad(IBEX)['2012::2017'],frequency = 256, start = c(2008, 1)),col=1)
lines(ts(ibex.des.m$datos$TD,frequency = 256, start = c(2008, 1)),col=2)
lines(ts(ibex.des.t$datos$TD,frequency = 256, start = c(2008, 1)),col=3)
lines(ts(ehlersSMA['2012::2016'],frequency = 256, start = c(2008, 1)), col = 4)
legend("bottomleft",lty=1,col=1:4,c("IBEX.sa","IBEX.per.21","IBEX.per.60","IBEX.Elher"))