4

I'm trying to model a time series with BSTS having a Poisson response variable. However, I just get an error message. Here is a reproducible example:

library(bsts)

holiday.list <- list(NamedHoliday("MemorialDay"),
                     NamedHoliday("IndependenceDay"),
                     NamedHoliday("LaborDay"),
                     NamedHoliday("Thanksgiving"),
                     NamedHoliday("Christmas"),
                     NamedHoliday("NewYearsDay"))

data <- ts(as.integer(EuStockMarkets))
ss <- AddLocalLinearTrend(list(), data)

ss <- AddRegressionHoliday(ss, data, holiday.list=holiday.list, time0=as.Date("1991-01-01"))

ss <- AddSeasonal(ss, data, nseasons=7) # weekly seasonal

bsts.poisson <- bsts(data, state.specification=ss, niter=500, family="poisson")

This exits with code 134 and prints the following message:

Abort trap: 6

It works without error when I remove the family="poisson" argument, but I need a Poisson response, not Gaussian. What am I doing wrong?

Edit: I know stock market data shouldn't really be a Poisson process. It is not what I'm modeling in my actual work. It's just a convenient substitute to provide a reproducible example.

Edit: Some version info - R version 3.6.0 (2019-04-26), Platform: x86_64-apple-darwin13.4.0 (64-bit), Running under: macOS Mojave 10.14.6

other attached packages:
[1] bsts_0.9.1          xts_0.11-2          zoo_1.8-6
[4] BoomSpikeSlab_1.1.1 Boom_0.9.1          MASS_7.3-51.4
Andy Carlson
  • 3,633
  • 24
  • 43
  • 1
    You seem to be modelling stock price time series data. What makes you think that is a Poisson process? – Maurits Evers Sep 12 '19 at 22:04
  • 1
    It's not what I'm actually modeling. It's just a reproducible example. I can't distribute my real data and `EuStockMarkets` was conveniently included in R in the global scope. – Andy Carlson Sep 13 '19 at 14:04
  • 3
    Your example is not actually reproducible. I get the error: "Error in bsts(data, state.specification = ss, niter = 500, family = "poisson"): Caught exception with the following error message: Cannot extract residual variance parameter." – ClancyStats Sep 19 '19 at 15:42
  • 1
    @AndyCarlson Did removing the `time0` fixed the issue? Did you encounter other problems after the removal? – David Sep 23 '19 at 16:11

2 Answers2

1

The problem is not with family="poisson". The problem is caused due to the input to AddRegressionHoliday the variable time0 causes the problem.

Running the following without time0:

library(bsts)

holiday.list <- list(NamedHoliday(holiday.name = "MemorialDay"),
                     NamedHoliday(holiday.name = "IndependenceDay"),
                     NamedHoliday(holiday.name = "LaborDay"),
                     NamedHoliday(holiday.name = "Thanksgiving"),
                     NamedHoliday(holiday.name = "Christmas"),
                     NamedHoliday(holiday.name = "NewYearsDay"))

data <- ts(as.integer(EuStockMarkets))
ss <- AddLocalLinearTrend(list(), data)
#ss <- AddRegressionHoliday(ss, data, holiday.list=holiday.list, time0=as.Date("1991-01-01"))
ss <- AddRegressionHoliday(ss, data, holiday.list=holiday.list)
ss <- AddSeasonal(ss, data, nseasons=7) # weekly seasonal

bsts.poisson <- bsts(data, ss, niter=500, family = "poisson")

resulted in the following output:

=-=-=-=-= Iteration 0 Sun Sep 22 22:35:51 2019 =-=-=-=-=
=-=-=-=-= Iteration 50 Sun Sep 22 22:35:57 2019 =-=-=-=-=
=-=-=-=-= Iteration 100 Sun Sep 22 22:36:03 2019 =-=-=-=-=
=-=-=-=-= Iteration 150 Sun Sep 22 22:36:10 2019 =-=-=-=-=
=-=-=-=-= Iteration 200 Sun Sep 22 22:36:16 2019 =-=-=-=-=
=-=-=-=-= Iteration 250 Sun Sep 22 22:36:23 2019 =-=-=-=-=
=-=-=-=-= Iteration 300 Sun Sep 22 22:36:29 2019 =-=-=-=-=
=-=-=-=-= Iteration 350 Sun Sep 22 22:36:36 2019 =-=-=-=-=
=-=-=-=-= Iteration 400 Sun Sep 22 22:36:42 2019 =-=-=-=-=
=-=-=-=-= Iteration 450 Sun Sep 22 22:36:49 2019 =-=-=-=-=

But in this case it will take data[1] as the time0, this might have an effect and might not depend on the data.

David
  • 8,113
  • 2
  • 17
  • 36
  • The problem is that by omitting the `time0` argument, the model is ignoring the holidays (I assume because it doesn't know the real dates for any of the time points). When I inspect the trained model and plot the components, it isn't using the holidays at all. – Andy Carlson Sep 24 '19 at 21:04
  • @AndyCarlson I understand, I've tried manipulating it by using `zoo` or `xts` objects with no success; might worth to raise an issue in the developer *GitHub* – David Sep 25 '19 at 08:07
1

I got an email response from Steven Scott, the author of bsts:

"Yes, some of the holiday models assume Gaussian errors because there is a shared residual variance parameter. That could be reworked, but it is low on my priority list."

So, it looks like using holidays with a Poisson response is currently unsupported.

Andy Carlson
  • 3,633
  • 24
  • 43