# NOT RUN {
#Simulate a price series containing both a flash crash and a flash-rally
flashCrashSim = function(iT, dSigma, dPhi, dMu){
vSim = numeric(iT)
vEps = rnorm(iT , sd =dSigma)
vEpsy = rnorm(iT)
vEps[30001:32000] = rnorm(2000 ,sd =seq(from = dSigma ,
to = 2*dSigma , length.out = 2000))
vEps[32001:34000] = rnorm(2000 ,sd =seq(from = 2*dSigma ,
to = dSigma , length.out = 2000))
vEpsy[30001:32000] = -rnorm(2000 ,mean =seq(from = 0,
to = 0.3 , length.out = 2000))
vEpsy[32001:34000] = -rnorm(2000 ,mean =seq(from = 0.3,
to = 0 , length.out = 2000))
vEps[60001:63000] = rnorm(3000,sd = seq(from = dSigma ,
to = 2*dSigma , length.out = 3000))
vEps[63001:66000] = rnorm(3000, sd = seq(from = 2*dSigma ,
to = dSigma, length.out = 3000))
vEpsy[60001:63000] = rnorm(3000 ,mean =seq(from = 0,
to = 0.2 , length.out = 3000))
vEpsy[63001:66000] = rnorm(3000 ,mean =seq(from = 0.2,
to = 0 , length.out = 3000))
vSim[1] = dMu + dPhi *rnorm(1 , mean = dMu , sd = dSigma /sqrt(1-dPhi^2))
for (i in 2:iT) {
vSim[i] = dMu + dPhi * (vSim[(i-1)] - dMu) + vEps[i]
}
vY = exp(vSim/2) * vEpsy
return(vY)
}
#Set parameter values of the simulation
iT = 66500; dSigma = 0.3; dPhi = 0.98; dMu = -10;
#set seed for reproducibility
set.seed(123)
#Simulate the series
y = 500+cumsum(flashCrashSim(iT, dSigma, dPhi, dMu))
#insert an outlier to illustrate robustness.
y[50000] = 500
#Here, the observations are equidistant, but the code can handle unevenly spaced observations.
timestamps = seq(34200 , 57600 , length.out = iT)
testTimes = seq(34260, 57600, 60)
logPrices = log(y)
library("DriftBurstHypothesis")
#calculating drift burst test statistic
DBH = driftBursts(timestamps, logPrices,
testTimes, preAverage = 5, ACLag = -1L,
meanBandwidth = 300L, varianceBandwidth = 900L,
parallelize = FALSE)
#plot test statistic
plot(DBH)
#plot both test statistic and price
plot(DBH, price = y, timestamps = timestamps)
#Plot the mu series
plot(DBH, which = "Mu")
#plot the sigma series
plot(DBH, which = "Sigma")
#plot both
plot(DBH, which = c("Mu", "Sigma"))
#Means of the tstat, sigma, and mu series.
mean(getDB(DBH))
mean(getSigma(DBH))
mean(getMu(DBH))
# }
# NOT RUN {
################## same example with xts object:
library("xts")
#Set parameter values of the simulation
iT = 66500; dSigma = 0.3; dPhi = 0.98; dMu = -10;
#set seed for reproducibility
set.seed(123)
#Simulate the series
y = 500+cumsum(flashCrashSim(iT, dSigma, dPhi, dMu))
#insert an outlier to illustrate robustness.
y[50000] = 500
#Here, the observations are equidistant, but the code can handle unevenly spaced observations.
timestamps = seq(34200 , 57600 , length.out = iT)
startTime = strptime("1970-01-01 00:00:00.0000", "%Y-%m-%d %H:%M:%OS", tz = "GMT")
testTimes = seq(34260, 57600, 60)
price = xts(vY, startTime + timestamps)
DBHxts = driftBursts(timestamps = NULL, log(price),
testTimes, preAverage = 5, ACLag = -1L,
meanBandwidth = 300L, varianceBandwidth = 900L,
parallelize = FALSE)
#check for equality
all.equal(as.numeric(getDB(DBH)), as.numeric(getDB(DBHxts)))
# }
Run the code above in your browser using DataLab