Текущая дисперсия, когда временное окно непостоянно

Я пытаюсь рассчитать скользящую дисперсию с окном, скажем, 4 года, для каждого из names A, B и C. Данные еженедельные:

> head(data1, 17)
         date name       value
1  1985-01-01    A -0.44008233
2  1985-01-01    B          NA #Observe that there are some NA's
3  1985-01-01    C  0.38682496
4  1985-01-08    A  0.41806540
5  1985-01-08    B -0.05460831
6  1985-01-08    C -0.52051435
7  1985-01-15    A  1.25769395
8  1985-01-15    B  0.80272053
9  1985-01-15    C -0.34501742
10 1985-01-22    A -0.43401839
11 1985-01-22    B  0.91113966
12 1985-01-22    C  1.07131717
13 1985-01-29    A -1.55395857
14 1985-01-29    B -0.43281709
15 1985-01-29    C  0.98034779
16 1985-02-05    A  1.70557396
17 1985-02-05    B  0.44688788

Мой подход до сих пор заключается в dcast данных, а затем запускать столбцы rollapply() (zoo) с движущимся окном 192 = 4 * 12 * 4:

v <- dcast(data1, date ~ name, value.var = "value")
var <- rollapply(v[-1], width=4*12*4, var, fill=NA, by.column = T)
var <- cbind(v$date, var)
var[,1] <- as.Date(var[,1])

Однако я понял, что для некоторых месяцев у меня есть четыре наблюдения (например, 7, 14, 21, 28 февраля), а для некоторых у меня пять еженедельных наблюдений (например, 1 , 8, 15, 22 и 29 января), поэтому использование окна 4 years * 12 months * 4 weeks наблюдений не корректно. Я думал добавить эти дополнительные наблюдения во временное окно (width), но я не уверен, как (и возможно ли это вообще), поскольку они меняются в зависимости от того, сколько 5 недель в месяц и сколько 4- недель в месяц наблюдения находятся внутри временного окна.

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

Пример данных:

set.seed(486)
date <- rep(seq(as.Date("1985-01-01"), as.Date("2010-01-1"), by="weeks"), each=3)
N <- length(date)
name <- c("A","B","C")
value <- rnorm(N)
i<-which(value %in% sample(value, 25)) ;i
j<-which(value %in% sample(value, 150)) ;j
value[i] <- NA
value[j] <- 0
data1 <- data.frame(date, name, value)

r variance sliding-window
person Per    schedule 18.06.2015    source источник

Ответы (2)


arrow_upward
2
arrow_downward

4 года имеют 208 недель плюс 5 дней, поэтому они не делятся на недели поровну. Если мы используем 209 недель, то мы отстаем всего на 2 дня за 4 года, поэтому давайте попробуем это.

Сначала преобразуйте класс data1 в класс "zoo", разделив данные на отдельные столбцы в соответствии со значением второго столбца. z будет иметь по одному столбцу для каждого из A, B и C. Затем определите функцию дисперсии, исключающую нули, и используйте ее с rollapplyr.

library(zoo)
z <- read.zoo(data1, split = 2) # 1305 x 3 
var0 <- function(x) var(x[x != 0])
r <- rollapplyr(z, 209, var0)

Оставить его как объект зоопарка может быть достаточно, но это приведет к преобразованию его в фрейм данных с 4 столбцами со столбцами Index, A, B и C:

fortify.zoo(r)
person G. Grothendieck    schedule 18.06.2015
comment
Хорошо, это хорошее приближение. На данный момент я пытаюсь подумать, можем ли мы использовать endpoints() для дат, чтобы получить точное решение. - person Per; 19.06.2015
comment
Точного решения не существует, потому что, как объяснялось, 4 года не кратны неделям. - person G. Grothendieck; 19.06.2015
comment
Да, я это понимаю, поэтому я пытаюсь подумать, можем ли мы изменить подход и вместо использования ряда наблюдений для движущегося окна использовать конкретные индексы из endpoints() (или, что то же самое, их расстояние). Это означало бы, что другая пара точек определяет скользящее окно для каждой вычисляемой дисперсии, что немного усложняет задачу. - person Per; 19.06.2015
comment
Это совсем не помогает. 4 года по-прежнему не кратны 7 дням, поэтому либо последняя неделя будет частично отсутствовать, а частично находиться в 4-летнем промежутке (что и делает использование 209 недель), либо будет одна неделя без дней в 4-недельном промежутке. что и делает 208 недель. Указываете ли вы недели по конечным точкам или нет, не имеет значения. - person G. Grothendieck; 19.06.2015

arrow_upward
2
arrow_downward

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

library(data.table)
library(zoo)
setDT(data1)[,var := {
           v1 <- rollapplyr(value,width=4*12*4, var, fill=N)
           v2 <- rollapplyr(value,width=4*12*5, var, fill=N)
           (v1+v2)/2},  name]

PS: Здесь я использую data.table, потому что он подходит для операций разделения (на группу) и повторной привязки.

Редактировать

Вы также можете преобразовать свои еженедельные данные в ежедневные, тогда вы сможете более точно рассчитать рулон на этой основе. Идея состоит в том, чтобы создать ежедневный индекс и объединить его с исходными данными. Это создаст новую таблицу данных с отсутствующими значениями. Вы заменяете отсутствующие значения первыми неотсутствующими значениями, используя na.locf.

library(data.table)
library(zoo)
ID <- 
data.table(
  date = seq(as.Date("1985-01-01"), as.Date("2010-01-1"), by="days"))
setkey(ID,date)

setDT(data1)[,date:=as.Date(date)][, 
        {
          merge(ID,.SD,all.x=TRUE)[,value := na.locf(value)]
        },
        
        name]
person agstudy    schedule 18.06.2015
comment
Это интересно. Я думаю, мы могли бы сделать еще один шаг вперед, используя средневзвешенное значение вместо простого среднего, после подсчета того, сколько существует еженедельных наблюдений 5 в месяц и 4 в месяц. Хотя меня больше интересует вычисление точного решения. Я начинаю подозревать, что прямое использование индексов наблюдений — не лучший способ. Мы должны как-то использовать даты. Все еще не уверен. - person Per; 19.06.2015
comment
@Per очень легко рассчитать вес. Я могу показать вам это, но если вам нужно точное решение, я думаю, вам следует аппроксимировать (интерполировать) ежедневный временной ряд из еженедельного. - person agstudy; 19.06.2015
comment
Да, мне кажется, я тоже знаю, как рассчитать веса, спасибо ^^ Что касается создания ежедневных рядов, я немного беспокоюсь, что тогда у меня будет слишком много данных, чтобы эффективно обрабатывать их в R, поскольку мой набор данных уже большой. Без учета выходных он вырастет в 5 раз, верно? - person Per; 19.06.2015
comment
@Per Я добавляю код для преобразования серии weekly.daily. - person agstudy; 19.06.2015