时间戳

首页 > 解决方案 > 在时间上向前移动 n 个时间单位的窗口上计算的两个时间序列之间的相关系数

r - 在时间上向前移动 n 个时间单位的窗口上计算的两个时间序列之间的相关系数

问题描述

是否有一个包或一个简单的代码来生成 (1) 两个时间序列之间的相关系数图，这些时间序列在时间上向前移动 n 个时间单位 (2) 的窗口上计算，以及为每次移动计算的它们各自的 p 值？

library(zoo)

x = ts(rnorm(1:121), start = 1900, end = 2021)
y = ts(rnorm(1:121), start = 1900, end = 2021)
data = data.frame(x, y)

# 40-year moving window lagged forward by 15 years per example

rollapply(data, width=40, by = 15, 
          function(x) cor(x[,1],x[,2], method =  "pearson"),
          by.column=FALSE)

[1]  0.92514750  0.5545223 -0.207100231 -0.119647462 -0.125114237  0.041334073

**如果还计算 p 值会更好，Hmisc::rcorr但我没有设法将它集成到rollapply.

在这里的结果中，第一个系数 (0.9251...) 适用于 1900:1940，第二个系数适用于 1915:1955，依此类推。

所以问题是：有没有一种快速的方法可以将此结果整合到具有时间、r 和 p 值的阶梯图中？

输出如下所示：

时间 r 磷
1900 0.92 0.000001
1901 0.92 0.000001
... ... ...
1915年 0.55 0.00045
1916年 0.55 0.00045

标签： rcorrelationrollapply

解决方案

几点：

  • 从 1900 年到 2021 年包括 2021-1900+1 = 122 年，而不是 121
  • 40/15 参数不适用于 122 点，因此从 1907 开始

rcorr返回一个包含 3 个组件的列表，我们想要每个组件的 1,2 个元素。我们可以使用 na.locf 填充来自 rollapplyr 的缺失值。输入输出均为mts/ts系列。

library(zoo)
library(Hmisc)

set.seed(123)
tt <- ts(cbind(x = rnorm(115), y = rnorm(115)), start = 1907)

na.locf(rollapplyr(tt, width=40, by = 15, 
          function(x) sapply(rcorr(x), `[`, 1, 2),
          by.column = FALSE, fill = NA), fromLast = TRUE)

上面返回一个与输入 tt 具有相同行数但基于以下年份范围的计算 rcorr 的系列：

rollapplyr(1907:2021, 40, by = 15, range)
##      [,1] [,2]
## [1,] 1907 1946
## [2,] 1922 1961
## [3,] 1937 1976
## [4,] 1952 1991
## [5,] 1967 2006
## [6,] 1982 2021

推荐阅读