# Reproducing Pankratz's ARIMA models for U.S. savings rate

I'm wondering if anyone has been able to reproduce the U.S. savings rate models in Pankratz's "Forecasting with Dynamic Regression Models"?

pankratz table-1.6

The first link for me is his book, page 264, referring to Table 1.6. This doesn't exist in my library copy, and I haven't found anything in a google for an errata. So I skipped that, but in the next page, he presents an ARIMA(1,0,2) model for the data. With some help from the comments to this question, I found online data from which I extracted the 1959-1979 subset, but my coefficients differ significantly from his Pankratz's. From his textbook:

                                 t
Parameter       Value       Statistic
-----------   -----------   -----------
Constant        1.6138        3.35
AR{1}        0.7359        9.58
MA{2}        0.3437        3.31
Variance        0.663716^2=0.44


I am using arima class in Matlab's Econometrics toolbox. With the coefficient for the 1st MA lag term set to 0, arima(1,0,2) estimation yields:

                                  t
Parameter       Value        Statistic
-----------   -----------    -----------
Constant         2.4916          3.774
AR{1}       0.778541        13.8733
MA{2}       0.268419        3.45492
Variance        0.40444        8.23616


The AR1 term is close to Pankratz's, but the MA2 term is about 80% of Pankratz's. And the variance is is a bit off. The t values are roughly the same.

I noticed that even though the shape of the plot is the same, the scaling is not. Pankratz's peaks at 10% while the online data peaks at 15%. If I scale the latter to peak at 10%, the lag terms stay the same, but the variance & constant term changes:

                              Standard          t
Parameter       Value          Error       Statistic
-----------   -----------   ------------   -----------
Constant        1.66106      0.440134          3.774
AR{1}       0.778541     0.0561178        13.8733
MA{2}       0.268419     0.0776917        3.45492
Variance       0.179751     0.0218246        8.23616


The constant term is now close to Pankratz's, but the variances matches even less than before. Here is the code (I've zeroed out the subtraction of the mean because it only affects the constant term):

svRate=[ 9.3 9.5 10.0 9.9 10.6 11.1 11.3 11.5 11.1 11.6 11.4 10.8 11.3 11.2 11.6 11.5 10.7 10.7 9.6 10.1 10.3 9.7 10.2 10.0 10.9 10.9 11.7 11.6 11.5 11.3 11.0 10.6 10.7 10.6 10.3 11.0 11.0 11.8 11.3 12.0 11.1 11.1 11.9 11.3 10.9 10.9 10.9 11.6 12.3 11.8 12.2 12.4 11.8 11.9 10.5 10.6 9.8 10.1 11.6 11.5 11.8 12.5 13.1 13.1 13.2 13.6 13.4 12.9 12.2 11.4 11.7 13.2 12.2 13.0 13.0 14.1 13.6 12.5 12.2 13.3 12.4 15.0 12.5 12.3 11.7 11.3 11.1 10.4 9.4 10.0 10.5 10.8 10.9 9.9 10.1 10.1 10.6 9.8 9.4 9.5 ]';

fprintf('max(svRate) before scaling to 10: %g\n',max(svRate));
svRate=svRate*10/max(svRate);
subplot(2,1,1)
plot(svRate)
subplot(2,1,2)
mdl=estimate(arima(1,0,2),svRate-0*mean(svRate));
mdl=estimate(arima('ARLags',1,'MALags',2),svRate-0*mean(svRate));


I've posted this to

• Why are you missing the first four years of data? The series in Figure 7.1 in the book seems close to this series of personal saving as a percentage of disposable personal income for the period 1955-1979. – javlacalle Jun 5 '15 at 21:11
• Thanks, javlacalle. I modified the question after using your more complete data set. – StatSmartWannaB Jun 5 '15 at 22:33

This is what I get fitting the following model in R:

$$y_t - \mu = \phi (y_{t-1} - \mu) + \epsilon_t + \theta \epsilon_{t-2} \,,$$

