/ Published in: R
code inspired by the story presented in http://www-stat.wharton.upenn.edu/~steele/Courses/434/434Context/Co-integration/Murray93DrunkAndDog.pdf
Expand |
Embed | Plain Text
Copy this code and paste it in your HTML
library(tseries) # po.test library(urca) # ca.jo Nmoves <- 1e5 cfreq <- 0.01 # correction frequency cfactor <- c(0.1,0.3,0.6) # correction efficiency; 1,2 - length, 3 - angle; (0-1) ### drunk_path <- matrix(0,Nmoves,2) dog_path <- matrix(0,Nmoves,2) random_walk <- rnorm(Nmoves*2,mean=0,sd=1) for (i in 2:Nmoves) { if (runif(1)>cfreq) { drunk_path[i,] <- drunk_path[i-1,]+rnorm(2,mean=0,sd=1) dog_path[i,] <- dog_path[i-1,]+rnorm(2,mean=0,sd=1) } else { d <- dog_path[i-1,]-drunk_path[i-1,] # delta h <- sqrt(d[1]^2+d[2]^2) # opposite alpha <- atan2(d[2],d[1]) # arc d1 <- h*runif(1,min=cfactor[1],max=cfactor[2]) # reduced delta drunk d2 <- h*runif(1,min=cfactor[1],max=cfactor[2]) # reduced delta dog # a1 <- alpha*runif(1,min=cfactor[3],max=(2-cfactor[3])) # distorted alpha drunk # a2 <- alpha*runif(1,min=cfactor[3],max=(2-cfactor[3])) # distorted alpha dog a1 <- alpha+(1-cfactor[3])*runif(1,-pi,pi) a2 <- alpha+(1-cfactor[3])*runif(1,-pi,pi) d_drunk <- c(d1*cos(a1),d1*sin(a1)) d_dog <- c(d2*cos(a2),d2*sin(a2)) drunk_path[i,] <- drunk_path[i-1,]+d_drunk dog_path[i,] <- dog_path[i-1,]-d_dog # cat("d=",d,"a=",alpha,a1,a2,"\n") } } par(mfrow=c(1,1)) xscope <- c(min(drunk_path[,1],dog_path[,1]),max(drunk_path[,1],dog_path[,1])) yscope <- c(min(drunk_path[,2],dog_path[,2]),max(drunk_path[,2],dog_path[,2])) plot(drunk_path[,1],drunk_path[,2],type="l",xlim=xscope,ylim=yscope,xlab="x",ylab="y") lines(dog_path[,1],dog_path[,2],col="Red") points(drunk_path[Nmoves,1],drunk_path[Nmoves,2],type="p") points(dog_path[Nmoves,1],dog_path[Nmoves,2],type="p") abline(v=drunk_path[Nmoves,1],col="Green") abline(h=drunk_path[Nmoves,2],col="Green") abline(v=dog_path[Nmoves,1],col="Blue") abline(h=dog_path[Nmoves,2],col="Blue") ### distance between pathes dd_path <- drunk_path-dog_path distance <- sqrt(dd_path[,1]^2+dd_path[,2]^2) plot(distance,type="l") ### mean(distance);sd(distance) drunk_path[Nmoves,]; dog_path[Nmoves,] delta <- drunk_path[Nmoves,]-dog_path[Nmoves,] sqrt(delta[1]^2+delta[2]^2) sum(drunk_path-dog_path)/Nmoves ### par(mfrow=c(4,1)) plot(drunk_path[,1],type="l") plot(drunk_path[,2],type="l") plot(dog_path[,1],type="l") plot(dog_path[,2],type="l") cor(drunk_path[,1],dog_path[,1]); cor(drunk_path[,2],dog_path[,2]); cor(drunk_path[,1],drunk_path[,2]) ### Phillips-Ouliaris Cointegration Test x_pathes <- cbind(drunk_path[,1],dog_path[,1]) y_pathes <- cbind(drunk_path[,2],dog_path[,2]) d_pathes <- cbind(drunk_path[,1],drunk_path[,2]) po.test(x_pathes) po.test(y_pathes) po.test(d_pathes) ### Johansen test colnames(x_pathes)<-c("drunk X","dog X") colnames(y_pathes)<-c("drunk Y","dog Y") colnames(d_pathes)<-c("drunk X","drunk Y") summary(ca.jo(x_pathes,type="eigen")) summary(ca.jo(y_pathes,type="eigen")) summary(ca.jo(d_pathes,type="eigen"))
URL: http://www.reakkt.com/2011/10/cointegrated-drunk-and-her-dog.html