对非均匀间隔时间序列进行季节性分解，R 或 Python 中是否有任何完善的算法？

Question

R 中的 stats 包具有 stl()，但它需要由 ts() 创建的均匀间隔的时间序列。它不能处理动物园对象。

奇怪的是，它也不能处理缺失值，尽管 STL 方法声称可以用 LOESS 填充缺失值。 （请参阅简历上的this question。）

因此，例如，如果您有工作日数据，则不能通过将 NA 放在周末并调用 stl() 来将其设为日历日。

我还看到 Python statsmodel 的作者试图迁移 stl() 以使用 Pandas TimeSeries，但它似乎还没有。

谢谢

编辑：补充一点，我知道我可以做一个非常简单的模型，比如拟合谐波，但我想要一个完善的模型，至少可以提供基准。我有分月数据，所以 X12 不适用。

Answer 1

根据@Julius 在这个post 中，是否可以使用stl 和na.approx，来自zoo 包，使用stl(x, na.action = na.approx, ...)。这有点像interpolation。

不幸的是，stl 更喜欢常规时间序列。

所以我采用@Julius 方法来测试loess 在填充NAs 方面的表现，假设时间序列不规则。显然误差很小。

让我们模拟

以下函数采用3^(0:p) NAs 和naapprox=TRUE 获取数据并模拟场景，用于使用na.approx 和optspan=TRUE 优化span 参数。损失函数为 MAPE。

#The approach is based in @Julius stl with NAs
library(zoo)
library(plyr)
library(reshape)
library(ggplot2)
mape <- function(f, x) colMeans(abs(1 - f / x) * 100,na.rm=T)

loessCheck <- function(data,p=2, naapprox=TRUE, optspan=TRUE){
  set.seed(20130201)
  pos <- lapply(3^(0:p), function(x) sample(1:length(data), x))
  datasetsNA <- lapply(pos, function(x) {data[x] <- NA; data})
  original <- data.frame(y.predict=as.numeric(data))
  original$id <- "Original"
  datasetsNA <- lapply(datasetsNA, function(y){ 
    posna=which(is.na(y)) 
    tvo=1:length(y) 
    yo=y;tv=tvo
    if (any(posna%in%c(1,length(y)))) tv=tvo[-posna[which(posna%in%c(1,length(y)))]]
    if(naapprox) y=na.approx(y) #instead of inside loess
    if(optspan){ fseq=function(x, y)
      mape(matrix(predict(loess(y ~ dt, span=x, data.frame(dt=tv, y=y))),ncol=1),as.vector(y[!is.na(y)]))
    ospan <- optimize(fseq, c(0.1,1), maximum=FALSE,y=y)
    spanmim <- ospan$minimum
    } else spanmim <- 0.75
    y.loess <- loess(y ~ dt, span=spanmim, data.frame(dt=tv, y=y))
    y.predict <- predict(y.loess, data.frame(dt=tvo))
    y.predict[-posna] <- yo[-posna]
    data.frame(y.predict, 
               id = paste(length(posna), "NAs"), 
               stringsAsFactors = FALSE)
    })
  loessAll <- rbind.fill(c(list(original), datasetsNA))
  loessAll$Date <- time(data)
  results <- data.frame(y.predict = sapply(lapply(datasetsNA, '[', i = "y.predict"), mape, original[, "y.predict"]))
  results$id <- unique(loessAll$id[-(1:nrow(original))])
  results <- melt(results, id.var = "id")
  results$x <- min(loessAll$Date) + diff(range(loessAll$Date)) / 4
  results$y <- min(original[, "y.predict"],na.rm=T) + diff(range(original[, "y.predict"],na.rm=T)) / (4 * p) * (0:p)
  results$value <- round(results$value, 2)
  ggplot(loessAll, aes(x = Date, y = y.predict, colour = id, group = id)) + geom_line() + 
     theme_bw() +
    theme(legend.title = element_blank(), strip.background = element_rect(fill = "white")) + 
    labs(x = NULL, y = NULL) + scale_colour_brewer(palette = "Set1") +
    lapply(unique(results$id), function(z)
      geom_text(data = results, colour = "black", size = 3,
                aes(x = x, y = y, label = paste0("MAPE (", id, "): ", value, "%"))))
}
args(loessCheck) # first boolean for using na.approx, the second for optimize span parameter in sample


loessCheck(nottem,p=3) #T T

loessCheck(nottem,p=3,FALSE) #FT

loessCheck(nottem,p=3,FALSE,FALSE)

loessCheck(nottem,p=3,TRUE,FALSE)

使用 zoo 包中的 na.approx 时减少 MAPE，推荐使用 span= 0.75。填写完 NA 后，可以考虑使用其他建模替代方案。