arima regression Improving the model

Question

Given the following model which relates the full year home sales to the unemployment rate (observed or estimated) I get a projected increase of 14% for 2013 over 2012... last year the same approach over projected by 6% (the 2012 projection was for 41,992 and the actual is coming in about 39,535) So I think the model is over projecting but I'm at a dead end thinking of a good (valid) way to modify it? BTW the projection)s) is identical to what a simple straight line regression in a spreadsheet yields.

I want to also thank the fine people at Stack overflow who got me this far in my 1st foray into R http://stackoverflow.com/questions/14032768/csv-input-to-r-forecast-with-dates-via-r-studio#14032768

Pointers appreciated.

 # load the base data as presented in the question
 Workbook1 <- structure(list(Year = structure(1:10, .Label = c("31-Dec-04", 
"31-Dec-05", "31-Dec-06", "31-Dec-07", "31-Dec-08", "31-Dec-09", 
"31-Dec-10", "31-Dec-11", "31-Dec-12", "31-Dec-13"), class = "factor"), 
total = c(51439L, 59674L, 58664L, 55698L, 42235L, 37918L, 
36234L, 36965L, 39535L, NA), UnemplRt = c(5.7, 4.7, 3.8, 
3.7, 4.3, 8.5, 10.9, 10, 8.3, 7.1)), .Names = c("Year", "total", 
"UnemplRt"), class = "data.frame", row.names = c(NA, -10L))

# Make a time series out of the value
dependent <- ts(Workbook1[1:9,]$total, start=c(2004), frequency=1)

# load forecast package
require(forecast)
# load independent variables in variables.
unemployment <- ts(Workbook1[1:9,]$UnemplRt, start=c(2004), frequency=1)
    unemployment_future <- ts(Workbook1[10:10,]$UnemplRt, start=c(2004), frequency=1)

# make a model that fits the history
fit2 <- auto.arima(dependent, xreg=unemployment)

# generate a forecast with the already known unemployment rate for 2013.
fcast2 <- forecast(fit2,xreg=unemployment_future)
fcast2
     Point Forecast    Lo 80    Hi 80    Lo 95    Hi 95
2013       45168.02 38848.92 51487.12 35503.79 54832.25

This yield exactly the same result a a simple spreadsheet liniear regression. And I belive the suggested increase in home sales is too high, last year the projection turned out to be about 6% too high. So I'm trying to torture the numbers in some statically valid method to come up with a somewhat lower number for 2013

Answer 1

I don't think 9 observations are enough to get a very good model. But my first question is that why you are trying time series approach? Is there any prior knowledge about that? So let's look at the data first.

plot(1:9,dependent,type="l",lwd=2)

To me, this looks like a polynomial in time. So let's try a very simple model based on linear models:

> t=1:9
> fit3=lm(dependent~poly(t,d=4)+unemployment)

> summary(fit3)

Call:
lm(formula = dependent ~ poly(t, d = 4) + unemployment)

Residuals:
      1       2       3       4       5       6       7       8       9 
  467.5 -1060.9  -673.9  3285.1 -1908.3  -786.1   256.0   784.4  -363.6 

Coefficients:
                 Estimate Std. Error t value Pr(>|t|)  
(Intercept)      43961.04    9056.19   4.854   0.0167 *
poly(t, d = 4)1 -25246.64    8288.92  -3.046   0.0556 .
poly(t, d = 4)2    -65.03    3248.32  -0.020   0.9853  
poly(t, d = 4)3  15749.95    6348.52   2.481   0.0892 .
poly(t, d = 4)4  -3661.36    3586.58  -1.021   0.3825  
unemployment       379.18    1355.26   0.280   0.7978  
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 

Residual standard error: 2426 on 3 degrees of freedom
Multiple R-squared: 0.977,      Adjusted R-squared: 0.9386 
F-statistic: 25.46 on 5 and 3 DF,  p-value: 0.01162 

Then I will compare the fitted values of your model (fit2) with my model (fit3).

> plot(1:9,dependent,type="l",ylim=range(dependent,fitted(fit2),fitted(fit3)),lwd=2)
> points(1:9,fitted(fit2),type="l",col="blue")
> points(1:9,fitted(fit3),type="l",col="green",lwd=2)

The blue one is based on fit2 and the green one is based on fit3. You can see that the polynomial (green curve) is doing a good job compared to fit2. We can compare the two models simply by using AIC criteria:

> AIC(fit3)
[1] 169.9483
> AIC(fit2)
[1] 184.5996

AIC for the fit3 is less than fit2. So fit3 provides a better fit (at least based on this criteria). As I said, this is just a very simple model, you should check the residuals for any possible problem such as auto correlation before using the model for any predictions.

But if you want to get a prediction at this stage then here is how to do it:

> predict(fit3,newdata=data.frame(t=10,unemployment=unemployment_future))
       1 
46118.13