x <- structure(c(9.3, 9.5, 10, 9.9, 10.6, 11.1, 11.3, 11.5, 11.1, 11.6, 11.4, 10.8, 11.3, 11.2, 11.6, 11.5, 10.7, 10.7, 9.6, 10.1,
10.3, 9.7, 10.2, 10, 10.9, 10.9, 11.7, 11.6, 11.5, 11.3, 11, 10.6, 10.7, 10.6, 10.3, 11, 11, 11.8, 11.3, 12, 11.1, 11.1, 11.9, 11.3, 10.9, 10.9, 10.9, 11.6, 12.3, 11.8, 12.2, 12.4, 11.8, 11.9, 10.5, 10.6, 9.8, 10.1, 11.6, 11.5, 11.8, 12.5, 13.1, 13.1, 13.2, 13.6, 13.4, 12.9, 12.2, 11.4, 11.7, 13.2, 12.2, 13, 13, 14.1, 13.6, 12.5, 12.2, 13.3, 12.4, 15, 12.5, 12.3, 11.7, 11.3, 11.1,
10.4, 9.4, 10, 10.5, 10.8, 10.9, 9.9, 10.1, 10.1, 10.6, 9.8, 9.4, 9.5),
.Dim = c(100L, 1L), .Dimnames = list(NULL, "V1"), .Tsp = c(1959, 1983.75, 4), class = "ts")
arima(x, order=c(1,0,2), include.mean=TRUE, fixed=c(NA,0,NA,NA))
# Coefficients:
#          ar1  ma1     ma2  intercept
#       0.7643    0  0.2839    11.1904
# s.e.  0.0729    0  0.1092     0.3390
# sigma^2 estimated as 0.4087:  log likelihood = -97.83,  aic = 203.65


The term reported as intercept is actually the mean $\mu$. The constant $\alpha$ in the model below can be recovered from the relationship $E(y_t)\equiv \mu = \frac{\alpha}{1-\phi}$.

$$y_t = \alpha + \phi y_{t-1} + \epsilon_t + \theta \epsilon_{t-2} \,.$$

I get $\hat{\alpha} = 11.1904 \times (1 - 0.7643) = 2.6376$. The difference with the constant reported in the book may be due to the different scale of the data.

These results are close to what you get and not very far from the results in the book except for the MA(2) coefficient.

As regards the variance, be aware that the standard deviation $\hat{\sigma}_a$ is reported in the book (not the variance). I get $\sqrt{0.4087}=0.6393$, close to your result and the value reported in the book.

• You show yt=α+ϕyt−1+ϵt+θϵt−1 but refer to MA(2). – IrishStat Jun 6 '15 at 10:07
• @IrishStat I have fixed a typo in the equations, thanks! – javlacalle Jun 6 '15 at 10:32
• Thank you, javlacalle. Looks like the only mystery is the MA2 coefficient. Your correction about my original variance actually being standard deviation has cleared up that discrepancy. Now, the constant term -- as per the question, Pankratz's seems to match the variation of mine that scales the data to match the peak of 10%, as in his book. So maybe it's just an editing error. – StatSmartWannaB Jun 7 '15 at 2:35
• I find it interesting that R's "intercept" doesn't correspond to the constant term in the regression, but as long as it's explained in the help. Isn't it odd how our close results can still be off by so much? While I haven't done much linear regression, I was under the impression that the procedure was very standard. – StatSmartWannaB Jun 7 '15 at 2:35
• Slight differences may arise across software implementations due to different choices in the initialization of the Kalman filter (if used), optimization algorithm, possible transformation of parameters to ensure the AR polynomial is stationary,... But, in general, I wouldn't expect relevant discrepancies in the results. Although the results from Matlab and R are close to each other, it seems that there are some differences in the implementations that would be interesting to explore. – javlacalle Jun 7 '15 at 13:14