Взаимная корреляция различных временных рядов в матрице

В поисках решения моей проблемы я нашел старый пост (Взаимная корреляция различных значений данных временных рядов в R), который спрашивает именно то, что мне нужно, но, к сожалению, не получил никакого ответа, поэтому я спрошу еще раз, надеясь на какое-то руководство.

Я создал большую матрицу из большого количества временных рядов одинакового размера, каждый столбец представляет собой отдельный временной ряд (что-то похожее на следующее, но намного больше и гораздо больше значений, отличных от нуля):

      [,1]  [,2]  [,3]  [,4]  [,5]  [,6]  [,7]  [,8]  [,9] [,10] [,11] [,12] [,13] [,14] [,15] [,16] [,17] [,18] [,19]
[1,]    0   0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0    NA    NA    NA   0.0    NA   0.0   0.0   0.0   0.0
[2,]    0   6.0   0.0   9.0   0.0   0.0   0.0   0.0   0.0   0.0    NA     0    NA   0.0    NA   0.0   0.0   0.0   0.0
[3,]    0   0.0   0.0   5.0   0.0   0.0   0.0   0.0   0.0   0.0    NA     0    NA   0.0    NA   0.0   0.0   0.0   0.0
[4,]    0   0.0   0.0  10.0   0.0   0.0   0.0   0.0   0.0   0.0    NA     0    NA   0.0    NA   0.0   0.0   0.0   0.0
[5,]    0   0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0    NA     0    NA   0.0    NA   0.0   0.0   0.0   0.0
[6,]    0   0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0    NA     0    NA   0.0    NA   0.0   0.0   0.0   0.0
[7,]    0   0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0    NA     0    NA   0.0    NA   0.0   0.0   0.0   0.0
[8,]    0   0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0    NA     0    NA   0.0    NA   0.0   0.0   0.0   0.0
[9,]    0   0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0    NA     0    NA  10.0    NA   0.0   0.0   0.0   0.0
.
.
.

Я хочу определить корреляцию между всеми временными рядами, я поместил их в матрицу, потому что подумал, что это лучший способ выполнить процедуру кросс-корреляции, я могу ошибаться.

Итак, я также знаю о функциях "ccf" и "diss()":

  1. ccf() # в базовых пакетах
  2. diss(meter_daywise,METHOD = "CORT",deltamethod = "DTW")# в пакете TSclust

но, как и в старом посте, у меня те же проблемы:

  1. ccf не принимает полную матрицу в качестве входных данных
  2. diss() принимает входную матрицу и создает некоторую матрицу, но, наблюдая за значениями, я обнаруживаю, что это не матрица взаимной корреляции, поскольку значения не находятся между -1 и 1.

Итак, вопрос в том, как мы вычисляем и выполняем взаимную корреляцию между различными временными рядами в R?


r time-series cross-correlation
person Davido    schedule 18.05.2017    source источник
comment
Как вы ожидаете, что результат будет выглядеть?   -  person emilliman5    schedule 18.05.2017
comment
Я надеялся создать в результате матрицу размерности AxA, где A=количество временных рядов, и каждое значение в этой матрице с координатами [x,y] является корреляцией между временным рядом (x) и временным рядом (y ). Возможно ли это и есть ли смысл?   -  person Davido    schedule 18.05.2017
comment
Похоже, вам может понадобиться старая добрая корреляционная матрица. Для вашей структуры это будет cor(myMat).   -  person lmo    schedule 18.05.2017

Ответы (2)


arrow_upward
1
arrow_downward

ccf возвращает парную корреляцию при каждом смещении (т.е. отставании), но я думаю, что вам нужно получить максимальную (абс (корреляцию)). Поскольку у вас есть NA, вам нужно установить аргумент na.action.

mat <- matrix(rnorm(100000), ncol=100)
mat[sample(1:length(mat), 100)] <- NA 

res <- sapply(1:ncol(mat), function(x) {
  sapply(1:ncol(mat), function(z){
    resTmp <- ccf(x = mat[, x], y = mat[, z], plot=F, na.action = na.pass)
    resTmp$acf[which.max(abs(resTmp$acf))]
  })
})

Из справки ccf:

По умолчанию пропущенные значения не допускаются. Если функция na.action проходит через отсутствующие значения (как это делает na.pass), ковариации вычисляются из полных наблюдений. Это означает, что вычисленная оценка может не быть действительной последовательностью автокорреляции и может содержать пропущенные значения.

