Time series analysis to quantify trend when seasonal amplitude is decreasing

Question

Short version: How would one be able to quantify an intervention effect in time-series analysis when the intervention decreases seasonal amplitude variation but doesn't directly effect the median?

Here is a link to my raw data.

I have a complex time-series of daily incidence numbers for a population over 7 years, totaling 2557 observations. There is a strong weekly and yearly seasonality (high incidence in winter months and low incidence in summer months). There is a baseline negative trend which is orders of magnitude smaller than the seasonality. An intervention was introduced at time = 1700. This intervention should theoretically not cause a level shift. My aim is to detect whether the intervention increases the baseline negative trend.

I have attempted to fit a dynamic linear regression with ARIMA errors in R using auto.arima() in the forecast package. I modeled the weekly season using a dummy variable for each weekday and the weekend. I modeled the monthly seasonality with harmonics using fourier() function in the forecast package. An the intervention effect was coded in by specifying the time index and post-intervention times as independent variables using the methods described in Segmented regression analysis of interrupted time series studies in medication use research. With these variables specified auto.arima() suggests an ARMA(7,7) process. The coefficients for baseline trend and post-intervention trend are however non-significant.

I am concerned that by using fourier terms to model away the seasonality I am artificially removing any intervention effect, as visual analysis of the time series indicates that the intervention is specifically decreasing incidence during the winter months and therefore reducing the yearly seasonal variability.

Answer 1

I took your data and used AUTOBOX. In addition your specified level shift at 1700 the program identified daily and monthly indicators and an additional level shift while incorporating an ARIMA model (1,0,0)(1,0,0) and a number of one time pulses. The equation is presented here and and . The residual plot is here and appears correct . The ACF of the original series is here while the ACF of the residual series is here . Don't be concerned about the apparent significant structure as the sample size is quite large yielding spurious limits. The Actual and Forecast picture is here with the forecasts shown here .

As to why your dynamic regression model didn't appear to work out my guesses (slanted opinions) are as follows in terms of potential importance:

1. The deleterious impact of all the untreated outliers (pulses)
2. The (way) over specification of your ARIMA model based upon a serious mis-identification approach (AIC/BIC which very simply assumes structure that doesn't exist such as the non-importance of the level shift series and of course no pulses)
3. Using harmonics rather than simple monthly effects
4. Omitting a second level shift (down) at observation 277

With respect to your reflection: " This intervention should theoretically not cause a level shift. My aim is to detect whether the intervention increases the baseline negative trend." . All level/step shifts reflect a change in intercept. A series that has 1,1,1,1,1,2,2,2,2,2,2 has a change in intercept. A series 1,2,3,4,5,15,16,17,18,19,..... has a change in intercept . No trend change was found ( perhaps because of the inclusion of your level/step indicator) but I can imagine a slight tweaking of conditions ( i.e. not using your level/step indicator and/or locking-out/precluding empirically developed level/step shifts ) it is possible that a useful model might include one or more trends.