Wavemulcor 包 - wave.multiple.cross.correlation 函数 - 替换长度为零
Wavemulcor package - wave.multiple.cross.correlation function - replacement has length zero
我想使用 wavemulcor 包,尤其是 wave.multiple.cross.correlation 函数对我的数据执行小波多重互相关。
我正在按照包手册中的示例进行操作,但改用我的数据。该函数适用于示例数据,但当我尝试使用我的数据时会抛出错误。错误指的是 "replacement has length zero" 但我不确定这到底是什么意思。
我用谷歌搜索了这个错误,但有很多不同功能的相同问题的例子,通常它们都与代码中的循环有关。
然后我用谷歌搜索了如何解决问题并阅读了有关调试的信息。我尝试调试代码,但无法弄清楚它在哪里崩溃,我仍处于学习编码的早期阶段。我认为可能是 wave.multiple.cross.correlation 函数中的这段代码导致了问题:
xy.cor.vec <- matrix(unlist(xy.cor), l, dd)
xy.mulcor <- matrix(0, l, 2 * lm + 1)
YmaxR <- vector("numeric", l)
for (i in 1:l) {
r <- xy.cor.vec[i, ]
P <- diag(d)/2
P[lower.tri(P)] <- r
P <- P + t(P)
Pidiag <- diag(solve(P))
if (is.null(ymaxr)) {
YmaxR[i] <- Pimax <- which.max(Pidiag)
}
有没有其他方法可以确定为什么这不起作用?
实际错误如下:
> Lst <- wave.multiple.cross.correlation(xx, lag.max = NULL, ymaxr = NULL)
Error in YmaxR[i] <- Pimax <- which.max(Pidiag) :
replacement has length zero
我试图追踪原始代码中的所有变量,但我就是想不通。
这是我尝试使用的代码,为了完整起见,here 是一个 link 到 xx 的 dput(),它是我希望使用的变量列表,请参阅xx
详情请见下方代码
library(wavemulcor)
library(readxl)
rm(list = ls()) # clear objects
graphics.off() # close graphics windows
RNC20_30Hourly <- read_excel(RNC20_30Hourly.xlsx")
RNC20_30Hourly <- RNC20_30Hourly[-1]
RNC20_30TS <- ts(RNC20_30Hourly, start = 1, frequency = 23)
wf <- "d4"
J <- 6
lmax <- 36
n <- nrow(RNC20_30TS)
CK0158U09A3.modwt <- modwt(RNC20_30TS[,"CK0158U09A3"], wf, J)
CK0158U09A3.modwt.bw <- brick.wall(CK0158U09A3.modwt, wf)
CK0158U21A1.modwt <- modwt(RNC20_30TS[,"CK0158U21A1"], wf, J)
CK0158U21A1.modwt.bw <- brick.wall(CK0158U21A1.modwt, wf)
CK0158U21A2.modwt <- modwt(RNC20_30TS[,"CK0158U21A2"], wf, J)
CK0158U21A2.modwt.bw <- brick.wall(CK0158U21A2.modwt, wf)
CK0158U21A3.modwt <- modwt(RNC20_30TS[,"CK0158U21A3"], wf, J)
CK0158U21A3.modwt.bw <- brick.wall(CK0158U21A3.modwt, wf)
xx <- list(CK0158U09A3.modwt.bw,
CK0158U21A1.modwt.bw,
CK0158U21A2.modwt.bw,
CK0158U21A3.modwt.bw)
Lst <- wave.multiple.cross.correlation(xx, lag.max = 13, ymaxr = NULL)
CK0158.RTWP.cross.cor <- as.matrix(Lst$xy.mulcor[1:J,])
YmaxR <- Lst$YmaxR
cell.names <- c("CK0158U09A3", "CK0158U21A1", "CK0158U21A2", "CK0158U21A3")
rownames(CK0158.RTWP.cross.cor)<-rownames(CK0158.RTWP.cross.cor,
do.NULL = FALSE, prefix = "Level ")
lags <- length(-lmax:lmax)
lower.ci <- tanh(atanh(CK0158.RTWP.cross.cor) - qnorm(0.975) /
sqrt(matrix(trunc(n/2^(1:J)), nrow=J, ncol=lags)- 3))
upper.ci <- tanh(atanh(CK0158.RTWP.cross.cor) + qnorm(0.975) /
sqrt(matrix(trunc(n/2^(1:J)), nrow=J, ncol=lags)- 3))
par(mfrow=c(3,2), las=1, pty="m", mar=c(2,3,1,0)+.1, oma=c(1.2,1.2,0,0))
for(i in J:1) {
matplot((1:(2*lmax+1)),CK0158.RTWP.cross.cor[i,], type="l", lty=1, ylim=c(-1,1),
xaxt="n", xlab="", ylab="", main=rownames(CK0158.RTWP.cross.cor)[[i]][1])
if(i<3) {axis(side=1, at=seq(1, 2*lmax+1, by=12),
labels=seq(-lmax, lmax, by=12))}
#axis(side=2, at=c(-.2, 0, .5, 1))
lines(lower.ci[i,], lty=1, col=2) ##Add Connected Line Segments to a Plot
lines(upper.ci[i,], lty=1, col=2)
abline(h=0,v=lmax+1) ##Add Straight horiz and vert Lines to a Plot
text(1,1, labels=cell.names[YmaxR[i]], adj=0.25, cex=.8)
}
par(las=0)
mtext('Lag (hours)', side=1, outer=TRUE, adj=0.5)
mtext('Wavelet Multiple Cross-Correlation', side=2, outer=TRUE, adj=0.5)
任何对 resolve/troubleshoot 这个问题的帮助将不胜感激。
仔细检查代码后,我发现了代码中的错误。根据手册,用法如下:
wave.multiple.cross.correlation(xx, lag.max = NULL, ymaxr = NULL)
但在提供的示例中,作者使用了一个名为 lmax 的不同变量,它指定要使用的滞后。
Lst <- wave.multiple.cross.correlation(xx, lmax)
从我上面的示例中可以看出,我为滞后指定了两个不同的参数,lmax = 36 和 lag.max = 13
啊好吧,我只花了 100 个声望!
看来,鉴于您的时间序列的长度,您将最大小波级别 J 设置得太高(代码的第 14 行):
14 焦 <- 6
试试 J <- 4 和 运行
Lst <- wave.multiple.cross.correlation.prueba(xx, lag.max = lmax, ymaxr = NULL)
其中 lmax 是您之前设置为 36(或您选择的任何其他合理值)的最大延迟。
一般来说,最大小波水平的推荐值为J = t运行c(log2(T))-3(为了安全起见),其中T = 时间序列的长度.
事实上,报告的错误是您的原始 xx 列表的结果,其中级别 d6 和 s6 填充了 NaN。
请注意,这不是 wavemulcor 包的一部分,而是小波变换先前计算的一部分
在输入您选择的 wavemulcor 包函数之前。
希望对您有所帮助。
我想使用 wavemulcor 包,尤其是 wave.multiple.cross.correlation 函数对我的数据执行小波多重互相关。
我正在按照包手册中的示例进行操作,但改用我的数据。该函数适用于示例数据,但当我尝试使用我的数据时会抛出错误。错误指的是 "replacement has length zero" 但我不确定这到底是什么意思。
我用谷歌搜索了这个错误,但有很多不同功能的相同问题的例子,通常它们都与代码中的循环有关。
然后我用谷歌搜索了如何解决问题并阅读了有关调试的信息。我尝试调试代码,但无法弄清楚它在哪里崩溃,我仍处于学习编码的早期阶段。我认为可能是 wave.multiple.cross.correlation 函数中的这段代码导致了问题:
xy.cor.vec <- matrix(unlist(xy.cor), l, dd)
xy.mulcor <- matrix(0, l, 2 * lm + 1)
YmaxR <- vector("numeric", l)
for (i in 1:l) {
r <- xy.cor.vec[i, ]
P <- diag(d)/2
P[lower.tri(P)] <- r
P <- P + t(P)
Pidiag <- diag(solve(P))
if (is.null(ymaxr)) {
YmaxR[i] <- Pimax <- which.max(Pidiag)
}
有没有其他方法可以确定为什么这不起作用?
实际错误如下:
> Lst <- wave.multiple.cross.correlation(xx, lag.max = NULL, ymaxr = NULL)
Error in YmaxR[i] <- Pimax <- which.max(Pidiag) :
replacement has length zero
我试图追踪原始代码中的所有变量,但我就是想不通。
这是我尝试使用的代码,为了完整起见,here 是一个 link 到 xx 的 dput(),它是我希望使用的变量列表,请参阅xx
详情请见下方代码library(wavemulcor)
library(readxl)
rm(list = ls()) # clear objects
graphics.off() # close graphics windows
RNC20_30Hourly <- read_excel(RNC20_30Hourly.xlsx")
RNC20_30Hourly <- RNC20_30Hourly[-1]
RNC20_30TS <- ts(RNC20_30Hourly, start = 1, frequency = 23)
wf <- "d4"
J <- 6
lmax <- 36
n <- nrow(RNC20_30TS)
CK0158U09A3.modwt <- modwt(RNC20_30TS[,"CK0158U09A3"], wf, J)
CK0158U09A3.modwt.bw <- brick.wall(CK0158U09A3.modwt, wf)
CK0158U21A1.modwt <- modwt(RNC20_30TS[,"CK0158U21A1"], wf, J)
CK0158U21A1.modwt.bw <- brick.wall(CK0158U21A1.modwt, wf)
CK0158U21A2.modwt <- modwt(RNC20_30TS[,"CK0158U21A2"], wf, J)
CK0158U21A2.modwt.bw <- brick.wall(CK0158U21A2.modwt, wf)
CK0158U21A3.modwt <- modwt(RNC20_30TS[,"CK0158U21A3"], wf, J)
CK0158U21A3.modwt.bw <- brick.wall(CK0158U21A3.modwt, wf)
xx <- list(CK0158U09A3.modwt.bw,
CK0158U21A1.modwt.bw,
CK0158U21A2.modwt.bw,
CK0158U21A3.modwt.bw)
Lst <- wave.multiple.cross.correlation(xx, lag.max = 13, ymaxr = NULL)
CK0158.RTWP.cross.cor <- as.matrix(Lst$xy.mulcor[1:J,])
YmaxR <- Lst$YmaxR
cell.names <- c("CK0158U09A3", "CK0158U21A1", "CK0158U21A2", "CK0158U21A3")
rownames(CK0158.RTWP.cross.cor)<-rownames(CK0158.RTWP.cross.cor,
do.NULL = FALSE, prefix = "Level ")
lags <- length(-lmax:lmax)
lower.ci <- tanh(atanh(CK0158.RTWP.cross.cor) - qnorm(0.975) /
sqrt(matrix(trunc(n/2^(1:J)), nrow=J, ncol=lags)- 3))
upper.ci <- tanh(atanh(CK0158.RTWP.cross.cor) + qnorm(0.975) /
sqrt(matrix(trunc(n/2^(1:J)), nrow=J, ncol=lags)- 3))
par(mfrow=c(3,2), las=1, pty="m", mar=c(2,3,1,0)+.1, oma=c(1.2,1.2,0,0))
for(i in J:1) {
matplot((1:(2*lmax+1)),CK0158.RTWP.cross.cor[i,], type="l", lty=1, ylim=c(-1,1),
xaxt="n", xlab="", ylab="", main=rownames(CK0158.RTWP.cross.cor)[[i]][1])
if(i<3) {axis(side=1, at=seq(1, 2*lmax+1, by=12),
labels=seq(-lmax, lmax, by=12))}
#axis(side=2, at=c(-.2, 0, .5, 1))
lines(lower.ci[i,], lty=1, col=2) ##Add Connected Line Segments to a Plot
lines(upper.ci[i,], lty=1, col=2)
abline(h=0,v=lmax+1) ##Add Straight horiz and vert Lines to a Plot
text(1,1, labels=cell.names[YmaxR[i]], adj=0.25, cex=.8)
}
par(las=0)
mtext('Lag (hours)', side=1, outer=TRUE, adj=0.5)
mtext('Wavelet Multiple Cross-Correlation', side=2, outer=TRUE, adj=0.5)
任何对 resolve/troubleshoot 这个问题的帮助将不胜感激。
仔细检查代码后,我发现了代码中的错误。根据手册,用法如下:
wave.multiple.cross.correlation(xx, lag.max = NULL, ymaxr = NULL)
但在提供的示例中,作者使用了一个名为 lmax 的不同变量,它指定要使用的滞后。
Lst <- wave.multiple.cross.correlation(xx, lmax)
从我上面的示例中可以看出,我为滞后指定了两个不同的参数,lmax = 36 和 lag.max = 13
啊好吧,我只花了 100 个声望!
看来,鉴于您的时间序列的长度,您将最大小波级别 J 设置得太高(代码的第 14 行):
14 焦 <- 6
试试 J <- 4 和 运行
Lst <- wave.multiple.cross.correlation.prueba(xx, lag.max = lmax, ymaxr = NULL)
其中 lmax 是您之前设置为 36(或您选择的任何其他合理值)的最大延迟。
一般来说,最大小波水平的推荐值为J = t运行c(log2(T))-3(为了安全起见),其中T = 时间序列的长度.
事实上,报告的错误是您的原始 xx 列表的结果,其中级别 d6 和 s6 填充了 NaN。 请注意,这不是 wavemulcor 包的一部分,而是小波变换先前计算的一部分 在输入您选择的 wavemulcor 包函数之前。
希望对您有所帮助。