person emilliman5    schedule 18.05.2017
comment
спасибо, но я не могу заставить его работать, должен ли я использовать код, кроме первых двух строк, и использовать мою матрицу (matri123) вместо матрицы мата этого примера? Я пробовал так и не работает, не могли бы вы помочь мне? Дополнительный вопрос: почему лучше применять абсолютную максимальную корреляцию? - person Davido; 19.05.2017
comment
1) Да, используйте свою матрицу вместо mat. Вам нужно будет уточнить, почему / как это не работает. Сообщения об ошибках неверные результаты? 2) Не обязательно лучше использовать корреляцию абс-макс, это было обоснованное предположение с моей стороны. ccf возвращает вектор корреляции для пары временных рядов, и поскольку вам нужна матрица AxA, мы должны уменьшить каждую парную корреляцию до одного значения, в данном случае это будет максимальная корреляция (без учета знака). Если вы уточните свой вариант использования, мы можем придумать что-то еще. - person emilliman5; 19.05.2017
comment
В результате я получаю матрицу AxA с тем, что кажется корреляцией каждой пары, что в значительной степени то, что я искал. У меня есть 2 последних вопроса, и мне жаль, что я мало знаю об этой теме. 1. Многие элементы в результирующей матрице: числовые, 0 знаете ли вы, почему это происходит и как этого избежать? 2. Есть ли способ получить значение корреляции не для пары временных рядов, а для временного ряда с их группой (например, корреляция временного ряда A с группой временных рядов A, B, C, D, E )? Я очень ценю любое руководство# - person Davido; 19.05.2017
comment
1) вам нужно будет создать воспроизводимый пример для устранения неполадок, но я подозреваю, что это связано с большим количеством столбцов NA или со всеми 0. 2) Да, вы можете выполнить иерархическую кластеризацию, а затем использовать дендрограмму для приблизительного расстояния между одним временным рядом и группой временных рядов. Но вы не можете сделать это, как в парном сценарии. - person emilliman5; 19.05.2017

arrow_upward
1
arrow_downward

Одна из возможностей — запустить ccf для всех комбинаций ваших столбцов, используя combn. Следующий код был протестирован по вопросу в ссылке:

myResults <- combn(seq_len(nrow(meter_daywise)), 2,
                   FUN=function(x) ccf(meter_daywise[x[1],], meter_daywise[x[2],]),
                   simplify=FALSE)

и создает вложенный список, подобный этому

str(myResults)
List of 10
 $ :List of 6
  ..$ acf   : num [1:17, 1, 1] 0.0241 0.0895 0.1463 0.0583 -0.0613 ...
  ..$ type  : chr "correlation"
  ..$ n.used: int 15
  ..$ lag   : num [1:17, 1, 1] -8 -7 -6 -5 -4 -3 -2 -1 0 1 ...
  ..$ series: chr "X"
  ..$ snames: chr "meter_daywise[x[1], ] & meter_daywise[x[2], ]"
  ..- attr(*, "class")= chr "acf"
 $ :List of 6
  ..$ acf   : num [1:17, 1, 1] -0.445 -0.493 -0.239 0.465 0.49 ...
  ..$ type  : chr "correlation"
  ..$ n.used: int 15
  ..$ lag   : num [1:17, 1, 1] -8 -7 -6 -5 -4 -3 -2 -1 0 1 ...
  ..$ series: chr "X"
  ..$ snames: chr "meter_daywise[x[1], ] & meter_daywise[x[2], ]"
  ..- attr(*, "class")= chr "acf"

...

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

myResults <- combn(seq_len(ncol(myMat)), 2,
                   FUN=function(x) ccf(myMat[, x[1]], myMat[, x[2]]), simplify=FALSE)

где myMat — имя вашей матрицы. Вы можете увидеть пары с более простым вызовом combn, например

myPairs <- combn(seq_len(ncol(myMat)), 2)
person lmo    schedule 18.05.2017
comment
спасибо за ответ, к сожалению, я не могу заставить его работать. Моя матрица (маты) имеет размеры 25 столбцов x 13514 строк, я попробовал с вашей командой: myResults ‹- combn(seq_len(ncol(matts)), 2, FUN=function(x) ccf(matts[, x[1]] , matts[, x[2]]), simple=FALSE), и я получаю сообщение об ошибке: Ошибка в na.fail.default(as.ts(x)) : пропущенные значения в объекте Поэтому я подумал, что это может быть связано с Значения NA, поэтому замените их на другие значения только для проверки, работает ли это, а затем я получаю сообщение об ошибке: Ошибка в plot.window (...): необходимы окончательные значения «ylim». Я действительно потерялся здесь. что я могу сделать? - person Davido; 19.05.2017