I have a quadratic regression model with three statistically significant coefficients: intercept, time, and time^2. The residuals of the model show significant autocorrelation. I understand that you can find an ARMA process for the residuals simply by using the ar or arima function and getting the coefficients [i.e. res.ar=ar(resid(fit),method='mle')]... but how do you refit the regression model with the autocorrelated residuals (for purposes in forecasting)? I'm thinking of terms of auto.arima
You are correct that it makes sense to estimate the "main part" of the model simultaneously allowing for the error term to have an ARMA structure (as opposed to doing this in two steps, which would be less efficient). Simply add the main regressors (
time^2) via the argument
xreg in the function
auto.arima. It will estimate a regression with ARMA errors. You can read more about this and related techniques in an enlightening blog post "The ARIMAX model muddle" by Rob J. Hyndman.
Rather than a two stage approach fitting time ,time sq , time cubed there could be level shifts , pulses and untreated ARIMA structure that are giving you "fits" . Following @whuber "Fitting polynomials to data can be a deceptively poor approach: a tiny bit of overfitting can result in models that are grossly bad because higher degree polynomials can (and often do) vary so wildly in between the data values and will be horrible extrapolators". . In building a model that includes memory , level shifts , time trends (as needed .. perhaps a break in trend ! , anomalies , changes in seasonal pulses , parameters and error variance you have all you need under 1 roof , so to speak .