Р: Эффективное расположение сегментов временных рядов с максимальной взаимной корреляцией с входным сегментом?

Question

У меня есть длинные числовые данные временного ряда приблизительно 200 000 строк (давайте назовем это Z ).

В цикле я подвожу x (около 30) последовательных строк из Z за раз и рассматриваю их как точку запроса q .

Я хочу найти в Z y (~ 300) наиболее коррелированные сегменты временного ряда длины x (наиболее коррелированные с q ).

Какой эффективный способ сделать это?

Answer 1

Код ниже находит 300 сегментов, которые вы ищете, и запускается за 8 секунд на моем не слишком мощном ноутбуке с Windows, поэтому он должен быть достаточно быстрым для ваших целей.

Во-первых, он создает матрицу 30 на 199971 (Zmat), столбцы которой содержат все 30 отрезков временного ряда, которые вы хотите исследовать.Один вызов cor(), работающий с вектором q и матрицей Zmat, затем вычисляет все требуемые коэффициенты корреляции.Наконец, результирующий вектор проверяется для идентификации 300 последовательностей, имеющих самые высокие коэффициенты корреляции.

# Simulate data
nZ <- 200000
nq <- 30
Z <- rnorm(nZ)
q <- seq_len(nq)

# From Z, construct a 30 by 199971 matrix, in which each column is a
# "time series segment". Column 1 contains observations 1:30, column 2
# contains observations 2:31, and so on through the end of the series.
Zmat <- sapply(seq_len(nZ - nq + 1),  
               FUN = function(X) Z[seq(from = X, length.out = nq)])

# Calculate the correlation of q with every column/"time series segment.
Cors <- cor(q, Zmat)

# Extract the starting position of the 300 most highly correlated segments    
ids <- order(Cors, decreasing=TRUE)[1:300]

# Maybe try something like the following to confirm that you have
# selected the most highly correlated segments.
hist(Cors, breaks=100)
hist(Cors[ids], col="red", add=TRUE)

Answer 2

Наивное решение действительно очень медленное (хотя бы несколько минут - я недостаточно терпелив):

library(zoo)
n <- 2e5
k <- 30
z <- rnorm(n)
x <- rnorm(k) # We do not use the fact that x is a part of z
rollapply(z, k, function(u) cor(u,x), align="left")

Вы можете вычислить корреляцию вручную, начиная с первых моментов и комментариев, но это все равно занимает несколько минут.

y <- zoo(rnorm(n), 1:n)
x <- rnorm(k)
exy <- exx <- eyy <- ex <- ey <- zoo( rep(0,n), 1:n )
for(i in 1:k) {
  cat(i, "\n")
  exy <- exy + lag(y,i-1) * x[i]
  ey  <- ey  + lag(y,i-1) 
  eyy <- eyy + lag(y,i-1)^2 
  ex  <- ex  + x[i]    # Constant time series
  exx <- exx + x[i]^2  # Constant time series
}
exy <- exy/k
ex <- ex/k
ey <- ey/k
exx <- exx/k
eyy <- eyy/k
covxy <- exy - ex * ey
vx <- exx - ex^2
vy <- eyy - ey^2
corxy <- covxy / sqrt( vx * vy )

Как только у вас есть временной ряд корреляций, положение вершины 300 легко извлечь.

i <- order(corxy, decreasing=TRUE)[1:300]
corxy[i]