Message Boards Message Boards

Time series with multiple seasonality using Fourier external regressors?

Posted 8 years ago

I'm studying a time series with clear seasonality at a daily and weekly period. However, the series is aggregated in 15-min (.25-hour) timesteps, and thus the daily and weekly periods are 96 and 672, respectively. Due to the nature of this particular problem, further aggregation at a larger timescale is not desirable, so let's say we're stuck with these perhaps large periods.

Let me first mention that, ostensibly due to the large period 672, even obtaining the usual fit hangs, abort can't seem to intervene, and I end up quitting the kernel:

p=d=q=sp=sd=sq=Automatic;  
modtest=TimeSeriesModelFit[ data, {"SARIMA", {{p, d, q}, {sp, sd, sq}, #}}] & /@ {96,672}

although 96 alone works fine.(* SARMA {{4,1},{0,3},96} *). The AIC for this 96-model is significantly lower than that of mere TimeSeriesModelFit[data];, as well.

Therefore, both because the 672 period cannot be found with TimeSeriesModelFit and in order to model multiple seasonality in time series generally, I would like to implement the standard technique of including Fourier terms with period 96 and with 672 as so-called external regressors, as explained here and here.

To be honest, my initial approach was to use RLink[] to achieve this, but issues surrounding loading additional R packages on Mac stymied the operation. I'm now reaching out to see if anyone has idea on how to do this directly, within Mathematica. For example, is it possible to "force" TimeSeriesModelFit to fit to a known process plus some additional functions, e.g., ARMA[p,q]+f1+f2, where f1, f2 are Fourier terms analogous to external regressors in R.

Or, is it possible to somehow use ARMAProcess as a valid argument in LinearModelFit?

Any ideas are welcome. Thanks!

Attachments:
POSTED BY: Crescent Curtis
6 Replies

That's one way of doing it yes, you can then subtract this fitted model from you data and perform a TimeSeriesModelFit on it. This result is probably good enough already depending on your requirements. If you want to do better iterate, but using your estimated correlation between the data (obtained from TimeSeriesModelFit) to perform a NonLinearModelFit. This is needed because in your first run, you implictly assumed your data to be uncorrelated when doing LinearModelFit.

POSTED BY: Eduardo Serna

A possible workaround in Mathematica is to first perform (least squares) regression with the seasonal trend (i.e. detrend the data), and then perform TimeSeriesModelFit, the result can be iterated ad infinitum by using the fitted model to perform a generalized least squares regression, which if the problem is well posed (for example use the same regressors every iteration and non degenerate) should converge eventually, probably in only a few iterations.

POSTED BY: Eduardo Serna

Thanks for the reply; I'm relatively new to this and trying to parse your suggestion. When you say to perform "regression with the seasonal trend" do you mean by, e.g., LinearModelFit[...,{f1,f2},t] with f1, f2 Fourier terms? Thanks.

POSTED BY: Crescent Curtis

Is there a link to a representative data set?

POSTED BY: Daniel Lichtblau

I just added a sample dataset in .mx format. Please let me know if this doesn't work out. Thanks!

POSTED BY: Crescent Curtis

Thanks, but...does it open for you? I saved it to my /tmp directory and I get no data when I do

In[291]:= data = Get["/tmp/somedata.mx"];

In[288]:= Length[data]

Out[288]= 0

Opening it in an editor does indicate it was constructed using DumpSave.

[Later note] I can open it directly using Mathematica. Not sure why Get is giving me a hard time but what I have now is usable. Thanks again.

POSTED BY: Daniel Lichtblau
Reply to this discussion
Community posts can be styled and formatted using the Markdown syntax.
Reply Preview
Attachments
Remove
or Discard

Group Abstract Group Abstract