testCCF: Statistic Test for CCF

In chengfeifan/myCCM: calculate CCM

Description

Takes in five parameters, and test whether there is causality between two variables.

Usage

testCCF(x, y, lag.max, alpha, sampleTimes, ...)

Arguments

x	A vector represents a time sequence
y	A vector represents another sequence
lag.max	Max time lag
alpha	Significant level
sampleTimes	The times of permutation sample for second order test

Details

First order test is to test whether the ccfvalue is different from zero. Second Test is to test whether the ccfvalue satisfies the 3-sigma test

Value

When pho < 0, it will return a sentence for no causality. When first order test is unsatisfied, it will return a sentence. When the second order test is unsatisfied, it will return a sentence and list of Q and sampleQ vector. When all the test is passed, it will return Q

Examples

data("BMData")
i=8
u<-BMData[[i]]$u
y<-BMData[[i]]$y
resultw<-testCCF(u,y,lag.max=97,alpha = 0.05,sampleTimes=10000)
resultq<-testCCF(u,y,lag.max=97,alpha = 0.05,sampleTimes=1000)
resultb<-testCCF(u,y,lag.max=97,alpha = 0.05,sampleTimes=100)
datab<-cbind(ccf=resultb[[2]],sampleTime=rep(100,100))
dataq<-cbind(ccf=resultq[[2]],sampleTime=rep(1000,1000))
dataw<-cbind(ccf=resultw[[2]],sampleTime=rep(10000,10000))
data<-rbind(datab,dataq,dataw)
data<-as.data.frame(data)
data$sampleTime<-as.factor(data$sampleTime)
library(plyr)
cdata <- ddply(data, "sampleTime", summarise, ccf.mean=mean(ccf),
               ccf.3sigma=mean(ccf)+3*sd(ccf))
library(ggplot2)
ggplot(data) + geom_density(aes(x = ccf,
                                colour = sampleTime),
                            size=1)+labs(x=expression(paste("C"[pi],""^"CCF")))+
  labs(y=expression(paste("p","(C"[pi],""^"CCF",")")))+
  labs(title = "2nd order system")+
  geom_vline(data=cdata, aes(xintercept=ccf.3sigma,  colour=sampleTime),
             linetype="dashed", size=1)+
  geom_vline(aes(xintercept=result[[1]]),
             linetype="dashed", size=1)

chengfeifan/myCCM documentation built on May 13, 2019, 3:40 p.m.