15

我有 2 个时间序列,我正在使用ccf它们来查找它们之间的互相关。 ccf(ts1, ts2)列出所有时间滞后的互相关。如何在不手动查看数据的情况下找到导致最大相关性的滞后?

4

4 回答 4

21

发布答案http://r.789695.n4.nabble.com/ccf-function-td2288257.html

Find_Max_CCF<- function(a,b)
{
 d <- ccf(a, b, plot = FALSE)
 cor = d$acf[,,1]
 lag = d$lag[,,1]
 res = data.frame(cor,lag)
 res_max = res[which.max(res$cor),]
 return(res_max)
} 
于 2012-04-30T14:01:05.693 回答
12

我想我会重做上述函数,但让它找到返回原始相关性(正或负)的绝对最大相关性。我还(几乎)最大化了滞后的数量。

Find_Abs_Max_CCF<- function(a,b)
{
 d <- ccf(a, b, plot = FALSE, lag.max = length(a)-5)
 cor = d$acf[,,1]
 abscor = abs(d$acf[,,1])
 lag = d$lag[,,1]
 res = data.frame(cor,lag)
 absres = data.frame(abscor,lag)
 absres_max = res[which.max(absres$abscor),]
 return(absres_max)
}
于 2013-11-21T22:15:02.147 回答
2

因为 3 大于 4,所以我也尝试修改这个函数,这次是从这里实现一个想法:

ccfmax <- function(a, b, e=0)
{
 d <- ccf(a, b, plot = FALSE, lag.max = length(a)/2)
 cor = d$acf[,,1]
 abscor = abs(d$acf[,,1])
 lag = d$lag[,,1]
 res = data.frame(cor, lag)
 absres = data.frame(abscor, lag)
 maxcor = max(absres$abscor)
 absres_max = res[which(absres$abscor >= maxcor-maxcor*e &
                        absres$abscor <= maxcor+maxcor*e),]
 return(absres_max)
}

本质上添加了一个“错误”术语,因此如果有多个接近最大值的值,它们都会被返回,例如:

ayy <- jitter(cos((1:360)/5), 100)
bee <- jitter(sin((1:360)/5), 100)

ccfmax(ayy, bee, 0.02)
           cor lag
348  0.9778319  -8
349  0.9670333  -7
363 -0.9650827   7
364 -0.9763180   8

如果没有e给定值,则将其视为零,并且该函数的行为就像发布的一个nvogen一样。

于 2015-11-04T18:47:31.490 回答
1

我也修改了原始解决方案,以便循环函数并输出与索引 (x) 的字符向量相对应的值:

abs.max.ccf <- function(x,a,b) {
  d <- ccf(a, b, plot=FALSE, lag.max=length(a)-5)
  cor <- d$acf[,,1]
  abscor <- abs(d$acf[,,1])
  lag <- d$lag[,,1]
  abs.cor.max <- abscor[which.max(abscor)]
  abs.cor.max.lag <- lag[which.max(abscor)]
  return(c(x, abs.cor.max, abs.cor.max.lag))
}

我删除了data.frame函数中的部分,因为它不必要地慢。要遍历 a 中的每一列data.frame并将结果返回到 new data.frame,我使用以下方法:

max.ccf <- lapply(colnames(df), function(x) unlist(abs.max.ccf(x, df$y, df[x])))
max.ccf <- data.frame(do.call(rbind, max.ccf))
colnames(max.ccf) <- c('Index','Cor','Lag')
于 2015-10-06T01:12:30.983 回答