Advertisement
Not a member of Pastebin yet?
Sign Up,
it unlocks many cool features!
- # http://reakkt.com
- library(fBasics) # stableFit
- k <- 1000
- n <- 250
- ###
- assets <- c("pko","peo")
- market.data <- list()
- for (asset in assets)
- market.data[[asset]] <- read.csv(paste(asset,"_d.csv",sep=""),header=TRUE)
- available.data <- min(unlist(lapply(market.data, nrow)))
- # stableFit:
- # ceny aktywow
- x <- do.call(cbind, lapply(market.data, function(i) tail(i[,"Close"],available.data)))
- d <- cbind(log(x[-1,1]/x[-available.data,1]),log(x[-1,2]/x[-available.data,2]))
- sf1 <- stableFit(d[,1])
- sf2 <- stableFit(d[,2])
- ###
- nv <- c(5,10,25,50,100,250,500,1000)
- nl <- length(nv)
- cor.mc.norm <- matrix(NA,k,nl)
- cor.mc.stable <- matrix(NA,k,nl)
- colnames(cor.mc.norm) <- nv
- colnames(cor.mc.stable) <- nv
- distFun <- function(dist.type,ni) {
- if (dist.type=="norm") x <- rnorm(ni,0,6)
- if (dist.type=="sf1") x <- rstable(ni,sf1@fit$estimate["alpha"],sf1@fit$estimate["beta"],sf1@fit$estimate["gamma"],sf1@fit$estimate["delta"])
- if (dist.type=="sf2") x <- rstable(ni,sf2@fit$estimate["alpha"],sf2@fit$estimate["beta"],sf2@fit$estimate["gamma"],sf2@fit$estimate["delta"])
- return( x )
- }
- for (ni in nv) {
- for (i in 1:k) {
- x <- distFun("sf1",ni)
- y <- distFun("sf2",ni)
- xn <- distFun("norm",ni)
- yn <- distFun("norm",ni)
- cor.mc.stable[i,which(nv==ni)] <- cor(x,y)
- cor.mc.norm[i,which(nv==ni)] <- cor(xn,yn)
- }
- }
- boxplot(cor.mc.stable,main="Random correlations depending on series length",ylab="correlation",xlab="period")
- lines(apply(cor.mc.norm, 2, max),col="Blue",lty="dotted")
- lines(apply(cor.mc.norm, 2, min),col="Blue",lty="dotted")
- #####
- ### biezaca korelacja dla okresu n:
- # ceny aktywow
- x <- do.call(cbind, lapply(market.data, function(i) tail(i[,"Close"],available.data)))
- x <- tail(x,(n+1))
- # zmiany procentowe
- d <- cbind(log(x[-1,1]/x[-n,1]),log(x[-1,2]/x[-n,2]))
- (cor.current <- cor(d[,1],d[,2]))
- ###
- cor.sample <- numeric(k)
- s1 <- round(runif(k,(n+1),available.data))
- x <- do.call(cbind, lapply(market.data, function(i) tail(i[,"Close"],available.data)))
- for (i in 1:k) {
- xs <- x[(s1[i]-n):s1[i],]
- d <- cbind(log(xs[-1,1]/xs[-n,1]),log(xs[-1,2]/xs[-n,2]))
- cor.sample[i] <- cor(d[,1],d[,2])
- }
- density.sample <- density(cor.sample)
- density.random <- density(cor.mc.stable[,paste(n)])
- ylim=range(density.sample$y,density.random$y)
- xlim=range(density.sample$x,density.random$x)
- plot(density.sample,"Historical Correlation Density",xlim=xlim,ylim=ylim)
- lines(density.random,lty="dotted",col="Blue")
- abline(v=cor.current,lty="dashed",col="Red")
- abline(v=mean(cor.sample),lty="dotted",col="Orange")
- mtext(paste("period=",n,
- "mean=",format(mean(cor.sample),digits=4),
- "current=",format(cor.current,digits=4)),
- line=0.5,cex=0.8)
- legend("topleft",legend=c("mean","current"),text.col=c("Orange","Red"))
Advertisement
Add Comment
Please, Sign In to add comment
Advertisement