Подтвердить что ты не робот

Оптимизированные функции качения на нерегулярных временных рядах с временным окном

Есть ли способ использовать rollapply (из zoo package или что-то подобное) оптимизированные функции (rollmean, rollmedian и т.д.) для вычисления функций качения с использованием временного окна вместо одного, основанного на числе наблюдений? Я хочу просто: для каждого элемента в нерегулярном временном ряду я хочу вычислить функцию качения с окном N дней. То есть, окно должно включать все наблюдения за N дней до текущего наблюдения. Временные ряды также могут содержать дубликаты.

Здесь следует пример. Учитывая следующие временные ряды:

      date  value
 1/11/2011      5
 1/11/2011      4
 1/11/2011      2
 8/11/2011      1
13/11/2011      0
14/11/2011      0
15/11/2011      0
18/11/2011      1
21/11/2011      4
 5/12/2011      3

Скользящая медиана с 5-дневным окном, выровненным вправо, должна приводить к следующему вычислению:

> c(
    median(c(5)),
    median(c(5,4)),
    median(c(5,4,2)),
    median(c(1)),
    median(c(1,0)), 
    median(c(0,0)),
    median(c(0,0,0)),
    median(c(0,0,0,1)),
    median(c(1,4)),
    median(c(3))
   )

 [1] 5.0 4.5 4.0 1.0 0.5 0.0 0.0 0.0 2.5 3.0

Я уже нашел некоторые решения там, но они обычно сложны, что обычно означает медленное. Мне удалось реализовать собственный калькулятор функции. Проблема в том, что для очень длинных рядов оптимизированная версия медианного (rollmedian) может сделать огромную разницу во времени, поскольку она учитывает перекрытие между окнами. Я бы хотел избежать повторного его реализации. Я подозреваю, что есть некоторые трюки с параметрами rollapply, которые заставят его работать, но я не могу понять это. Заранее спасибо за помощь.

4b9b3361

Ответ 1

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

Аргумент ширины rollapply (пакет zoo) может быть списком (подробнее см. справку rollapply). На основании этого я написал функцию, которая создает список, который будет использоваться с rollapply в качестве параметра ширины. Функция извлекает индексы для нерегулярных объектов зоопарка, если движущееся окно должно быть временным, а не индексированным. Поэтому индекс объекта зоопарка должен быть фактическим.

# Create a zoo object where index represents time (e.g. in seconds) 

d <- zoo(c(1,1,1,1,1,2,2,2,2,2,16,25,27,27,27,27,27,31),     
         c(1:5,11:15,16,25:30,31))

# Create function 

createRollapplyWidth = function(zoodata, steps, window ){   

  mintime =  min(time(zoodata))     

  maxtime =  max(time(zoodata)) 

  spotstime = seq(from = mintime , to = maxtime, by = steps)

  spotsindex = list() 

    for (i in 1:length(spotstime)){
    spotsindex[[i]] =  as.numeric(which(spotstime[i]  <=  time(zoodata) & time(zoodata) < spotstime[i] + window))}

  rollapplywidth = list()
    for (i in 1:length(spotsindex)){
    if (!is.na(median(spotsindex[[i]])) ){ 
      rollapplywidth[[round(median(spotsindex[[i]]))]] = spotsindex[[i]] - round(median(spotsindex[[i]]))}
  }
  return(rollapplywidth)
  }


# Create width parameter for rollapply using function

rollwidth =  createRollapplyWidth(zoodata = d, steps = 5, window = 5) 

# Use parameter in rollapply 

result = rollapply(d, width = rollwidth , FUN =  sum, na.rm = T) 
result

Ограничение: не основано на дате, но в секундах. Параметр "частичный" rollapply не работает.

Ответ 2

Вот моя работа с проблемой. Если такой подход зависит от того, что вы хотели (я не знаю, удовлетворительно ли это с точки зрения скорости), я могу написать его как более подробный ответ (хотя он основан на идее @rbatt).

library(zoo)
library(dplyr)

# create a long time series
start <- as.Date("1800-01-01")
end <- as.Date(Sys.Date())

df <- data.frame(V1 = seq.Date(start, end, by = "day"))
df$V2 <- sample(1:10, nrow(df), replace = T)

# make it an irregular time series by sampling 10000 rows
# including allowing for duplicates (replace = T)
df2 <- df %>% 
  sample_n(10000, replace = T)

# create 'complete' time series & join the data & compute the rolling median
df_rollmed <- data.frame(V1 = seq.Date(min(df$V1), max(df$V1), by = "day")) %>% 
  left_join(., df2) %>% 
  mutate(rollmed = rollapply(V2, 5, median, na.rm = T, align = "right", partial = T)) %>% 
  filter(!is.na(V2)) # throw out the NAs from the complete dataset

Ответ 3

Не проверяйте скорость, но если дата не имеет более чем max.dup, то должно быть, что последние 5 * max.dup записей содержат последние 5 дней, поэтому приведенная ниже однострочная функция fn на rollapplyr сделает это:

k <- 5

dates <- as.numeric(DF$date)
values <- DF$value

max.dup <- max(table(dates))

fn <- function(ix, d = dates[ix], v = values[ix], n = length(ix)) median(v[d >= d[n]-k])

rollapplyr(1:nrow(DF), max.dup * k, fn, partial = TRUE)
## [1] 5.0 4.5 4.0 1.0 0.5 0.0 0.0 0.0 2.5 3.0

Примечание: Мы использовали это для DF:

 Lines <- "
      date  value
 1/11/2011      5
 1/11/2011      4
 1/11/2011      2
 8/11/2011      1
13/11/2011      0
14/11/2011      0
15/11/2011      0
18/11/2011      1
21/11/2011      4
 5/12/2011      3
"
DF <- read.table(text = Lines, header = TRUE)
DF$date <- as.Date(DF$date, format = "%d/%m/%Y")

Ответ 4

Мы можем сделать это, просто используя базу, следующим образом:

Сначала настройте данные (на основе примечания @g-grothendieck)

library(data.table)
Lines <- "
      date  value
1/11/2011      5
1/11/2011      4
1/11/2011      2
8/11/2011      1
13/11/2011      0
14/11/2011      0
15/11/2011      0
18/11/2011      1
21/11/2011      4
5/12/2011      3
"
DT <- as.data.table(read.table(text = Lines, header = TRUE))
DT$date <- as.Date(DF$date, format = "%d/%m/%Y")
DT$row <- 1:NROW(DF)
setkey(DT, row, date) #mark columns as sorted, for speed

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

roll.median.DT <- function(x){
  this.date <- as.Date(x[1])
  this.row <- as.numeric(x[3])
  median(DT[row <= this.row & date >= (this.date-5)]$value) #NB DT is not defined within function, so it is found from parent scope
}
apply(DT, FUN=roll.median.DT, MARGIN = 1)
[1] 5.0 4.5 4.0 1.0 0.5 0.0 0.0 0.0 2.5 3.0