基于BSTS的贝叶斯结构时间序列建模：R BSTS包复刻Facebook Prophet方法问询

Great question! Yes, researchers and R practitioners have absolutely explored replicating Facebook Prophet's core methodology using the bsts package. Since both frameworks focus on decomposable time series models, the translation is quite intuitive—let’s walk through the key components and parameters involved.

Prophet’s model breaks down time series into trend, seasonality, holiday effects, and noise. Here’s how to map each component to bsts:

• Trend Component: Prophet’s flexible piecewise linear (or logistic) trend is replicated using bsts's AddLocalLinearTrend or AddSemilocalLinearTrend. For logistic growth (when you have a known carrying capacity), you can combine a regression term for the capacity variable with a trend component to mimic Prophet’s logistic trend behavior.
• Seasonality Components: Prophet’s built-in yearly, weekly, and daily seasonality maps directly to bsts's AddSeasonal function. You can add multiple seasonal layers (one for each frequency) to capture periodic patterns at different scales.
• Holiday/Special Event Effects: Similar to Prophet’s holidays parameter, you can encode holidays as binary indicator variables and add them to the model using AddDynamicRegression (or AddRegression for fixed effects). This lets the model estimate the impact of one-off or recurring events.
• Noise/Error Term: bsts uses a Gaussian observation error by default, which aligns with Prophet’s noise model. For time series with autocorrelated errors, you can optionally add an AddAutoArma component to mimic Prophet’s implicit handling of residual autocorrelation.

To match Prophet’s default behavior (or adjust it similarly), focus on these critical parameters:

• Trend Flexibility:
  • prior.level.sd (in AddLocalLinearTrend): Controls how much the level of the trend can change over time—analogous to Prophet’s changepoint_prior_scale. A higher value makes the trend more responsive to changepoints.
  • prior.slope.sd: Governs the variability of the trend’s slope, adding another layer of control over trend smoothness.
• Seasonality Strength:
  • prior.sd (in AddSeasonal): Determines the magnitude of seasonal effects, matching Prophet’s seasonality_prior_scale. Smaller values produce smoother seasonal patterns, while larger values allow more extreme seasonal fluctuations.
  • frequency: Sets the periodicity of the seasonality (e.g., 7 for weekly, 365 for yearly), just like specifying yearly.seasonality or weekly.seasonality in Prophet.
• Holiday Impact:
  • prior.sd (in AddDynamicRegression): Controls the strength of holiday effects, equivalent to Prophet’s holidays_prior_scale. Lower values shrink holiday impacts toward zero, reducing overfitting to rare events.
• Model Fitting:
  • niter: The number of MCMC iterations for bsts (since it’s a Bayesian sampling framework). Unlike Prophet’s optimization-based fitting, this parameter dictates how many posterior samples you generate for inference and prediction.

Here’s a minimal snippet to build a Prophet-like model in bsts:

library(bsts)
library(dplyr)

# Simulate time series data with trend, seasonality, and a holiday
set.seed(123)
dates <- seq.Date(as.Date("2020-01-01"), as.Date("2022-12-31"), by = "day")
y <- 50 + 0.01*(1:length(dates)) + sin(seq(0, 2*pi*3, length.out = length(dates)))*10
# Add a one-time holiday spike
y[dates == as.Date("2021-07-04")] <- y[dates == as.Date("2021-07-04")] + 15

# Prepare data frame with holiday indicator
data <- tibble(ds = dates, y = y) %>%
  mutate(independence_day = ifelse(ds == as.Date("2021-07-04"), 1, 0))

# Build model specification
model_spec <- AddLocalLinearTrend(list(), data$y)
model_spec <- AddSeasonal(model_spec, data$y, frequency = 7, prior.sd = 0.1)  # Weekly seasonality
model_spec <- AddSeasonal(model_spec, data$y, frequency = 365, prior.sd = 0.2) # Yearly seasonality
model_spec <- AddDynamicRegression(model_spec, formula = ~ independence_day, data = data, prior.sd = 0.5)

# Fit the model
bsts_model <- bsts(model_spec, data = data$y, niter = 1000)

# Generate predictions (with uncertainty intervals)
predictions <- predict(bsts_model, horizon = 30)

It’s worth noting that while this replicates Prophet’s core structure, bsts provides full posterior uncertainty estimates (a key Bayesian advantage) whereas Prophet uses approximate intervals. The two frameworks also differ in fitting methods—bsts uses MCMC sampling, while Prophet uses Stan-based maximum a posteriori (MAP) optimization—but the end result is a functionally similar decomposable time series model.

内容的提问来源于stack exchange，提问作者John Zicker