1

我试图理解 R 中的交叉小波函数,但无法弄清楚如何使用双小波包将相位滞后箭头转换为时间滞后。例如:

require(gamair)
data(cairo)
data_1 <- within(cairo, Date <- as.Date(paste(year, month, day.of.month, sep = "-")))
data_1 <- data_1[,c('Date','temp')]
data_2 <- data_1

# add a lag
n <- nrow(data_1)
nn <- n - 49
data_1 <- data_1[1:nn,]
data_2 <- data_2[50:nrow(data_2),]
data_2[,1] <- data_1[,1]

require(biwavelet)
d1 <- data_1[,c('Date','temp')]
d2 <- data_2[,c('Date','temp')]
xt1 <- xwt(d1,d2)
plot(xt1, plot.phase = TRUE)

在此处输入图像描述 在此处输入图像描述

这是我的两个时间序列。两者相同,但一个落后于另一个。箭头表示 45 度的相位角——显然向下或向上表示 90 度(同相或异相),所以我的解释是我看到的是 45 度的滞后。

我现在如何将其转换为时间滞后,即如何计算这些信号之间的时间滞后?

我在网上读到这只能针对特定波长进行(我认为这意味着一段时间?)。那么,假设我们对 365 周期感兴趣,并且信号之间的时间步长是一天,那么如何计算时间滞后呢?

4

2 回答 2

0

所以我相信你在问如何确定两个时间序列的滞后时间(在这种情况下,你人为地添加了 49 天的滞后时间)。

我不知道有任何软件包可以使此过程成为一步,但由于我们本质上是在处理正弦波,因此一种选择是将波“归零”,然后找到零交叉点。然后,您可以计算第 1 波和第 2 波的零交叉点之间的平均距离。如果您知道测量之间的时间步长,则可以轻松计算滞后时间(在这种情况下,测量步骤之间的时间为一天)。

这是我用来完成此操作的代码:

#smooth the data to get rid of the noise that would introduce excess zero crossings)
#subtracted 70 from the temp to introduce a "zero" approximately in the middle of the wave
spline1 <- smooth.spline(data_1$Date, y = (data_1$temp - 70), df = 30)
plot(spline1)
#add the smoothed y back into the original data just in case you need it
data_1$temp_smoothed <- spline1$y

#do the same for wave 2
spline2 <- smooth.spline(data_2$Date, y = (data_2$temp - 70), df = 30)
plot(spline2)
data_2$temp_smoothed <- spline2$y

#function for finding zero crossing points, borrowed from the msProcess package
zeroCross <- function(x, slope="positive")
{
  checkVectorType(x,"numeric")
  checkScalarType(slope,"character")
  slope <- match.arg(slope,c("positive","negative"))
  slope <- match.arg(lowerCase(slope), c("positive","negative"))

  ipost  <- ifelse1(slope == "negative", sort(which(c(x, 0) < 0 & c(0, x) > 0)),
  sort(which(c(x, 0) > 0 & c(0, x) < 0)))
  offset <- apply(matrix(abs(x[c(ipost-1, ipost)]), nrow=2, byrow=TRUE), MARGIN=2, order)[1,] - 2
  ipost + offset
}

#find zero crossing points for the two waves
zcross1 <- zeroCross(data_1$temp_smoothed, slope = 'positive')
length(zcross1)
[1] 10

zcross2 <- zeroCross(data_2$temp_smoothed, slope = 'positive')
length(zcross2)
[1] 11

#join the two vectors as a data.frame (using only the first 10 crossing points for wave2 to avoid any issues of mismatched lengths)
zcrossings <- as.data.frame(cbind(zcross1, zcross2[1:10]))

#calculate the mean of the crossing point differences
mean(zcrossings$zcross1 - zcrossings$V2)
[1] 49

我确信有更多雄辩的方法可以解决这个问题,但它应该为您提供所需的信息。

于 2015-08-05T21:10:03.917 回答
0

就我而言,对于半日潮,90 度等于 3 小时(90*12.5 小时/360 = 3.125 小时)。12.5小时是半昼夜时段。因此,对于 45 度等于 -> 45*12.5/360 = 1.56 小时。

因此,在您的情况下:90 度 -> 90*365/360 = 91.25 小时。45 度 -> 45*365/360= 45.625 小时。

于 2018-12-07T07:39:32.047 回答