Initialize ARIMA simulations with different time-series

I have a fairly long time-series of annual abundances ($N_t$) of a wildlife species (73 years of abundances). To forecast the population’s trajectory, I have used ARIMA modeling. Examination of the ACF and PACF of the first-order differenced time-series suggested a 10-year cycle exists. So I used a span 10 seasonal difference to account for this periodic pattern. Therefore, the response variable was: $$Y_t=(\sqrt{N_t}-\sqrt{N_{t-1}})-(\sqrt{N_{t-10}}-\sqrt{N_{t-11}})$$ Typically, I would have used a logarithmic transformation but it resulted in heteroscedastic residuals. Examination of the ACF and PACF of $Y_t$ indicated a multiplicative seasonal structure so I fit the model: $$ARIMA(0,1,1)(0,1,1)_{10}$$ using the Forecast Package in R....library(forecast).

Example code for fitting the model:

m1=Arima(y,order=c(0,1,1),seasonal=list(order=c(0,1,1),period=10),include.mean=FALSE)

The residuals of this model were normally distributed, not autocorrelated, and homoscedastic.

I have been using the fitted model from above for some additional simulation work using the simulate.Arima function. However, I would like to initialize the simulation with a different time-series. The arima.sim function allows this but the arima.sim function doesn't seem to handle seasonal ARIMA models. With the simulate.Arima function one can use the future=TRUE option to simulate values that are "future to and conditional on the data" in the model m1. Can the data in the model object m1 simply be replaced to create a simulation that is conditional on different data?

For example:

# Create a new model object for simulation.
m.sim=m1
# Replace the data in the model object with the new data.
m.sim\$x=new
# Simulation conditional on the new data.
sim.forecasts=replicate(1000,simulate.Arima(m.sim,future=TRUE,bootstrap=TRUE))
• Just curious, why do you want to initialize the model with a different time series? – Michael Chernick Aug 29 '12 at 15:31
• @MichaelChernick I can think of few reasons why one might want to initialize the simulation with a different time-series. 1). One might have a population that is very similar to the one modeled but few time points to build a new model from (of course the big assumption would be the new population exhibits population dynamics just as the original one did which is probably not the best assumption but could still be informative about what to expect from the new population). – RioRaider Aug 29 '12 at 16:16
• 2). Also, one might want to simulate potential impact of a rare but extreme short-term event that causes high mortality in the population (i.e., hurricane, chemical spill, disease outbreak, etc). Of course one would have make some assumptions about how those events might change the population’s growth rate but could use the simulation results to inform wildlife managers as to potential effects by magnitude and duration of declines in population growth rates on the population’s trajectory. Those results could then be used to help inform planning and management decisions for the population. – RioRaider Aug 29 '12 at 16:17

You can "fit" the model to different data and then simulate:

m2 <- Arima(z,model=m1)
simulate.Arima(m2,future=TRUE,bootstrap=TRUE)

m2 will have the same parameters as m1 (they are not re-estimated), but the residuals, etc., are computed on the new data.

However, I am concerned with your model. Seasonal models are for when the seasonality is fixed and known. With animal population data, you almost certainly have aperiodic population cycling. This is a well-known phenomenon and can easily be handled with non-seasonal ARIMA models. Look at the literature on the Canadian lynx data for discussion.

By all means, use the square root, but then I would use a non-seasonal ARIMA model. Provided the AR order is greater than 1, it is possible to have cycles. See

You can do all this in one step:

m1 <- auto.arima(y, lambda=0.5)

Then proceed with your simulations as above.

• Wouldn't using a non-seasonal ARIMA model with higher AR orders (p > 1) results in the cyclic pattern being dampened as the forecasts increase a cycle or two beyond the last data point? – RioRaider Aug 29 '12 at 19:02
• Yes, but that's because the cycles are not strictly periodic. They are almost periodic, and as you go further out in the forecast horizon, it becomes harder to predict which part of the cycle you will be in. The point forecasts are means, and so they will naturally flatten out, indicating this uncertainty. – Rob Hyndman Aug 30 '12 at 4:24
• I still think this is not a very sensible thing to do. It is a rather strong assumption to think that a different series will have the same model form and coefficients as another series. When you don't have sufficient data to forecast I think it is best not to. – Michael Chernick Aug 30 '12 at 12:19
• @RobHyndman Thank you for your insight. Is the link you were refering to Research tips - Cyclic and seasonal time series? – RioRaider Aug 30 '12 at 14:05
• That's relevant, but I was referring to papers such as publish.csiro.au/?paper=ZO9530163 – Rob Hyndman Aug 30 '12 at 20:08

On a related note, you can accomplish the same objective if your ARIMA model has external regressors. This has been helpful for me on occasion.

For instance, say your first model was created as follows:

fit.arimax <- Arima(response, order=c(1, 0, 1), xreg=xreg)

Then suppose that after creating your model, you observe additional values in your response and external regression variables, and would like to forecast or simulate future outcomes given these new observations. E.g., say you are predicting electricity demand, and you observe another hour of demand (i.e. response) and temperature (i.e. external regression) data.

Then, you may fit the original model to the updated time series as follows, where response.new and xreg.new are your updated response and regression variables.

fit.arimax.new <- Arima(response.new, model=fit.arimax, xreg=xreg.new)

You can use this new model to forecast or simulate future outcomes, conditional on all observed data. Note that you must provide forecast external regressors for each. E.g.,

forecast.Arima(fit.arimax.new, h=length(xreg.forecast), xreg=xreg.forecast)

simulate.Arima(fit.arimax.new, n=length(xreg.forecast), xreg=xreg.forecast)

Another way to accomplish all of this is to make an entirely new model using the updated data. But the method described above is appropriate in real-time applications, in which case fitting a new ARIMA model would take too long.