[This article was first published on R snippets, and kindly contributed to R-bloggers]. (You can report issue about the content on this page here)
Want to share your content on R-bloggers? click here if you have a blog, or here if you don't.
During ESSA2013 conference I had a discussion about Cont model I have commented a year ago.Want to share your content on R-bloggers? click here if you have a blog, or here if you don't.
In original paper Cont highlights that his model produces distribution of returns characterized by positive excess kurtosis. In this post I want to investigate this assertion.
Cont simulation has three key parameters (check out the paper for details):
- d: standard deviation of new public information incoming to the investors
- l: market depth (this implies that maximum percentage price change is in range [-1/l,1/l])
- s: percentage of investors updating their decision threshold
We want to simulate the model with 1000 agents for 11000 periods (dropping first 1000 as warm up) and calculate kurtosis of observed returns distribution. The selected ranges of d, l and s parameters follow Cont recommendations. Here is the code that generates the data:
library(e1071)< o:p>
library(mgcv)< o:p>
cont.run <- function(burn.in, reps, n, d, l ,s) {< o:p>
tr <- rep(0, n)< o:p>
sig <- rnorm(reps, 0, d)< o:p>
r <- rep(0, reps)< o:p>
for (i in 1:reps) {< o:p>
r[i] <- (sum(sig[i] > tr) – sum(sig[i] < (-tr))) / (l * n)< o:p>
tr[runif(n) < s] <- abs(r[i])< o:p>
}< o:p>
return(r[burn.in:reps])< o:p>
}< o:p>
sim.points <- 60000< o:p>
d <- runif(sim.points,0.001,0.01)< o:p>
l <- runif(sim.points,5,20)< o:p>
s <- runif(sim.points,0.01,0.1)< o:p>
k <- mapply(function(d, l, s) {< o:p>
kurtosis(cont.run(1000, 10000, 1000, d, l ,s))< o:p>
}, d, l, s)Warning! The simulation takes a lot of time so I have uploaded here the file sim_output.txt containing the results. It contains four columns: d, l, s, k, where last column holds calculated kurtosis.
Now let us turn to the code that allows us to visually explore the data. I have saved the generated plots to GIF file using animation package. They show the probability of encountering excess kurtosis in generated data conditional on d, l and s.
library(ggplot2)< o:p>
data.set <- read.table(“sim_output.txt”, head = T,< o:p>
colClasses = rep(“numeric”, 4))< o:p>
data.set$cd <- cut(data.set$d, seq(0.001, 0.01, len = 10))< o:p>
data.set$cl <- cut(data.set$l, seq(5, 20, len = 16))< o:p>
data.set$cs <- cut(data.set$s, seq(0.01, 0.1, len = 10))< o:p>
data.set$p.excess <- as.numeric(data.set$k > 0)< o:p>
sum.data <- aggregate(p.excess ~ cd + cl + cs, data = data.set, mean)< o:p>
for (i in levels(sum.data$cs)[c(1:9, 8:2)]) {< o:p>
print(ggplot() +< o:p>
geom_point(data = sum.data[sum.data$cs == i,],< o:p>
aes(x = cl, y = cd, colour = p.excess),< o:p>
shape = 15, size = 10) +< o:p>
scale_colour_gradient(low = “blue”, high = “red”) +< o:p>
theme(panel.background =element_blank(),< o:p>
axis.title.x = element_blank(),< o:p>
axis.title.y = element_blank(),< o:p>
text=element_text(colour = “black”, size = 14),< o:p>
axis.text.x=element_text(angle = –90)) +< o:p>
ggtitle(paste(“cs:”, i)))< o:p>
}And here is the result (d is on y-axis and l on x-axis):
We can see that s is the least important parameter and there is a non linear interaction between d and l. But most importantly – in the given range parameters excess kurtosis is not guaranteed to appear. In particular it does not happen when both d and l are large.
To leave a comment for the author, please follow the link and comment on their blog: R snippets.
R-bloggers.com offers daily e-mail updates about R news and tutorials about learning R and many other topics. Click here if you're looking to post or find an R/data-science job.
Want to share your content on R-bloggers? click here if you have a blog, or here if you don't.