ARIMAX and VAR Models

The Federal Funds Rate as a Regressor and as a Variable in a System

Published

December 3, 2023

Modified

September 15, 2026

Purpose

This section adds the Federal Funds Rate to the univariate models of the previous section, in two ways. The first two models are regressions with ARMA errors, the ARIMAX form: a target series is regressed on the Federal Funds Rate and the regression residuals are modelled as an ARMA process, so that the coefficient on the policy rate is estimated with the autocorrelation of the target accounted for. The targets are monthly CPI inflation and the monthly change in the 30-year Treasury yield. The last two models are vector autoregressions, in which the Federal Funds Rate is one variable in a system and every variable is explained by the lagged values of all of them. The first VAR covers the Treasury yield curve; the second covers output, prices, unemployment and the mortgage rate.

The two approaches differ in what they assume about the policy rate. A regression treats it as exogenous, set without regard to the target. That is defensible for the 30-year yield, which the Federal Reserve does not target, and it is not defensible for inflation, which the Federal Reserve responds to. The inflation regression therefore uses the policy rate lagged twelve months, so that the regressor is fixed before the target is observed, and its coefficient is read as a forecasting relationship rather than a causal one. The VARs make no exogeneity assumption; their Granger-causality tests and impulse responses are the results in this section that speak to the effect of the policy rate. Those two tools extend the course’s VAR material, which covers lag selection, estimation, residual diagnostics and forecasting; each is explained where it is first used.

Each regression is built in the order the course sets out: the regression is fitted by least squares, the autocorrelation and partial autocorrelation of its residuals are inspected, auto.arima is run on those residuals to suggest an error structure, and the regression and the error model are then estimated together with Arima(..., xreg = ). A grid of ARMA orders ranked by AIC, AICc and BIC is reported alongside the auto.arima choice (Hyndman and Khandakar 2008), followed by coefficient tables, Ljung-Box tests and rolling-origin cross-validation. For the ARIMAX models the cross-validation compares each model with the same ARMA model without the regressor, which is the test of whether the policy rate adds forecasting information. The VARs are estimated with the vars package (Pfaff 2008).

ARIMAX: Inflation and the Federal Funds Rate

Data

Inflation is the monthly change in the log of CPIAUCSL, multiplied by 100, so a value of 0.30 is a 0.30% rise in the index over the month, about 3.6% at an annual rate. The regressor is the effective Federal Funds Rate, FEDFUNDS, lagged twelve months. Both series are monthly; after the lag the sample runs from July 1955 to August 2023, 818 observations.

Stationarity

Monthly inflation                  ADF statistic  -3.73   p-value 0.023

Federal Funds Rate, level          ADF statistic  -3.12   p-value 0.105

Federal Funds Rate, change         ADF statistic  -8.37   p-value < 0.01

Monthly inflation                  KPSS statistic   1.33   p-value < 0.01

Monthly inflation rejects the unit root, with an ADF statistic of −3.73 and p = 0.023, and is modelled in levels. It is persistent: the ACF is 0.61 at lag one and stays near 0.40 through lag eleven, and the KPSS test rejects level stationarity, which is the usual verdict on a series with a slowly shifting mean. The Federal Funds Rate does not reject the unit root in levels on this sample (statistic −3.12, p = 0.105) and is used as a lagged level in the regression, where the question is whether its past value helps predict inflation; the change in the rate is used as an alternative regressor below. The ACF of inflation also shows a spike at lag twelve, 0.30, which is what remains of the seasonal pattern after the Bureau of Labor Statistics’ adjustment; seasonal autoregressive terms are included in the grid for that reason.

Regression and Its Errors

The model is built in the order the course sets out: fit the regression first, look at what its residuals do, and let that choose the ARMA structure for the errors. The regression of monthly inflation on the Federal Funds Rate lagged twelve months is estimated by ordinary least squares, and its residuals are the series whose autocorrelation has to be modelled.

Least-Squares Regression of Monthly Inflation on the FFR Lagged 12 Months, R² = 0.111
Coefficient	Estimate	Std. Error	t value	Pr(>\|t\|)
(Intercept)	0.1673	0.0166	10.1044	0
ffr_l12	0.0285	0.0028	10.0689	0

auto.arima on the regression residuals:

Series: res1 
ARIMA(2,0,4)(2,0,0)[12] with zero mean 

Coefficients:
         ar1      ar2      ma1     ma2     ma3     ma4    sar1    sar2
      1.8141  -0.8236  -1.3852  0.2574  0.0712  0.1146  -0.142  -0.110
s.e.  0.0855   0.0832   0.0907  0.0775  0.0632  0.0365   0.037   0.037

sigma^2 = 0.05295:  log likelihood = 44.49
AIC=-70.98   AICc=-70.75   BIC=-28.61

The regression coefficient is 0.029 with an R² of 0.11, but its residuals are far from white noise: the ACF starts at 0.57 and decays slowly, and the PACF has a large first spike. Ordinary least squares is therefore not a valid summary of this relationship on its own, which is the reason for the ARMA error term. Applied to these residuals, auto.arima selects an ARMA(2,4)(2,0,0)[12] structure; the two seasonal autoregressive terms confirm that the annual pattern survives the Bureau’s seasonal adjustment. The next step estimates the regression and an error structure together, searching over orders rather than taking this one as given.

Model Selection

The grid covers ARMA error orders p and q from 0 to 3 and seasonal autoregressive order P from 0 to 2, 48 models, each with the twelve-month lag of the Federal Funds Rate as regressor.

Six Lowest-AICc Error Structures, Regressor FFR Lagged 12 Months
p	q	P	AIC	AICc	BIC
3	2	2	-73.20	-72.93	-26.14
1	3	2	-68.61	-68.39	-26.25
2	1	2	-68.18	-68.01	-30.53
3	1	2	-67.79	-67.57	-25.43
2	3	2	-67.83	-67.56	-20.76
2	2	2	-67.29	-67.07	-24.93

Three Lowest-BIC Error Structures
p	q	P	AIC	AICc	BIC
2	1	2	-68.18	-68.01	-30.53
1	2	2	-66.60	-66.42	-28.94
2	1	1	-60.12	-59.98	-27.18

auto.arima, same regressor:

Series: infl 
Regression with ARIMA(3,0,2)(2,0,0)[12] errors 

Coefficients:
         ar1      ar2     ar3      ma1     ma2     sar1     sar2  intercept
      2.1906  -1.5482  0.3521  -1.7627  0.8077  -0.1495  -0.1131     0.2553
s.e.  0.0814   0.1176  0.0454   0.0754  0.0675   0.0373   0.0371     0.0669
      ffr_l12
       0.0094
s.e.   0.0099

sigma^2 = 0.05273:  log likelihood = 46.6
AIC=-73.2   AICc=-72.93   BIC=-26.14

AICc selects ARMA(3,2)(2,0,0)[12] errors; BIC selects ARMA(2,1)(2,0,0)[12], 4.4 points behind on AIC with two fewer parameters. Every model in the top six has two seasonal autoregressive terms. auto.arima with the same regressor selects the same model as AICc. The BIC model is retained; the AICc model is shown for comparison.

Diagnostics

AIC -68.18   AICc -68.01   BIC -30.53   residual s.d. 0.2306
Ljung-Box p-value at lag 8: 0.314   Ljung-Box p-value at lag 12: 0.009   Ljung-Box p-value at lag 24: 0.001

Coefficient	Estimate	Std. Error	z	p
ar1	1.2927	0.0443	29.15	0.0000
ar2	-0.3048	0.0419	-7.27	0.0000
ma1	-0.8726	0.0251	-34.75	0.0000
sar1	-0.1244	0.0356	-3.49	0.0005
sar2	-0.1151	0.0361	-3.19	0.0014
intercept	0.2652	0.0760	3.49	0.0005
ffr_l12	0.0079	0.0092	0.86	0.3878

AIC -73.20   AICc -72.93   BIC -26.14   residual s.d. 0.2296
Ljung-Box p-value at lag 8: 0.129   Ljung-Box p-value at lag 12: 0.183   Ljung-Box p-value at lag 24: 0.013

Coefficient	Estimate	Std. Error	z	p
ar1	2.1906	0.0814	26.91	0.0000
ar2	-1.5482	0.1176	-13.17	0.0000
ar3	0.3521	0.0454	7.76	0.0000
ma1	-1.7627	0.0754	-23.39	0.0000
ma2	0.8077	0.0675	11.96	0.0000
sar1	-0.1495	0.0373	-4.01	0.0001
sar2	-0.1131	0.0371	-3.05	0.0023
intercept	0.2553	0.0669	3.81	0.0001
ffr_l12	0.0094	0.0099	0.96	0.3387

AIC -96.06   AICc -95.88   BIC -58.41   residual s.d. 0.2268
Ljung-Box p-value at lag 8: 0.259   Ljung-Box p-value at lag 12: 0.014   Ljung-Box p-value at lag 24: 0.002

Coefficient	Estimate	Std. Error	z	p
ar1	1.2740	0.0438	29.11	0.0000
ar2	-0.2900	0.0406	-7.14	0.0000
ma1	-0.8854	0.0245	-36.21	0.0000
sar1	-0.1291	0.0357	-3.62	0.0003
sar2	-0.1137	0.0362	-3.14	0.0017
intercept	0.1093	0.0545	2.00	0.0450
ffr_l1	0.0416	0.0072	5.82	0.0000

The coefficient on the Federal Funds Rate lagged twelve months is 0.008 with a standard error of 0.009 in the retained model and 0.009 with a standard error of 0.010 in the larger one: the policy rate a year earlier has no measurable relationship with this month’s inflation once inflation’s own dynamics are in the model. The third table replaces the twelve-month lag with a one-month lag. That coefficient is 0.042 with a standard error of 0.007, positive and significant, and the AIC falls by 28 points. The sign is the opposite of what monetary transmission predicts. It is the reaction function: the Federal Reserve raises the rate when inflation is high, so last month’s rate is high when this month’s inflation is high. A single-equation regression cannot separate the two directions, which is why the VAR later in this section is the model used for the causal question.

The ARMA terms are all significant. In the retained model the Ljung-Box p-values are 0.31, 0.009 and 0.001 at lags 8, 12 and 24; in the AICc model they are 0.13, 0.18 and 0.013. Neither set of residuals is white at the annual lags, though the autocorrelations involved are small. The residual standard deviation is 0.23 points of monthly inflation, against an unconditional standard deviation of 0.31.

Cross-Validation

Initial training set: 572 observations. Forecast origins: 246.

Model	MAE h = 1	RMSE h = 1	RMSE h = 3	RMSE h = 6	RMSE h = 12
ARMA(2,1)(2,0,0)[12], no regressor	0.1953	0.2877	0.3287	0.3271	0.3347
Regressor FFR lagged 12 months	0.1954	0.2879	0.3291	0.3280	0.3357
Regressor FFR lagged 1 month	0.2011	0.2906	0.3305	0.3324	0.3418
ARMA(2,1), no seasonal terms	0.1990	0.2900	0.3307	0.3303	0.3412

The initial training set ends in February 2003 and the models are re-estimated at each of the 246 later months, forecasting one to twelve months ahead. The twelve-month lag of the policy rate is known at every horizon, so its forecasts are genuine; the one-month lag is known only one step ahead, and its longer-horizon forecasts use the realised rate. Adding the twelve-month lag changes the one-step RMSE from 0.2877 to 0.2879 and the twelve-step RMSE from 0.335 to 0.336: no gain. The one-month lag, despite its in-sample significance, makes the forecasts worse at every horizon, 0.291 one step ahead and 0.342 at twelve. The seasonal terms are worth about 0.002 of RMSE. For forecasting monthly inflation, the policy rate adds nothing to what inflation’s own history provides.

Forecast

Because the regressor is lagged twelve months, the values it takes over the forecast period are the Federal Funds Rates of September 2022 to August 2023, which are known, so no assumption about future policy is needed. The forecast runs from September 2023 to August 2024. The point forecasts average 0.38% a month, 4.5% at an annual rate, against 3.6% over the twelve months to August 2023; the model carries forward the persistence of the 2021–2023 readings. The 80% interval at twelve months is 0.02 to 0.73% a month, an annual range of roughly 0 to 9%.

ARIMAX: The 30-Year Yield and the Federal Funds Rate

Data

The target is the monthly change in the 30-year constant-maturity yield, DGS30, with the daily series averaged within each month. The regressor is the monthly change in the Federal Funds Rate. The sample runs from March 1977 to August 2023, 557 monthly changes. Both are in % points.

Stationarity

30-year yield, level               ADF statistic  -3.65   p-value 0.028

30-year yield, monthly change      ADF statistic  -7.38   p-value < 0.01

Federal Funds Rate, monthly change ADF statistic  -6.34   p-value < 0.01

Correlation of levels: 0.904.  Correlation of monthly changes: 0.299.  Correlation of the yield change with the previous month's FFR change: 0.066.

The two levels are I(1) and move together over the sample, with a correlation of 0.90 that reflects the common decline from 1981 to 2020 rather than a month-to-month relationship. The monthly changes are stationary, with ADF statistics of −7.38 and −6.34, and their correlation is 0.30. The change in the yield has a lag-one autocorrelation of 0.32 and a PACF that is significant at lags 1 and 2, so an ARMA error structure is needed.

Regression and Its Errors

Least-Squares Regression of the Change in the 30-Year Yield on the Change in the FFR, R² = 0.090
Coefficient	Estimate	Std. Error	t value	Pr(>\|t\|)
(Intercept)	-0.0065	0.0107	-0.6048	0.5456
dffr	0.1494	0.0202	7.3886	0.0000

auto.arima on the regression residuals:

Series: res2 
ARIMA(0,0,2) with zero mean 

Coefficients:
         ma1      ma2
      0.3195  -0.1435
s.e.  0.0422   0.0418

sigma^2 = 0.05692:  log likelihood = 8.76
AIC=-11.53   AICc=-11.48   BIC=1.44

The least-squares coefficient is 0.149 with an R² of 0.09. Its residuals have a first-order autocorrelation of about 0.3 and a PACF that cuts off after two lags, so an ARMA error term is again needed; auto.arima on the residuals selects MA(2). Both series are already differenced, which is what the course requires before regressing one non-stationary series on another: fitting the levels against each other would risk a spurious regression.

Model Selection

Six Lowest-AICc Error Structures, Regressor Change in FFR
p	q	AIC	AICc	BIC
0	2	-8.65	-8.54	12.96
1	1	-7.33	-7.22	14.29
1	2	-6.75	-6.60	19.18
0	3	-6.74	-6.59	19.19
2	1	-6.49	-6.34	19.44
3	3	-6.67	-6.34	32.23

auto.arima, same regressor:

Series: dy30 
Regression with ARIMA(0,0,2) errors 

Coefficients:
         ma1      ma2    dffr
      0.3319  -0.1327  0.1262
s.e.  0.0437   0.0429  0.0217

sigma^2 = 0.05694:  log likelihood = 9.18
AIC=-10.37   AICc=-10.29   BIC=6.92

AIC, AICc and BIC all select MA(2) errors, and auto.arima agrees. Adding the previous month’s change in the Federal Funds Rate as a second regressor raises AICc by 1.9 points and its coefficient is −0.006 with p = 0.78, so it is not retained.

Diagnostics

AIC -8.65   AICc -8.54   BIC 12.96   residual s.d. 0.2388
Ljung-Box p-value at lag 8: 0.575   Ljung-Box p-value at lag 12: 0.849   Ljung-Box p-value at lag 24: 0.234

Coefficient	Estimate	Std. Error	z	p
ma1	0.3312	0.0438	7.57	0.0000
ma2	-0.1334	0.0429	-3.11	0.0019
intercept	-0.0065	0.0121	-0.54	0.5925
dffr	0.1264	0.0217	5.82	0.0000

AIC 21.65   AICc 21.73   BIC 38.95   residual s.d. 0.2456
Ljung-Box p-value at lag 8: 0.196   Ljung-Box p-value at lag 12: 0.295   Ljung-Box p-value at lag 24: 0.171

Coefficient	Estimate	Std. Error	z	p
ma1	0.4006	0.0428	9.37	0.0000
ma2	-0.0715	0.0426	-1.68	0.0929
intercept	-0.0063	0.0138	-0.46	0.6481

The coefficient on the change in the Federal Funds Rate is 0.126 with a standard error of 0.022: a one-point rise in the policy rate within a month is associated with a 0.13-point rise in the 30-year yield in the same month, holding the yield’s own dynamics fixed. The 30-year yield is priced from expectations of short rates over three decades and a term premium, so a pass-through of about an eighth is the expected order of magnitude. The moving-average terms are 0.33 and −0.13, both significant, and the residuals pass the Ljung-Box test at lags 8, 12 and 24 (p = 0.58, 0.85 and 0.23). Removing the regressor raises the AIC by 30 points. The residual standard deviation is 0.24 points a month.

Cross-Validation

Initial training set: 389 observations. Forecast origins: 168.

Model	MAE h = 1	RMSE h = 1	RMSE h = 3	RMSE h = 6	RMSE h = 12
MA(2), no regressor	0.1282	0.1628	0.1763	0.1768	0.1764
Regressor change in FFR	0.1268	0.1610	0.1721	0.1725	0.1720
Regressors change in FFR and its lag	0.1268	0.1610	0.1722	0.1725	0.1720
No change (random walk in the level)	0.1342	0.1757	0.1763	0.1768	0.1764

The initial training set ends in July 2009 and there are 168 forecast origins. The regressor’s forecasts are conditional on the realised change in the policy rate, so they measure how much knowing the Federal Reserve’s move would have helped, not what a forecaster could have achieved. On that basis the regressor lowers the one-step RMSE from 0.163 to 0.161 and the twelve-step RMSE from 0.176 to 0.172, a gain of 1 to 2%. The random walk, which forecasts no change, has a one-step RMSE of 0.176. The lagged regressor adds nothing. The improvement from the policy rate is real and small, consistent with the coefficient of 0.13 and a correlation of 0.30 between the two monthly changes.

Forecast

A forecast of the yield change requires an assumption about the policy rate. The one shown holds the Federal Funds Rate at its August 2023 level, so the regressor is zero throughout. Under that assumption the forecast change is within 0.01 of zero from the third month, and the 80% interval is about ±0.3 points a month: the model’s information about the direction of long yields is exhausted within two months, after which the forecast is the level as of August 2023, 4.28%.

VAR: The Federal Funds Rate and the Treasury Curve

Data

Seven monthly series from January 1982 to August 2023, 500 months: the effective Federal Funds Rate and Treasury yields at 3 months, 6 months, 1 year, 5 years, 10 years and 30 years. The four short-maturity yields are reconstructed from the spread series described on the data page by adding the Federal Funds Rate back to the monthly average spread; since the monthly FEDFUNDS value is itself the average of the daily effective rate, the sum is the monthly average yield up to the difference between trading-day and calendar-day averaging. The 10- and 30-year yields are monthly averages of the daily constant-maturity series. All are in %.

ADF Statistics and Correlation with the Federal Funds Rate
Series	Level	First difference	Correlation with FFR
FFR	-2.13	-6.65	1.000
T3M	-1.95	-7.06	0.996
T6M	-2.30	-7.01	0.994
T1Y	-2.40	-6.80	0.990
T5Y	-2.83	-7.47	0.955
T10Y	-3.37	-7.64	0.927
T30Y	-3.66	-7.71	0.896

Every series fails to reject the unit root in levels at the 5% level except the 30-year yield (statistic −3.66), and every first difference rejects it. A VAR is fitted to the levels rather than the differences. This follows the treatment of the mortality and pollution example in the course, where the model is fitted to data that need not be stationary and the deterministic terms are estimated inside the VAR rather than by de-trending each series first. The deterministic term here is a constant only, not the constant-plus-trend used in that example: interest rates have no deterministic trend, and a linear trend term would extrapolate the 1982–2020 decline into the forecast.

Lag Selection

Lag Order Selected by Each Criterion, Maximum 12
Criterion	Selected lag
AIC	2
HQ	2
SC (BIC)	2
FPE	2

AIC, HQ and SC by Lag Order, Lags 1 to 6
	1	2	3	4	5	6
AIC(n)	-33.23	-33.73	-33.68	-33.64	-33.58	-33.59
HQ(n)	-33.04	-33.38	-33.16	-32.95	-32.73	-32.58
SC(n)	-32.75	-32.83	-32.36	-31.90	-31.42	-31.01

All four criteria select two lags. VAR(2) has 7 equations with 15 coefficients each. The largest root of its companion matrix is 0.978, inside the unit circle, and the equations’ R² are all above 0.99, which is the persistence of the levels rather than explanatory power.

VAR(2): 7 equations, 15 coefficients each, 498 observations. Largest root 0.978.

Portmanteau test for residual autocorrelation to lag 24: statistic 1335.4, p < 0.001.

The portmanteau test rejects the absence of residual autocorrelation. With 49 cross-equation autocorrelations at each lag the test has power against small residual dependence, and the lag-selection criteria did not favour a longer lag; the result is noted and the model is used as selected.

Residual Diagnostics

Correlation Matrix of VAR Residuals
	FFR	T3M	T6M	T1Y	T5Y	T10Y	T30Y
FFR	1.000	0.693	0.659	0.596	0.338	0.292	0.252
T3M	0.693	1.000	0.935	0.850	0.597	0.493	0.430
T6M	0.659	0.935	1.000	0.961	0.733	0.626	0.538
T1Y	0.596	0.850	0.961	1.000	0.835	0.726	0.627
T5Y	0.338	0.597	0.733	0.835	1.000	0.947	0.848
T10Y	0.292	0.493	0.626	0.726	0.947	1.000	0.948
T30Y	0.252	0.430	0.538	0.627	0.848	0.948	1.000

The correlation matrix shows how much of each month’s movement the model leaves unexplained in common across equations. The residual correlations between neighbouring maturities are high, 0.95 between the 10- and 30-year yields and 0.94 between the 3- and 6-month, and they fall steadily with the distance between maturities, to 0.25 between the policy rate and the 30-year yield. The model does not capture what the curve does within a month, which is a limitation of a monthly VAR rather than a fault in the lag order, and it is also why the impulse responses below need an ordering assumption. In the plot the diagonal panels are the autocorrelation of each equation’s residuals and the off-diagonal panels the cross-correlations between pairs, with the second-named series leading. The diagonals are close to zero beyond lag one; the large spikes at lag zero in the off-diagonals are the contemporaneous correlations just described.

Granger Causality

The Granger test asks a forecasting question, not a causal one: does adding the past of one variable to a regression improve the prediction of another, beyond that variable’s own past? If it does, the first variable is said to Granger-cause the second. The test is an F-test on the lags of the candidate variable in one equation of the VAR, so it uses the model already estimated above. It extends the course’s VAR material, which stops at lag selection, residual diagnostics and forecasting.

The table reports two things. The first line is the joint test from the full VAR: all the other variables at once. The table below it fits a two-variable VAR for each pair at the same lag order and tests each direction separately, which is how the direction of the relationship becomes visible.

Joint test in the full VAR(2), FFR does not Granger-cause the other variables: F = 6.86, p < 0.001

Bivariate Granger Tests, VAR(2), p-values
Variable	FFR causes variable	Variable causes FFR
T3M	< 0.001	< 0.001
T6M	< 0.001	< 0.001
T1Y	0.012	< 0.001
T5Y	0.128	< 0.001
T10Y	0.326	< 0.001
T30Y	0.461	0.003

The joint test rejects: the Federal Funds Rate helps predict the yield curve. The bivariate tests show where. The policy rate Granger-causes the 3-month, 6-month and 1-year yields, with p-values below 0.02, and does not Granger-cause the 5-, 10- or 30-year yields, with p-values of 0.13, 0.33 and 0.46. In the other direction, every yield Granger-causes the Federal Funds Rate, with p-values below 0.003 at all six maturities. The pattern is what the term structure implies. Short yields are averages of expected policy rates over the next months, so the policy rate and the short yields predict each other. Long yields embody expectations of policy years ahead, so they anticipate the policy rate, but last month’s policy rate carries no information about the 30-year yield beyond what the 30-year yield already contains.

Impulse Responses

An impulse response traces what the model implies happens to every variable over the following months when one variable receives an unexpected one-period push. It is computed by feeding a one-unit disturbance into the estimated VAR and iterating the equations forward with no further disturbances, so it uses only the coefficients already estimated; the grey bands come from re-estimating the model on 200 bootstrap resamples. This tool is not part of the course’s VAR material and is included because it is the only way in this study to put a number on the size and timing of a policy-rate effect.

Because the variables move together within a month, a disturbance to one is not independent of the others. The responses below resolve that by ordering the variables and letting each affect only those after it within the month, which is the standard convention; the ordering is stated with each figure and is an assumption, not a result.

Response to a One-Standard-Deviation FFR Shock, % Points, and Share of Forecast-Error Variance Due to FFR Shocks
Series	h = 0	h = 3	h = 6	h = 12	h = 24	FEVD share at h = 12, %
FFR	0.208	0.291	0.272	0.229	0.171	53.5
T3M	0.163	0.276	0.251	0.210	0.157	43.8
T6M	0.155	0.263	0.249	0.213	0.161	39.4
T1Y	0.145	0.244	0.239	0.212	0.162	34.2
T5Y	0.086	0.159	0.176	0.173	0.143	21.2
T10Y	0.069	0.120	0.140	0.143	0.123	16.4
T30Y	0.053	0.093	0.110	0.113	0.102	12.6

The shock is orthogonalised with the Federal Funds Rate ordered first, so it is the part of the month’s movement in the policy rate not explained by the system’s past, and the yields are allowed to respond within the month. Its size is 0.21 points. The bands are 95% bootstrap intervals from 200 replications.

The response is largest at the short end and declines with maturity at every horizon. On impact the 3-month yield rises 0.16 points, three-quarters of the shock, the 5-year 0.09 and the 30-year 0.05. The responses build for several months as the policy rate itself continues to rise, since the Federal Reserve moves in sequences: the policy rate peaks at 0.29 points above baseline after three months. At twelve months the 3-month yield is 0.21 points above baseline, 92% of the policy rate’s own response, and the 30-year is 0.11, 49%. The 30-year response is still significant at 24 months. Policy-rate shocks account for 54% of the twelve-month forecast-error variance of the policy rate, 44% for the 3-month yield, and 13% for the 30-year, whose variance is dominated by its own shocks.

Cross-Validation

Initial training set: 350 observations. Forecast origins: 150.

One-Step-Ahead RMSE by Variable
Model	1	2	3	4	5	6	7
Random walk	0.163	0.171	0.170	0.166	0.189	0.192	0.176
VAR(1)	0.130	0.157	0.163	0.166	0.203	0.201	0.183
VAR(2)	0.115	0.156	0.156	0.159	0.199	0.195	0.178
VAR(4)	0.117	0.151	0.154	0.162	0.206	0.203	0.185

RMSE Averaged Across Variables, by Horizon
Model	h = 1	h = 3	h = 6	h = 12
Random walk	0.175	0.407	0.698	1.181
VAR(1)	0.172	0.393	0.647	1.117
VAR(2)	0.165	0.398	0.666	1.154
VAR(4)	0.168	0.401	0.665	1.160

The initial training set ends in February 2011 and there are 150 forecast origins. VAR(2) beats the random walk one step ahead for the policy rate, 0.115 against 0.163, and for the yields out to one year, 0.156 to 0.159 against 0.166 to 0.171; for the 5-, 10- and 30-year yields it is no better than the random walk, 0.18 to 0.20 for both. Averaged across variables the VAR’s advantage over the random walk is 6% one step ahead and 2% at twelve months. VAR(4) is within 0.01 of VAR(2) at every horizon; VAR(1) is slightly worse one step ahead and slightly better from three months on, 1.117 against 1.154 at twelve. The gains are concentrated where the Granger tests placed them, at the short end.

Forecast

From the August 2023 values, 5.33% for the policy rate and 4.28% for the 30-year, the VAR(2) forecasts the policy rate at 5.59% after twelve months and 5.03% after 24, the 3-month yield at 5.25% and 4.73%, and the 30-year at 5.03% and 5.15%. The 95% interval for the policy rate at twelve months runs from 3.2% to 8.0%. The forecast has the curve flattening: the short yields fall and the long yields rise, closing most of the inversion in the August 2023 data.

VAR: The Federal Funds Rate and the Macroeconomy

Data

Five quarterly series from the second quarter of 1971 to the second quarter of 2023, 209 quarters, the sample set by the start of the mortgage-rate series. Real GDP, GDPC1, and the CPI, CPIAUCSL, enter as 100 times the natural log, so that changes are in %; the Federal Funds Rate, the unemployment rate and the 30-year mortgage rate are quarterly averages of the monthly or weekly data, in %. Nothing is interpolated: the monthly series are averaged to the quarterly frequency of GDP. The variables are ordered output, prices, unemployment, policy rate, mortgage rate.

Summary Statistics and ADF Statistics
Series	Mean	Std. dev.	Level	First difference
LGDP	937.60	41.63	-1.48	-6.04
LCPI	494.35	53.93	-3.74	-3.65
UNRATE	6.17	1.70	-3.57	-6.05
FFR	4.87	3.94	-4.64	-5.17
MORT	7.74	3.26	-2.82	-5.71

Log GDP and the mortgage rate do not reject the unit root in levels; log CPI, unemployment and the policy rate reject at the 5% level on this quarterly sample, log CPI because the ADF regression includes a trend. The system is estimated in levels with a constant, for the reasons given for the yield-curve VAR, and because differencing log GDP and log CPI would discard the level relationships among output, prices and rates that the impulse responses are meant to trace.

Lag Selection

Lag Order Selected by Each Criterion, Maximum 8
Criterion	Selected lag
AIC	2
HQ	2
SC (BIC)	2
FPE	2

AIC, HQ and SC by Lag Order
	1	2	3	4	5	6	7	8
AIC(n)	-4.94	-5.46	-5.43	-5.35	-5.20	-5.33	-5.22	-5.13
HQ(n)	-4.75	-5.09	-4.90	-4.65	-4.33	-4.30	-4.02	-3.77
SC(n)	-4.45	-4.55	-4.12	-3.62	-3.06	-2.78	-2.26	-1.76


VAR(2): 5 equations, 11 coefficients each, 207 observations. Largest root 0.995. Portmanteau test to lag 16: p = 0.324.

All four criteria select two lags. The largest root is 0.995, close to but inside the unit circle, and the portmanteau test does not reject residual whiteness (p = 0.32). The Federal Funds Rate equation has significant coefficients on its own two lags, 1.17 and −0.36, and on the second lag of the mortgage rate; the mortgage-rate equation has significant coefficients on the policy rate at both lags, 0.17 and −0.12, and on the price level at both lags.

Residual Diagnostics

Correlation Matrix of VAR Residuals
	LGDP	LCPI	UNRATE	FFR	MORT
LGDP	1.000	0.304	-0.801	0.278	0.155
LCPI	0.304	1.000	-0.286	0.317	0.345
UNRATE	-0.801	-0.286	1.000	-0.248	-0.094
FFR	0.278	0.317	-0.248	1.000	0.528
MORT	0.155	0.345	-0.094	0.528	1.000

The residual correlations here are smaller than in the yield-curve system, since quarterly macroeconomic series move less in lockstep within a period than interest rates do within a month. The largest is −0.80 between output and unemployment, the within-quarter form of the relationship between growth and joblessness, followed by 0.53 between the policy rate and the mortgage rate.

Granger Causality

The same test is applied to the macroeconomic system: whether the past of the policy rate improves the forecast of each other variable, and whether each improves the forecast of the policy rate.

Joint test in the full VAR(2), FFR does not Granger-cause the other variables: F = 5.11, p < 0.001

Bivariate Granger Tests, VAR(2), p-values
Variable	FFR causes variable	Variable causes FFR
LGDP	< 0.001	0.023
LCPI	0.002	0.052
UNRATE	0.004	0.647
MORT	< 0.001	0.163

The policy rate Granger-causes every other variable in the system, with bivariate p-values below 0.001 for output and the mortgage rate, 0.002 for prices and 0.004 for unemployment. In the other direction, output Granger-causes the policy rate (p = 0.023) and prices do so at the 10% level (p = 0.052); unemployment and the mortgage rate do not. This is the reaction function seen from the other side: the Federal Reserve’s rate responds to lagged output and inflation, and the rate in turn predicts all four.

Impulse Responses

The same exercise for the macroeconomic system. The variance decomposition in the last column of the table reports the share of each variable’s forecast error at that horizon attributable to policy-rate disturbances.

Response to a One-Standard-Deviation FFR Shock and Share of Forecast-Error Variance Due to FFR Shocks
Series	h = 0	h = 2	h = 4	h = 8	h = 12	h = 20	FEVD share at h = 8, %
LGDP	0.000	0.002	-0.275	-0.489	-0.530	-0.440	11.6
LCPI	0.000	0.299	0.487	0.628	0.633	0.512	17.3
UNRATE	0.000	-0.063	0.068	0.159	0.163	0.100	3.1
FFR	0.812	0.809	0.573	0.315	0.180	0.062	63.4
MORT	0.189	0.371	0.371	0.311	0.244	0.140	41.7

The shock is orthogonalised with the policy rate ordered after output, prices and unemployment and before the mortgage rate: the Federal Reserve is allowed to respond to the quarter’s macroeconomic data, the mortgage rate is allowed to respond to the quarter’s policy rate, and output, prices and unemployment respond only from the next quarter. The shock is 0.88 points, and the policy rate stays 0.57 points above baseline after a year and 0.18 after three years.

Three of the four responses have the sign monetary theory predicts. Real GDP is unchanged for two quarters, then falls, to 0.49% below baseline after two years and 0.53% after three, with the 95% band excluding zero from the third quarter. Unemployment rises by 0.16 points at two to three years, with the band excluding zero from the sixth quarter. The mortgage rate rises 0.19 points on impact and 0.37 at two quarters, about 40% of the policy move, then decays with it. The price level is the exception: it rises, by 0.63% after two years, with a band that excludes zero. A positive price response to a tightening is a known result in this literature, called the price puzzle, and the usual explanation is that the policy rate rises when the Federal Reserve expects inflation that a five-variable system does not otherwise see, so the estimated shock is partly anticipation of inflation rather than policy alone. A measure of expected inflation is the standard remedy; this data set does not contain one. Policy-rate shocks account for 63% of the two-year forecast-error variance of the policy rate, 42% for the mortgage rate, 17% for prices, 12% for output and 3% for unemployment.

Cross-Validation

Initial training set: 146 observations. Forecast origins: 63.

One-Step-Ahead RMSE by Variable
Model	1	2	3	4	5
Random walk	1.613	0.934	1.344	0.462	0.360
VAR(1)	1.596	0.726	1.339	0.497	0.386
VAR(2)	2.196	0.703	1.690	1.022	0.386
VAR(3)	2.597	0.895	1.851	1.427	0.410
VAR(4)	2.647	0.983	1.903	1.549	0.456

RMSE Averaged Across Variables, by Horizon
Model	h = 1	h = 2	h = 4	h = 8
Random walk	0.943	1.408	2.136	3.199
VAR(1)	0.909	1.316	1.883	2.566
VAR(2)	1.199	1.829	2.526	2.873
VAR(3)	1.436	2.039	2.585	2.893
VAR(4)	1.507	2.060	2.644	3.044

The initial training set ends in the third quarter of 2007 and there are 63 forecast origins, so the test period is the financial crisis, the zero-rate years, the pandemic and the 2022 tightening. On this period VAR(1) forecasts better than VAR(2) on average at every horizon and for every variable but the price level one step ahead, and the random walk forecasts the policy rate better than any VAR: one step ahead its RMSE is 0.46 points against 0.50 for VAR(1) and 1.02 for VAR(2). The VARs’ advantage is in output and prices at longer horizons, where VAR(1) has a two-year RMSE of 3.0% for log GDP against 4.7% for the random walk, and 3.4% against 5.4% for log CPI. The pattern says that on a test period dominated by a rate held at zero for seven years and a one-quarter fall in GDP of 8%, the lag-2 dynamics estimated on 1971–2007 do not carry over, and a simpler model is safer. VAR(2) is retained as the model the selection criteria chose on the full sample, and the cross-validation is the caveat on its forecasts.

Forecast

From the second quarter of 2023 the VAR(2) forecasts the policy rate falling from 4.99% to 4.06% after four quarters and 2.91% after eight, unemployment rising from 3.6% to 4.8% and 5.6%, the mortgage rate easing from 6.51% to 5.76%, real GDP flat over the two years, and the price level rising 2.5% a year. It is the profile of a slowdown in which the Federal Reserve cuts rates as unemployment rises. The cross-validation above shows this model’s policy-rate forecasts have been poor over the last fifteen years, and the intervals are wide: the 95% band for the policy rate at eight quarters runs from −1.7% to 7.5%.

Summary

Retained Models and Main Findings
Section	Frequency	Sample	Retained model	Effect of the FFR
ARIMAX, inflation	Monthly	1955–2023	ARMA(2,1)(2,0,0)[12] errors, FFR lagged 12 months	Coefficient 0.008 (s.e. 0.009); no forecast gain
ARIMAX, 30-year yield	Monthly	1977–2023	MA(2) errors, change in FFR	Coefficient 0.126 (s.e. 0.022); RMSE 1–2% lower
VAR, Treasury curve	Monthly	1982–2023	VAR(2), levels, constant	Granger-causes yields to 1 year; not 5 to 30 years
VAR, macroeconomy	Quarterly	1971–2023	VAR(2), levels, constant	Granger-causes all four; GDP −0.5%, unemployment +0.16 at 2–3 years

What the Models Found

The four models agree on where the policy rate’s influence is direct and where it is not. In the Treasury market, the effect runs through the short end. The policy rate Granger-causes the 3-month, 6-month and 1-year yields and not the longer ones; the impulse responses show the 3-month yield moving three-quarters of the way with the policy rate on impact and nine-tenths after a year, while the 30-year moves a quarter on impact and a half after a year; and the single-equation regression puts the same-month pass-through to the 30-year yield at 0.13. Every yield Granger-causes the policy rate, which is the market anticipating the Federal Reserve. In the macroeconomy the policy rate Granger-causes output, prices, unemployment and the mortgage rate; a 0.9-point tightening shock lowers real GDP by half a percent over two to three years and raises unemployment by 0.16 points, with the mortgage rate absorbing about 40% of the move within two quarters.

The inflation regression is the model that did not find an effect, and the reason is instructive. The policy rate lagged a year has a coefficient indistinguishable from zero and adds nothing to inflation forecasts; the policy rate lagged a month has a positive, significant coefficient and makes the forecasts worse. Both results are what a reaction function produces. The Federal Reserve raises the rate when inflation is high, so the rate is correlated with inflation in the wrong direction for a causal reading and in a direction that does not survive out of sample. The VAR, which lets both variables respond to each other, finds that the policy rate Granger-causes prices, but its impulse response has the price level rising after a tightening, the price puzzle. The macroeconomic effect of policy on inflation is the one question in this section the data set cannot answer cleanly, because the system lacks the variable that would let it separate the Federal Reserve’s information from its action.

Model Choice

The regression form suits the 30-year yield and does not suit inflation. The Federal Reserve does not set the policy rate in response to the 30-year yield, so treating the policy rate as exogenous to it is a reasonable approximation at monthly frequency, and the estimated pass-through is the kind of quantity a regression can deliver. The Federal Reserve does set the policy rate in response to inflation, and a regression of inflation on the policy rate estimates that response as much as it estimates transmission. The twelve-month lag makes the regressor predetermined but cannot make it exogenous, since the rate a year ago responded to inflation a year ago, which is correlated with inflation today. For questions of effect, the VAR is the right instrument in this data set. Two of the tools applied to the VARs, the Granger-causality test and the impulse response function, go beyond the vector-autoregression material covered in the course, which ends at lag selection, residual diagnostics and forecasting; both are explained where they appear, and the impulse responses rest on an ordering assumption that is stated with them.

The two VARs are estimated in levels with a constant, following the course’s treatment of a VAR on non-stationary data. Both lag selections settled on two lags, and in both the residuals carry large contemporaneous cross-correlations and either fail the Portmanteau test or lose to a simpler model out of sample, which is the usual position of a small VAR on fifty years of data that include several policy regimes. The macroeconomic VAR would be improved by a measure of expected inflation, which is the variable most often added when the price response comes out positive.

Forecasting

None of the four models forecasts well far ahead, and none should be expected to. The inflation model forecasts 4.5% for the year from September 2023, carrying forward the 2021–2023 persistence, against 3.6% in the year before. The 30-year model has no information beyond two months once the policy rate is held fixed. The yield-curve VAR beats a random walk by 6% one step ahead and 2% at a year, entirely at the short end. The macroeconomic VAR’s policy-rate forecasts were worse than a random walk over 2007–2023, and its eight-quarter projection, a policy rate of 2.9%, unemployment at 5.6% and flat real GDP, should be read as the pattern of a typical easing cycle in the sample rather than a prediction.

Limitations

The findings rest on Granger causality and orthogonalised impulse responses, both of which identify effects from the ordering of variables in time and within the period, not from any measure of the Federal Reserve’s intentions. The yield-curve VAR uses the effective rate, which moved with the target during the 1982–2008 period and sat at the floor of the target range from 2009 to 2015 and 2020 to 2022; a VAR on that sample treats a rate the Federal Reserve could not lower as an ordinary observation. The short-maturity yields are reconstructed from spreads and differ from the published constant-maturity yields by the averaging convention. The macroeconomic VAR has five variables and 209 quarters; with two lags each equation has 11 coefficients, which is adequate, but the bootstrap bands are wide at long horizons, and the cross-validation on the post-2007 period shows the estimated dynamics do not hold at the zero bound.

References

Hyndman, Rob J, and Yeasmin Khandakar. 2008. “Automatic Time Series Forecasting: The Forecast Package for R.” Journal of Statistical Software 27 (3): 1–22. https://doi.org/10.18637/jss.v027.i03.

Pfaff, Bernhard. 2008. “VAR, SVAR and SVEC Models: Implementation Within R Package vars.” Journal of Statistical Software 27 (4): 1–32. https://doi.org/10.18637/jss.v027.i04.

--- title: "ARIMAX and VAR Models" subtitle: "The Federal Funds Rate as a Regressor and as a Variable in a System" date: "2023-12-3" date-modified: today date-format: long format: html: smooth-scroll: true toc: true code-fold: true code-tools: true lightbox: true embed-resources: true bibliography: citation.bib execute: echo: false warning: false message: false --- ```{r setup} suppressMessages({ library(readxl); library(dplyr); library(lubridate); library(ggplot2) library(forecast); library(tseries); library(vars); library(knitr) }) knitr::opts_chunk$set(dev.args = list(bg = "#f4f4f4")) theme_set(theme_minimal(base_size = 11) + theme(panel.grid.minor = element_blank(), panel.grid.major = element_line(colour = "#e0e0e0", linewidth = 0.3), plot.background = element_rect(fill = "#f4f4f4", colour = NA), panel.background = element_rect(fill = NA, colour = NA), legend.background = element_rect(fill = NA, colour = NA), legend.key = element_rect(fill = NA, colour = NA), strip.background = element_rect(fill = NA, colour = NA), plot.title = element_text(size = 11, face = "plain"), axis.title = element_text(size = 9.5, colour = "#555555"), axis.text = element_text(size = 9, colour = "#555555"))) par(bg = "#f4f4f4") rd <- function(f) { d <- na.omit(read_excel(file.path("data/data", f))); data.frame(Date = as.Date(d[[1]]), v = as.numeric(d[[2]])) } mavg <- function(f) { d <- rd(f); d$Date <- floor_date(d$Date, "month"); d |> group_by(Date) |> summarize(v = mean(v), .groups = "drop") } qavg <- function(f) { d <- rd(f); d$Date <- floor_date(d$Date, "quarter"); d |> group_by(Date) |> summarize(v = mean(v), .groups = "drop") } fmt_p <- function(p) ifelse(p < 0.001, "< 0.001", formatC(p, format = "f", digits = 3)) adf_line <- function(x, label) { a <- suppressWarnings(adf.test(x)) cat(sprintf("%-34s ADF statistic %6.2f p-value %s\n", label, a$statistic, ifelse(a$p.value <= 0.01, "< 0.01", ifelse(a$p.value >= 0.99, "> 0.99", sprintf("%.3f", a$p.value))))) } # Coefficient table, information criteria and Ljung-Box p-values for an Arima fit. report <- function(fit, lags = c(8, 12, 24)) { se <- sqrt(diag(fit$var.coef)); z <- coef(fit) / se k <- sum(fit$arma[c(1, 2, 3, 4)]) tab <- data.frame(Coefficient = names(coef(fit)), Estimate = round(coef(fit), 4), `Std. Error` = round(se, 4), z = round(z, 2), p = formatC(2 * pnorm(-abs(z)), format = "f", digits = 4), check.names = FALSE) lb <- sapply(lags, function(L) Box.test(residuals(fit), lag = L, type = "Ljung-Box", fitdf = k)$p.value) cat(sprintf("AIC %.2f AICc %.2f BIC %.2f residual s.d. %.4f\n", fit$aic, fit$aicc, fit$bic, sqrt(fit$sigma2))) cat(sprintf("Ljung-Box p-value at lag %d: %.3f", lags, lb), sep = " "); cat("\n") kable(tab, row.names = FALSE, align = "lrrrr") } # Grid of ARMA(p,q)(P,0,0)[s] error structures for a regression on xreg. grid_x <- function(y, x, pmax = 3, qmax = 3, Pmax = 0) { g <- expand.grid(p = 0:pmax, q = 0:qmax, P = 0:Pmax) ic <- t(apply(g, 1, function(o) { m <- tryCatch(Arima(y, order = c(o["p"], 0, o["q"]), seasonal = c(o["P"], 0, 0), xreg = x), error = function(e) NULL) if (is.null(m)) c(NA, NA, NA) else c(m$aic, m$aicc, m$bic) })) g <- cbind(g, AIC = round(ic[, 1], 2), AICc = round(ic[, 2], 2), BIC = round(ic[, 3], 2)) if (Pmax == 0) g$P <- NULL g[order(g$AICc), ] } # Rolling-origin cross-validation for regression-with-ARMA-errors models. # From each origin the model is re-estimated on all earlier data and forecast # 1 to H steps ahead; the regressor's future values are taken as realised. cv_x <- function(y, X, specs, frac = 0.7, H = 12) { n <- length(y); k <- floor(n * frac); origins <- k:(n - 1) err <- lapply(specs, function(s) matrix(NA, length(origins), H)) for (ii in seq_along(origins)) { o <- origins[ii]; hh <- min(H, n - o) ytr <- ts(y[1:o], frequency = frequency(y)); yte <- y[(o + 1):(o + hh)] for (j in seq_along(specs)) { s <- specs[[j]] xtr <- if (is.null(s$x)) NULL else X[1:o, s$x, drop = FALSE] xte <- if (is.null(s$x)) NULL else X[(o + 1):(o + hh), s$x, drop = FALSE] f <- tryCatch(Arima(ytr, order = s$order, seasonal = s$seasonal, xreg = xtr), error = function(e) NULL) if (is.null(f)) next f$call$xreg <- xtr err[[j]][ii, 1:hh] <- forecast(f, h = hh, xreg = xte)$mean - yte } } nm <- sapply(specs, `[[`, "name") tab <- data.frame(Model = nm, `MAE h = 1` = sapply(err, function(e) round(mean(abs(e[, 1]), na.rm = TRUE), 4)), t(sapply(err, function(e) round(sqrt(colMeans(e^2, na.rm = TRUE))[c(1, 3, 6, 12)], 4))), check.names = FALSE) names(tab)[3:6] <- paste0("RMSE h = ", c(1, 3, 6, 12)) cat(sprintf("\n\nInitial training set: %d observations. Forecast origins: %d.\n\n", k, length(origins))) print(kable(tab, row.names = FALSE)); cat("\n\n") } # Rolling-origin cross-validation for VARs against a random walk. cv_var <- function(Y, ps, frac = 0.7, H = 12, hs = c(1, 3, 6, 12)) { n <- nrow(Y); k <- floor(n * frac); origins <- k:(n - 1) res <- lapply(ps, function(p) array(NA, c(length(origins), H, ncol(Y)))) for (ii in seq_along(origins)) { o <- origins[ii]; hh <- min(H, n - o) Ytr <- ts(Y[1:o, ], frequency = frequency(Y)); Yte <- Y[(o + 1):(o + hh), , drop = FALSE] for (j in seq_along(ps)) { if (ps[j] == 0) fc <- matrix(rep(Y[o, ], each = hh), nrow = hh) else { pr <- predict(VAR(Ytr, p = ps[j], type = "const"), n.ahead = hh) fc <- sapply(pr$fcst, function(x) x[1:hh, "fcst"]); if (hh == 1) fc <- matrix(fc, nrow = 1) } res[[j]][ii, 1:hh, ] <- fc - Yte } } nm <- ifelse(ps == 0, "Random walk", paste0("VAR(", ps, ")")) rm <- function(e, h) sqrt(colMeans(e[, h, ]^2, na.rm = TRUE)) tab1 <- data.frame(Model = nm, t(sapply(res, function(e) round(rm(e, 1), 3))), check.names = FALSE) tabH <- data.frame(Model = nm, t(sapply(res, function(e) round(sapply(hs, function(h) mean(rm(e, h))), 3))), check.names = FALSE) names(tabH)[-1] <- paste0("h = ", hs) cat(sprintf("\n\nInitial training set: %d observations. Forecast origins: %d.\n\n", k, length(origins))) print(kable(tab1, row.names = FALSE, caption = "One-Step-Ahead RMSE by Variable")); cat("\n\n") print(kable(tabH, row.names = FALSE, caption = "RMSE Averaged Across Variables, by Horizon")); cat("\n\n") } # Granger-causality table: the joint test from the full VAR and the two # directions of each bivariate VAR at the same lag order. granger_tab <- function(Y, p, cause = "FFR") { full <- VAR(Y, p = p, type = "const") rows <- lapply(setdiff(colnames(Y), cause), function(v) { bv <- VAR(Y[, c(cause, v)], p = p, type = "const") data.frame(Variable = v, `FFR causes variable` = fmt_p(causality(bv, cause = cause)$Granger$p.value), `Variable causes FFR` = fmt_p(causality(bv, cause = v)$Granger$p.value), check.names = FALSE) }) cat(sprintf("Joint test in the full VAR(%d), FFR does not Granger-cause the other variables: F = %.2f, p %s\n\n", p, causality(full, cause = cause)$Granger$statistic, fmt_p(causality(full, cause = cause)$Granger$p.value))) kable(do.call(rbind, rows), row.names = FALSE, caption = sprintf("Bivariate Granger Tests, VAR(%d), p-values", p)) } # Impulse responses to an orthogonalised one-standard-deviation FFR shock. irf_plot <- function(fit, n.ahead, unit, xlab) { ir <- irf(fit, impulse = "FFR", n.ahead = n.ahead, ortho = TRUE, boot = TRUE, runs = 200, seed = 1) long <- do.call(rbind, lapply(colnames(ir$irf$FFR), function(v) data.frame(Variable = factor(v, levels = colnames(ir$irf$FFR)), h = 0:n.ahead, irf = ir$irf$FFR[, v], lo = ir$Lower$FFR[, v], hi = ir$Upper$FFR[, v]))) p <- ggplot(long, aes(h, irf)) + geom_ribbon(aes(ymin = lo, ymax = hi), fill = "grey80") + geom_line() + geom_hline(yintercept = 0, linetype = 2) + facet_wrap(~ Variable, scales = "free_y") + ggtitle("Response to a One-Standard-Deviation Federal Funds Rate Shock") + xlab(xlab) + ylab(unit) print(p) invisible(ir) } # VAR forecast with 80% and 95% bands, panels in variable order, recent history only. var_fc_plot <- function(fit, h, from, title, unit) { Y <- fit$y; tt <- time(Y); vars_ <- colnames(Y) pr80 <- predict(fit, n.ahead = h, ci = 0.80); pr95 <- predict(fit, n.ahead = h, ci = 0.95) step <- 1 / frequency(Y); tf <- max(tt) + step * (1:h) hist <- do.call(rbind, lapply(vars_, function(v) data.frame(Variable = v, t = as.numeric(tt), y = as.numeric(Y[, v])))) fc <- do.call(rbind, lapply(vars_, function(v) data.frame(Variable = v, t = tf, y = pr95$fcst[[v]][, "fcst"], lo95 = pr95$fcst[[v]][, "lower"], hi95 = pr95$fcst[[v]][, "upper"], lo80 = pr80$fcst[[v]][, "lower"], hi80 = pr80$fcst[[v]][, "upper"]))) hist$Variable <- factor(hist$Variable, levels = vars_); fc$Variable <- factor(fc$Variable, levels = vars_) ggplot() + geom_ribbon(data = fc, aes(t, ymin = lo95, ymax = hi95), fill = "grey85") + geom_ribbon(data = fc, aes(t, ymin = lo80, ymax = hi80), fill = "grey65") + geom_line(data = hist[hist$t >= from, ], aes(t, y)) + geom_line(data = fc, aes(t, y), colour = "steelblue4") + facet_wrap(~ Variable, scales = "free_y") + scale_x_continuous(breaks = function(l) seq(ceiling(l[1]), floor(l[2]), by = 2)) + ggtitle(title) + xlab("Year") + ylab(unit) } # Residual correlation matrix and cross-correlation plot for a VAR, as taught: # diagonals are the ACF of each residual series, off-diagonals the CCFs. var_resid_diag <- function(fit, lag.max) { e <- resid(fit) print(kable(round(cor(e), 3), caption = "Correlation Matrix of VAR Residuals")) acf(e, lag.max = lag.max, main = "") } ``` *** ### Purpose This section adds the Federal Funds Rate to the univariate models of the previous section, in two ways. The first two models are regressions with ARMA errors, the ARIMAX form: a target series is regressed on the Federal Funds Rate and the regression residuals are modelled as an ARMA process, so that the coefficient on the policy rate is estimated with the autocorrelation of the target accounted for. The targets are monthly CPI inflation and the monthly change in the 30-year Treasury yield. The last two models are vector autoregressions, in which the Federal Funds Rate is one variable in a system and every variable is explained by the lagged values of all of them. The first VAR covers the Treasury yield curve; the second covers output, prices, unemployment and the mortgage rate. The two approaches differ in what they assume about the policy rate. A regression treats it as exogenous, set without regard to the target. That is defensible for the 30-year yield, which the Federal Reserve does not target, and it is not defensible for inflation, which the Federal Reserve responds to. The inflation regression therefore uses the policy rate lagged twelve months, so that the regressor is fixed before the target is observed, and its coefficient is read as a forecasting relationship rather than a causal one. The VARs make no exogeneity assumption; their Granger-causality tests and impulse responses are the results in this section that speak to the effect of the policy rate. Those two tools extend the course's VAR material, which covers lag selection, estimation, residual diagnostics and forecasting; each is explained where it is first used. Each regression is built in the order the course sets out: the regression is fitted by least squares, the autocorrelation and partial autocorrelation of its residuals are inspected, `auto.arima` is run on those residuals to suggest an error structure, and the regression and the error model are then estimated together with `Arima(..., xreg = )`. A grid of ARMA orders ranked by AIC, AICc and BIC is reported alongside the `auto.arima` choice [@hyndman2008], followed by coefficient tables, Ljung-Box tests and rolling-origin cross-validation. For the ARIMAX models the cross-validation compares each model with the same ARMA model without the regressor, which is the test of whether the policy rate adds forecasting information. The VARs are estimated with the `vars` package [@pfaff2008]. ## ARIMAX: Inflation and the Federal Funds Rate ### Data Inflation is the monthly change in the log of `CPIAUCSL`, multiplied by 100, so a value of 0.30 is a 0.30% rise in the index over the month, about 3.6% at an annual rate. The regressor is the effective Federal Funds Rate, `FEDFUNDS`, lagged twelve months. Both series are monthly; after the lag the sample runs from July 1955 to August 2023, 818 observations. ```{r} ff <- rd("fedfunds.xlsx"); cpi <- rd("cpi-u.xlsx") m <- merge(ff, cpi, by = "Date"); names(m) <- c("Date", "ffr", "cpi") m$infl <- c(NA, 100 * diff(log(m$cpi))); m$dffr <- c(NA, diff(m$ffr)) m$ffr_l12 <- dplyr::lag(m$ffr, 12); m$ffr_l1 <- dplyr::lag(m$ffr, 1) d1 <- m[complete.cases(m), ] infl <- ts(d1$infl, start = c(year(d1$Date[1]), month(d1$Date[1])), frequency = 12) X1 <- as.matrix(d1[, c("ffr_l12", "ffr_l1")]) autoplot(ts(d1[, c("infl", "ffr")], start = start(infl), frequency = 12), facets = TRUE) + ggtitle("Monthly CPI Inflation and the Federal Funds Rate") + xlab("Year") + ylab("% per month (top), % (bottom)") ``` ### Stationarity ```{r} adf_line(infl, "Monthly inflation"); adf_line(d1$ffr, "Federal Funds Rate, level"); adf_line(d1$dffr, "Federal Funds Rate, change") k <- suppressWarnings(kpss.test(infl)); cat(sprintf("%-34s KPSS statistic %6.2f p-value %s\n", "Monthly inflation", k$statistic, ifelse(k$p.value <= 0.01, "< 0.01", sprintf("%.3f", k$p.value)))) ``` Monthly inflation rejects the unit root, with an ADF statistic of −3.73 and p = 0.023, and is modelled in levels. It is persistent: the ACF is 0.61 at lag one and stays near 0.40 through lag eleven, and the KPSS test rejects level stationarity, which is the usual verdict on a series with a slowly shifting mean. The Federal Funds Rate does not reject the unit root in levels on this sample (statistic −3.12, p = 0.105) and is used as a lagged level in the regression, where the question is whether its past value helps predict inflation; the change in the rate is used as an alternative regressor below. The ACF of inflation also shows a spike at lag twelve, 0.30, which is what remains of the seasonal pattern after the Bureau of Labor Statistics' adjustment; seasonal autoregressive terms are included in the grid for that reason. ```{r fig.width=9, fig.height=3} gridExtra::grid.arrange(ggAcf(infl, 36) + ggtitle("ACF, Monthly Inflation"), ggPacf(infl, 36) + ggtitle("PACF, Monthly Inflation"), nrow = 1) ``` ### Regression and Its Errors The model is built in the order the course sets out: fit the regression first, look at what its residuals do, and let that choose the ARMA structure for the errors. The regression of monthly inflation on the Federal Funds Rate lagged twelve months is estimated by ordinary least squares, and its residuals are the series whose autocorrelation has to be modelled. ```{r fig.width=9, fig.height=3} ols1 <- lm(infl ~ ffr_l12, data = d1) kable(data.frame(Coefficient = rownames(summary(ols1)$coefficients), round(summary(ols1)$coefficients, 4), check.names = FALSE), row.names = FALSE, caption = sprintf("Least-Squares Regression of Monthly Inflation on the FFR Lagged 12 Months, R² = %.3f", summary(ols1)$r.squared)) res1 <- ts(residuals(ols1), start = start(infl), frequency = 12) gridExtra::grid.arrange(ggAcf(res1, 36) + ggtitle("ACF, Regression Residuals"), ggPacf(res1, 36) + ggtitle("PACF, Regression Residuals"), nrow = 1) cat("auto.arima on the regression residuals:\n"); auto.arima(res1, seasonal = TRUE) ``` The regression coefficient is 0.029 with an R² of 0.11, but its residuals are far from white noise: the ACF starts at 0.57 and decays slowly, and the PACF has a large first spike. Ordinary least squares is therefore not a valid summary of this relationship on its own, which is the reason for the ARMA error term. Applied to these residuals, `auto.arima` selects an `ARMA(2,4)(2,0,0)[12]` structure; the two seasonal autoregressive terms confirm that the annual pattern survives the Bureau's seasonal adjustment. The next step estimates the regression and an error structure together, searching over orders rather than taking this one as given. ### Model Selection The grid covers ARMA error orders p and q from 0 to 3 and seasonal autoregressive order P from 0 to 2, 48 models, each with the twelve-month lag of the Federal Funds Rate as regressor. ```{r} g1 <- grid_x(infl, X1[, "ffr_l12", drop = FALSE], Pmax = 2) kable(head(g1, 6), row.names = FALSE, caption = "Six Lowest-AICc Error Structures, Regressor FFR Lagged 12 Months") kable(head(g1[order(g1$BIC), ], 3), row.names = FALSE, caption = "Three Lowest-BIC Error Structures") cat("auto.arima, same regressor:\n"); auto.arima(infl, xreg = X1[, "ffr_l12", drop = FALSE]) ``` AICc selects `ARMA(3,2)(2,0,0)[12]` errors; BIC selects `ARMA(2,1)(2,0,0)[12]`, 4.4 points behind on AIC with two fewer parameters. Every model in the top six has two seasonal autoregressive terms. `auto.arima` with the same regressor selects the same model as AICc. The BIC model is retained; the AICc model is shown for comparison. ### Diagnostics ```{r} fit1 <- Arima(infl, order = c(2, 0, 1), seasonal = c(2, 0, 0), xreg = X1[, "ffr_l12", drop = FALSE]) report(fit1) ``` ```{r} report(Arima(infl, order = c(3, 0, 2), seasonal = c(2, 0, 0), xreg = X1[, "ffr_l12", drop = FALSE])) ``` ```{r} report(Arima(infl, order = c(2, 0, 1), seasonal = c(2, 0, 0), xreg = X1[, "ffr_l1", drop = FALSE])) ``` The coefficient on the Federal Funds Rate lagged twelve months is 0.008 with a standard error of 0.009 in the retained model and 0.009 with a standard error of 0.010 in the larger one: the policy rate a year earlier has no measurable relationship with this month's inflation once inflation's own dynamics are in the model. The third table replaces the twelve-month lag with a one-month lag. That coefficient is 0.042 with a standard error of 0.007, positive and significant, and the AIC falls by 28 points. The sign is the opposite of what monetary transmission predicts. It is the reaction function: the Federal Reserve raises the rate when inflation is high, so last month's rate is high when this month's inflation is high. A single-equation regression cannot separate the two directions, which is why the VAR later in this section is the model used for the causal question. The ARMA terms are all significant. In the retained model the Ljung-Box p-values are 0.31, 0.009 and 0.001 at lags 8, 12 and 24; in the AICc model they are 0.13, 0.18 and 0.013. Neither set of residuals is white at the annual lags, though the autocorrelations involved are small. The residual standard deviation is 0.23 points of monthly inflation, against an unconditional standard deviation of 0.31. ### Cross-Validation ```{r} #| results: asis cv_x(infl, X1, list( list(name = "ARMA(2,1)(2,0,0)[12], no regressor", order = c(2, 0, 1), seasonal = c(2, 0, 0), x = NULL), list(name = "Regressor FFR lagged 12 months", order = c(2, 0, 1), seasonal = c(2, 0, 0), x = "ffr_l12"), list(name = "Regressor FFR lagged 1 month", order = c(2, 0, 1), seasonal = c(2, 0, 0), x = "ffr_l1"), list(name = "ARMA(2,1), no seasonal terms", order = c(2, 0, 1), seasonal = c(0, 0, 0), x = NULL))) ``` The initial training set ends in February 2003 and the models are re-estimated at each of the 246 later months, forecasting one to twelve months ahead. The twelve-month lag of the policy rate is known at every horizon, so its forecasts are genuine; the one-month lag is known only one step ahead, and its longer-horizon forecasts use the realised rate. Adding the twelve-month lag changes the one-step RMSE from 0.2877 to 0.2879 and the twelve-step RMSE from 0.335 to 0.336: no gain. The one-month lag, despite its in-sample significance, makes the forecasts worse at every horizon, 0.291 one step ahead and 0.342 at twelve. The seasonal terms are worth about 0.002 of RMSE. For forecasting monthly inflation, the policy rate adds nothing to what inflation's own history provides. ### Forecast ```{r} xf <- matrix(d1$ffr[(nrow(d1) - 11):nrow(d1)], ncol = 1, dimnames = list(NULL, "ffr_l12")) fc1 <- forecast(fit1, h = 12, xreg = xf) autoplot(fc1) + ggtitle("Monthly CPI Inflation, 12-Month Forecast") + xlab("Year") + ylab("% per month") ``` Because the regressor is lagged twelve months, the values it takes over the forecast period are the Federal Funds Rates of September 2022 to August 2023, which are known, so no assumption about future policy is needed. The forecast runs from September 2023 to August 2024. The point forecasts average 0.38% a month, 4.5% at an annual rate, against 3.6% over the twelve months to August 2023; the model carries forward the persistence of the 2021–2023 readings. The 80% interval at twelve months is 0.02 to 0.73% a month, an annual range of roughly 0 to 9%. ## ARIMAX: The 30-Year Yield and the Federal Funds Rate ### Data The target is the monthly change in the 30-year constant-maturity yield, `DGS30`, with the daily series averaged within each month. The regressor is the monthly change in the Federal Funds Rate. The sample runs from March 1977 to August 2023, 557 monthly changes. Both are in % points. ```{r} y30 <- mavg("30-yrs-t-bill-yield.xlsx"); names(y30)[2] <- "y30" m2 <- merge(ff, y30, by = "Date"); names(m2)[2] <- "ffr" m2$dy <- c(NA, diff(m2$y30)); m2$dffr <- c(NA, diff(m2$ffr)); m2$dffr_l1 <- dplyr::lag(m2$dffr, 1) d2 <- m2[complete.cases(m2), ] dy30 <- ts(d2$dy, start = c(year(d2$Date[1]), month(d2$Date[1])), frequency = 12) X2 <- as.matrix(d2[, c("dffr", "dffr_l1")]) autoplot(ts(m2[, c("y30", "ffr")], start = c(year(m2$Date[1]), month(m2$Date[1])), frequency = 12)) + ggtitle("30-Year Treasury Yield and the Federal Funds Rate, Monthly") + xlab("Year") + ylab("%") + guides(colour = guide_legend(title = "Series")) ``` ### Stationarity ```{r} adf_line(m2$y30, "30-year yield, level"); adf_line(dy30, "30-year yield, monthly change"); adf_line(d2$dffr, "Federal Funds Rate, monthly change") cat(sprintf("Correlation of levels: %.3f. Correlation of monthly changes: %.3f. Correlation of the yield change with the previous month's FFR change: %.3f.\n", cor(m2$y30, m2$ffr), cor(d2$dy, d2$dffr), cor(d2$dy, d2$dffr_l1))) ``` The two levels are I(1) and move together over the sample, with a correlation of 0.90 that reflects the common decline from 1981 to 2020 rather than a month-to-month relationship. The monthly changes are stationary, with ADF statistics of −7.38 and −6.34, and their correlation is 0.30. The change in the yield has a lag-one autocorrelation of 0.32 and a PACF that is significant at lags 1 and 2, so an ARMA error structure is needed. ```{r fig.width=9, fig.height=3} gridExtra::grid.arrange(ggAcf(dy30, 24) + ggtitle("ACF, Change in 30-Year Yield"), ggPacf(dy30, 24) + ggtitle("PACF, Change in 30-Year Yield"), nrow = 1) ``` ### Regression and Its Errors ```{r fig.width=9, fig.height=3} ols2 <- lm(dy ~ dffr, data = d2) kable(data.frame(Coefficient = rownames(summary(ols2)$coefficients), round(summary(ols2)$coefficients, 4), check.names = FALSE), row.names = FALSE, caption = sprintf("Least-Squares Regression of the Change in the 30-Year Yield on the Change in the FFR, R² = %.3f", summary(ols2)$r.squared)) res2 <- ts(residuals(ols2), start = start(dy30), frequency = 12) gridExtra::grid.arrange(ggAcf(res2, 24) + ggtitle("ACF, Regression Residuals"), ggPacf(res2, 24) + ggtitle("PACF, Regression Residuals"), nrow = 1) cat("auto.arima on the regression residuals:\n"); auto.arima(res2, seasonal = FALSE) ``` The least-squares coefficient is 0.149 with an R² of 0.09. Its residuals have a first-order autocorrelation of about 0.3 and a PACF that cuts off after two lags, so an ARMA error term is again needed; `auto.arima` on the residuals selects `MA(2)`. Both series are already differenced, which is what the course requires before regressing one non-stationary series on another: fitting the levels against each other would risk a spurious regression. ### Model Selection ```{r} g2 <- grid_x(dy30, X2[, "dffr", drop = FALSE]) kable(head(g2, 6), row.names = FALSE, caption = "Six Lowest-AICc Error Structures, Regressor Change in FFR") cat("auto.arima, same regressor:\n"); auto.arima(dy30, xreg = X2[, "dffr", drop = FALSE]) ``` AIC, AICc and BIC all select `MA(2)` errors, and `auto.arima` agrees. Adding the previous month's change in the Federal Funds Rate as a second regressor raises AICc by 1.9 points and its coefficient is −0.006 with p = 0.78, so it is not retained. ### Diagnostics ```{r} fit2 <- Arima(dy30, order = c(0, 0, 2), xreg = X2[, "dffr", drop = FALSE]) report(fit2) ``` ```{r} report(Arima(dy30, order = c(0, 0, 2))) ``` The coefficient on the change in the Federal Funds Rate is 0.126 with a standard error of 0.022: a one-point rise in the policy rate within a month is associated with a 0.13-point rise in the 30-year yield in the same month, holding the yield's own dynamics fixed. The 30-year yield is priced from expectations of short rates over three decades and a term premium, so a pass-through of about an eighth is the expected order of magnitude. The moving-average terms are 0.33 and −0.13, both significant, and the residuals pass the Ljung-Box test at lags 8, 12 and 24 (p = 0.58, 0.85 and 0.23). Removing the regressor raises the AIC by 30 points. The residual standard deviation is 0.24 points a month. ### Cross-Validation ```{r} #| results: asis cv_x(dy30, X2, list( list(name = "MA(2), no regressor", order = c(0, 0, 2), seasonal = c(0, 0, 0), x = NULL), list(name = "Regressor change in FFR", order = c(0, 0, 2), seasonal = c(0, 0, 0), x = "dffr"), list(name = "Regressors change in FFR and its lag", order = c(0, 0, 2), seasonal = c(0, 0, 0), x = c("dffr", "dffr_l1")), list(name = "No change (random walk in the level)", order = c(0, 0, 0), seasonal = c(0, 0, 0), x = NULL))) ``` The initial training set ends in July 2009 and there are 168 forecast origins. The regressor's forecasts are conditional on the realised change in the policy rate, so they measure how much knowing the Federal Reserve's move would have helped, not what a forecaster could have achieved. On that basis the regressor lowers the one-step RMSE from 0.163 to 0.161 and the twelve-step RMSE from 0.176 to 0.172, a gain of 1 to 2%. The random walk, which forecasts no change, has a one-step RMSE of 0.176. The lagged regressor adds nothing. The improvement from the policy rate is real and small, consistent with the coefficient of 0.13 and a correlation of 0.30 between the two monthly changes. ### Forecast ```{r} xf2 <- matrix(0, nrow = 12, ncol = 1, dimnames = list(NULL, "dffr")) fc2 <- forecast(fit2, h = 12, xreg = xf2) autoplot(fc2) + ggtitle("Monthly Change in the 30-Year Yield, 12-Month Forecast, FFR Unchanged") + xlab("Year") + ylab("% points") ``` A forecast of the yield change requires an assumption about the policy rate. The one shown holds the Federal Funds Rate at its August 2023 level, so the regressor is zero throughout. Under that assumption the forecast change is within 0.01 of zero from the third month, and the 80% interval is about ±0.3 points a month: the model's information about the direction of long yields is exhausted within two months, after which the forecast is the level as of August 2023, 4.28%. ## VAR: The Federal Funds Rate and the Treasury Curve ### Data Seven monthly series from January 1982 to August 2023, 500 months: the effective Federal Funds Rate and Treasury yields at 3 months, 6 months, 1 year, 5 years, 10 years and 30 years. The four short-maturity yields are reconstructed from the spread series described on the data page by adding the Federal Funds Rate back to the monthly average spread; since the monthly `FEDFUNDS` value is itself the average of the daily effective rate, the sum is the monthly average yield up to the difference between trading-day and calendar-day averaging. The 10- and 30-year yields are monthly averages of the daily constant-maturity series. All are in %. ```{r} sp <- c(T3M = "3-mon-t-bill-rate.xlsx", T6M = "6-mon-t-bill-rate.xlsx", T1Y = "1-yrs-t-bill-rate.xlsx", T5Y = "5-yrs-t-bill-rate.xlsx") yl <- c(T10Y = "10-yrs-t-bill-yield.xlsx", T30Y = "30-yrs-t-bill-yield.xlsx") mv <- ff; names(mv)[2] <- "FFR" for (nm in names(sp)) { d <- mavg(sp[[nm]]); names(d)[2] <- nm; mv <- merge(mv, d, by = "Date") } for (nm in names(yl)) { d <- mavg(yl[[nm]]); names(d)[2] <- nm; mv <- merge(mv, d, by = "Date") } for (nm in names(sp)) mv[[nm]] <- mv[[nm]] + mv$FFR Yc <- ts(as.matrix(mv[, -1]), start = c(year(mv$Date[1]), month(mv$Date[1])), frequency = 12) autoplot(Yc) + ggtitle("Federal Funds Rate and Treasury Yields, Monthly") + xlab("Year") + ylab("%") + guides(colour = guide_legend(title = "Series")) ``` ```{r} adf_c <- t(sapply(colnames(Yc), function(v) c(Level = unname(suppressWarnings(adf.test(Yc[, v]))$statistic), `First difference` = unname(suppressWarnings(adf.test(diff(Yc[, v])))$statistic)))) kable(data.frame(Series = colnames(Yc), round(adf_c, 2), `Correlation with FFR` = round(cor(Yc)[, "FFR"], 3), check.names = FALSE), row.names = FALSE, caption = "ADF Statistics and Correlation with the Federal Funds Rate") ``` Every series fails to reject the unit root in levels at the 5% level except the 30-year yield (statistic −3.66), and every first difference rejects it. A VAR is fitted to the levels rather than the differences. This follows the treatment of the mortality and pollution example in the course, where the model is fitted to data that need not be stationary and the deterministic terms are estimated inside the VAR rather than by de-trending each series first. The deterministic term here is a constant only, not the constant-plus-trend used in that example: interest rates have no deterministic trend, and a linear trend term would extrapolate the 1982–2020 decline into the forecast. ### Lag Selection ```{r} vs_c <- VARselect(Yc, lag.max = 12, type = "const") kable(data.frame(Criterion = c("AIC", "HQ", "SC (BIC)", "FPE"), `Selected lag` = as.integer(vs_c$selection), check.names = FALSE), row.names = FALSE, caption = "Lag Order Selected by Each Criterion, Maximum 12") kable(round(vs_c$criteria[1:3, 1:6], 2), caption = "AIC, HQ and SC by Lag Order, Lags 1 to 6") ``` All four criteria select two lags. `VAR(2)` has 7 equations with 15 coefficients each. The largest root of its companion matrix is 0.978, inside the unit circle, and the equations' R² are all above 0.99, which is the persistence of the levels rather than explanatory power. ```{r} var_c <- VAR(Yc, p = 2, type = "const") cat(sprintf("VAR(2): %d equations, %d coefficients each, %d observations. Largest root %.3f.\n", ncol(Yc), ncol(Yc) * 2 + 1, nobs(var_c$varresult$FFR), max(vars::roots(var_c)))) st <- serial.test(var_c, lags.pt = 24, type = "PT.asymptotic") cat(sprintf("Portmanteau test for residual autocorrelation to lag 24: statistic %.1f, p %s.\n", st$serial$statistic, fmt_p(st$serial$p.value))) ``` The portmanteau test rejects the absence of residual autocorrelation. With 49 cross-equation autocorrelations at each lag the test has power against small residual dependence, and the lag-selection criteria did not favour a longer lag; the result is noted and the model is used as selected. ### Residual Diagnostics ```{r fig.width=9, fig.height=7} #| results: asis cat("\n\n"); var_resid_diag(var_c, 24); cat("\n\n") ``` The correlation matrix shows how much of each month's movement the model leaves unexplained in common across equations. The residual correlations between neighbouring maturities are high, 0.95 between the 10- and 30-year yields and 0.94 between the 3- and 6-month, and they fall steadily with the distance between maturities, to 0.25 between the policy rate and the 30-year yield. The model does not capture what the curve does within a month, which is a limitation of a monthly VAR rather than a fault in the lag order, and it is also why the impulse responses below need an ordering assumption. In the plot the diagonal panels are the autocorrelation of each equation's residuals and the off-diagonal panels the cross-correlations between pairs, with the second-named series leading. The diagonals are close to zero beyond lag one; the large spikes at lag zero in the off-diagonals are the contemporaneous correlations just described. ### Granger Causality The Granger test asks a forecasting question, not a causal one: does adding the past of one variable to a regression improve the prediction of another, beyond that variable's own past? If it does, the first variable is said to Granger-cause the second. The test is an F-test on the lags of the candidate variable in one equation of the VAR, so it uses the model already estimated above. It extends the course's VAR material, which stops at lag selection, residual diagnostics and forecasting. The table reports two things. The first line is the joint test from the full VAR: all the other variables at once. The table below it fits a two-variable VAR for each pair at the same lag order and tests each direction separately, which is how the direction of the relationship becomes visible. ```{r} granger_tab(Yc, 2) ``` The joint test rejects: the Federal Funds Rate helps predict the yield curve. The bivariate tests show where. The policy rate Granger-causes the 3-month, 6-month and 1-year yields, with p-values below 0.02, and does not Granger-cause the 5-, 10- or 30-year yields, with p-values of 0.13, 0.33 and 0.46. In the other direction, every yield Granger-causes the Federal Funds Rate, with p-values below 0.003 at all six maturities. The pattern is what the term structure implies. Short yields are averages of expected policy rates over the next months, so the policy rate and the short yields predict each other. Long yields embody expectations of policy years ahead, so they anticipate the policy rate, but last month's policy rate carries no information about the 30-year yield beyond what the 30-year yield already contains. ### Impulse Responses An impulse response traces what the model implies happens to every variable over the following months when one variable receives an unexpected one-period push. It is computed by feeding a one-unit disturbance into the estimated VAR and iterating the equations forward with no further disturbances, so it uses only the coefficients already estimated; the grey bands come from re-estimating the model on 200 bootstrap resamples. This tool is not part of the course's VAR material and is included because it is the only way in this study to put a number on the size and timing of a policy-rate effect. Because the variables move together within a month, a disturbance to one is not independent of the others. The responses below resolve that by ordering the variables and letting each affect only those after it within the month, which is the standard convention; the ordering is stated with each figure and is an assumption, not a result. ```{r fig.width=9, fig.height=6} #| results: asis ir_c <- irf_plot(var_c, 24, "% points", "Months after shock") irt <- data.frame(Series = colnames(Yc), t(sapply(colnames(Yc), function(v) round(ir_c$irf$FFR[c(1, 4, 7, 13, 25), v], 3))), check.names = FALSE) names(irt)[-1] <- paste0("h = ", c(0, 3, 6, 12, 24)) fe_c <- fevd(var_c, n.ahead = 24) irt$`FEVD share at h = 12, %` <- sapply(fe_c, function(x) round(100 * x[12, "FFR"], 1)) cat("\n\n"); print(kable(irt, row.names = FALSE, caption = "Response to a One-Standard-Deviation FFR Shock, % Points, and Share of Forecast-Error Variance Due to FFR Shocks")); cat("\n\n") ``` The shock is orthogonalised with the Federal Funds Rate ordered first, so it is the part of the month's movement in the policy rate not explained by the system's past, and the yields are allowed to respond within the month. Its size is 0.21 points. The bands are 95% bootstrap intervals from 200 replications. The response is largest at the short end and declines with maturity at every horizon. On impact the 3-month yield rises 0.16 points, three-quarters of the shock, the 5-year 0.09 and the 30-year 0.05. The responses build for several months as the policy rate itself continues to rise, since the Federal Reserve moves in sequences: the policy rate peaks at 0.29 points above baseline after three months. At twelve months the 3-month yield is 0.21 points above baseline, 92% of the policy rate's own response, and the 30-year is 0.11, 49%. The 30-year response is still significant at 24 months. Policy-rate shocks account for 54% of the twelve-month forecast-error variance of the policy rate, 44% for the 3-month yield, and 13% for the 30-year, whose variance is dominated by its own shocks. ### Cross-Validation ```{r} #| results: asis cv_var(Yc, c(0, 1, 2, 4), H = 12) ``` The initial training set ends in February 2011 and there are 150 forecast origins. `VAR(2)` beats the random walk one step ahead for the policy rate, 0.115 against 0.163, and for the yields out to one year, 0.156 to 0.159 against 0.166 to 0.171; for the 5-, 10- and 30-year yields it is no better than the random walk, 0.18 to 0.20 for both. Averaged across variables the VAR's advantage over the random walk is 6% one step ahead and 2% at twelve months. `VAR(4)` is within 0.01 of `VAR(2)` at every horizon; `VAR(1)` is slightly worse one step ahead and slightly better from three months on, 1.117 against 1.154 at twelve. The gains are concentrated where the Granger tests placed them, at the short end. ### Forecast ```{r fig.width=9, fig.height=7} var_fc_plot(var_c, 24, 2018, "VAR(2) Forecast, 24 Months, with 80% and 95% Bands", "%") ``` From the August 2023 values, 5.33% for the policy rate and 4.28% for the 30-year, the `VAR(2)` forecasts the policy rate at 5.59% after twelve months and 5.03% after 24, the 3-month yield at 5.25% and 4.73%, and the 30-year at 5.03% and 5.15%. The 95% interval for the policy rate at twelve months runs from 3.2% to 8.0%. The forecast has the curve flattening: the short yields fall and the long yields rise, closing most of the inversion in the August 2023 data. ## VAR: The Federal Funds Rate and the Macroeconomy ### Data Five quarterly series from the second quarter of 1971 to the second quarter of 2023, 209 quarters, the sample set by the start of the mortgage-rate series. Real GDP, `GDPC1`, and the CPI, `CPIAUCSL`, enter as 100 times the natural log, so that changes are in %; the Federal Funds Rate, the unemployment rate and the 30-year mortgage rate are quarterly averages of the monthly or weekly data, in %. Nothing is interpolated: the monthly series are averaged to the quarterly frequency of GDP. The variables are ordered output, prices, unemployment, policy rate, mortgage rate. ```{r} mq <- rd("GDP.xlsx"); names(mq) <- c("Date", "GDP") qfiles <- c(FFR = "fedfunds.xlsx", CPI = "cpi-u.xlsx", UNRATE = "unemploy.xlsx", MORT = "30-yrs-mortgage-rate.xlsx") for (nm in names(qfiles)) { d <- qavg(qfiles[[nm]]); names(d)[2] <- nm; mq <- merge(mq, d, by = "Date") } mq$LGDP <- 100 * log(mq$GDP); mq$LCPI <- 100 * log(mq$CPI) Ym <- ts(as.matrix(mq[, c("LGDP", "LCPI", "UNRATE", "FFR", "MORT")]), start = c(year(mq$Date[1]), quarter(mq$Date[1])), frequency = 4) autoplot(Ym, facets = TRUE) + ggtitle("Macroeconomic Series, Quarterly") + xlab("Year") + ylab("100 × log (GDP, CPI); % (others)") ``` ```{r} adf_m <- t(sapply(colnames(Ym), function(v) c(Level = unname(suppressWarnings(adf.test(Ym[, v]))$statistic), `First difference` = unname(suppressWarnings(adf.test(diff(Ym[, v])))$statistic)))) kable(data.frame(Series = colnames(Ym), Mean = round(colMeans(Ym), 2), `Std. dev.` = round(apply(Ym, 2, sd), 2), round(adf_m, 2), check.names = FALSE), row.names = FALSE, caption = "Summary Statistics and ADF Statistics") ``` Log GDP and the mortgage rate do not reject the unit root in levels; log CPI, unemployment and the policy rate reject at the 5% level on this quarterly sample, log CPI because the ADF regression includes a trend. The system is estimated in levels with a constant, for the reasons given for the yield-curve VAR, and because differencing log GDP and log CPI would discard the level relationships among output, prices and rates that the impulse responses are meant to trace. ### Lag Selection ```{r} vs_m <- VARselect(Ym, lag.max = 8, type = "const") kable(data.frame(Criterion = c("AIC", "HQ", "SC (BIC)", "FPE"), `Selected lag` = as.integer(vs_m$selection), check.names = FALSE), row.names = FALSE, caption = "Lag Order Selected by Each Criterion, Maximum 8") kable(round(vs_m$criteria[1:3, ], 2), caption = "AIC, HQ and SC by Lag Order") var_m <- VAR(Ym, p = 2, type = "const") st_m <- serial.test(var_m, lags.pt = 16, type = "PT.asymptotic") cat(sprintf("\nVAR(2): 5 equations, 11 coefficients each, %d observations. Largest root %.3f. Portmanteau test to lag 16: p = %s.\n", nobs(var_m$varresult$FFR), max(vars::roots(var_m)), fmt_p(st_m$serial$p.value))) ``` All four criteria select two lags. The largest root is 0.995, close to but inside the unit circle, and the portmanteau test does not reject residual whiteness (p = 0.32). The Federal Funds Rate equation has significant coefficients on its own two lags, 1.17 and −0.36, and on the second lag of the mortgage rate; the mortgage-rate equation has significant coefficients on the policy rate at both lags, 0.17 and −0.12, and on the price level at both lags. ### Residual Diagnostics ```{r fig.width=9, fig.height=6} #| results: asis cat("\n\n"); var_resid_diag(var_m, 16); cat("\n\n") ``` The residual correlations here are smaller than in the yield-curve system, since quarterly macroeconomic series move less in lockstep within a period than interest rates do within a month. The largest is −0.80 between output and unemployment, the within-quarter form of the relationship between growth and joblessness, followed by 0.53 between the policy rate and the mortgage rate. ### Granger Causality The same test is applied to the macroeconomic system: whether the past of the policy rate improves the forecast of each other variable, and whether each improves the forecast of the policy rate. ```{r} granger_tab(Ym, 2) ``` The policy rate Granger-causes every other variable in the system, with bivariate p-values below 0.001 for output and the mortgage rate, 0.002 for prices and 0.004 for unemployment. In the other direction, output Granger-causes the policy rate (p = 0.023) and prices do so at the 10% level (p = 0.052); unemployment and the mortgage rate do not. This is the reaction function seen from the other side: the Federal Reserve's rate responds to lagged output and inflation, and the rate in turn predicts all four. ### Impulse Responses The same exercise for the macroeconomic system. The variance decomposition in the last column of the table reports the share of each variable's forecast error at that horizon attributable to policy-rate disturbances. ```{r fig.width=9, fig.height=6} #| results: asis ir_m <- irf_plot(var_m, 20, "% (GDP, CPI: % of level)", "Quarters after shock") irtm <- data.frame(Series = colnames(Ym), t(sapply(colnames(Ym), function(v) round(ir_m$irf$FFR[c(1, 3, 5, 9, 13, 21), v], 3))), check.names = FALSE) names(irtm)[-1] <- paste0("h = ", c(0, 2, 4, 8, 12, 20)) fe_m <- fevd(var_m, n.ahead = 20) irtm$`FEVD share at h = 8, %` <- sapply(fe_m, function(x) round(100 * x[8, "FFR"], 1)) cat("\n\n"); print(kable(irtm, row.names = FALSE, caption = "Response to a One-Standard-Deviation FFR Shock and Share of Forecast-Error Variance Due to FFR Shocks")); cat("\n\n") ``` The shock is orthogonalised with the policy rate ordered after output, prices and unemployment and before the mortgage rate: the Federal Reserve is allowed to respond to the quarter's macroeconomic data, the mortgage rate is allowed to respond to the quarter's policy rate, and output, prices and unemployment respond only from the next quarter. The shock is 0.88 points, and the policy rate stays 0.57 points above baseline after a year and 0.18 after three years. Three of the four responses have the sign monetary theory predicts. Real GDP is unchanged for two quarters, then falls, to 0.49% below baseline after two years and 0.53% after three, with the 95% band excluding zero from the third quarter. Unemployment rises by 0.16 points at two to three years, with the band excluding zero from the sixth quarter. The mortgage rate rises 0.19 points on impact and 0.37 at two quarters, about 40% of the policy move, then decays with it. The price level is the exception: it rises, by 0.63% after two years, with a band that excludes zero. A positive price response to a tightening is a known result in this literature, called the price puzzle, and the usual explanation is that the policy rate rises when the Federal Reserve expects inflation that a five-variable system does not otherwise see, so the estimated shock is partly anticipation of inflation rather than policy alone. A measure of expected inflation is the standard remedy; this data set does not contain one. Policy-rate shocks account for 63% of the two-year forecast-error variance of the policy rate, 42% for the mortgage rate, 17% for prices, 12% for output and 3% for unemployment. ### Cross-Validation ```{r} #| results: asis cv_var(Ym, c(0, 1, 2, 3, 4), H = 8, hs = c(1, 2, 4, 8)) ``` The initial training set ends in the third quarter of 2007 and there are 63 forecast origins, so the test period is the financial crisis, the zero-rate years, the pandemic and the 2022 tightening. On this period `VAR(1)` forecasts better than `VAR(2)` on average at every horizon and for every variable but the price level one step ahead, and the random walk forecasts the policy rate better than any VAR: one step ahead its RMSE is 0.46 points against 0.50 for `VAR(1)` and 1.02 for `VAR(2)`. The VARs' advantage is in output and prices at longer horizons, where `VAR(1)` has a two-year RMSE of 3.0% for log GDP against 4.7% for the random walk, and 3.4% against 5.4% for log CPI. The pattern says that on a test period dominated by a rate held at zero for seven years and a one-quarter fall in GDP of 8%, the lag-2 dynamics estimated on 1971–2007 do not carry over, and a simpler model is safer. `VAR(2)` is retained as the model the selection criteria chose on the full sample, and the cross-validation is the caveat on its forecasts. ### Forecast ```{r fig.width=9, fig.height=6} var_fc_plot(var_m, 8, 2015, "VAR(2) Forecast, 8 Quarters, with 80% and 95% Bands", "100 × log (GDP, CPI); % (others)") ``` From the second quarter of 2023 the `VAR(2)` forecasts the policy rate falling from 4.99% to 4.06% after four quarters and 2.91% after eight, unemployment rising from 3.6% to 4.8% and 5.6%, the mortgage rate easing from 6.51% to 5.76%, real GDP flat over the two years, and the price level rising 2.5% a year. It is the profile of a slowdown in which the Federal Reserve cuts rates as unemployment rises. The cross-validation above shows this model's policy-rate forecasts have been poor over the last fifteen years, and the intervals are wide: the 95% band for the policy rate at eight quarters runs from −1.7% to 7.5%. ## Summary ```{r} summary_tab <- data.frame( Section = c("ARIMAX, inflation", "ARIMAX, 30-year yield", "VAR, Treasury curve", "VAR, macroeconomy"), Frequency = c("Monthly", "Monthly", "Monthly", "Quarterly"), Sample = c("1955–2023", "1977–2023", "1982–2023", "1971–2023"), `Retained model` = c("ARMA(2,1)(2,0,0)[12] errors, FFR lagged 12 months", "MA(2) errors, change in FFR", "VAR(2), levels, constant", "VAR(2), levels, constant"), `Effect of the FFR` = c("Coefficient 0.008 (s.e. 0.009); no forecast gain", "Coefficient 0.126 (s.e. 0.022); RMSE 1–2% lower", "Granger-causes yields to 1 year; not 5 to 30 years", "Granger-causes all four; GDP −0.5%, unemployment +0.16 at 2–3 years"), check.names = FALSE) kable(summary_tab, row.names = FALSE, caption = "Retained Models and Main Findings") ``` ### What the Models Found The four models agree on where the policy rate's influence is direct and where it is not. In the Treasury market, the effect runs through the short end. The policy rate Granger-causes the 3-month, 6-month and 1-year yields and not the longer ones; the impulse responses show the 3-month yield moving three-quarters of the way with the policy rate on impact and nine-tenths after a year, while the 30-year moves a quarter on impact and a half after a year; and the single-equation regression puts the same-month pass-through to the 30-year yield at 0.13. Every yield Granger-causes the policy rate, which is the market anticipating the Federal Reserve. In the macroeconomy the policy rate Granger-causes output, prices, unemployment and the mortgage rate; a 0.9-point tightening shock lowers real GDP by half a percent over two to three years and raises unemployment by 0.16 points, with the mortgage rate absorbing about 40% of the move within two quarters. The inflation regression is the model that did not find an effect, and the reason is instructive. The policy rate lagged a year has a coefficient indistinguishable from zero and adds nothing to inflation forecasts; the policy rate lagged a month has a positive, significant coefficient and makes the forecasts worse. Both results are what a reaction function produces. The Federal Reserve raises the rate when inflation is high, so the rate is correlated with inflation in the wrong direction for a causal reading and in a direction that does not survive out of sample. The VAR, which lets both variables respond to each other, finds that the policy rate Granger-causes prices, but its impulse response has the price level rising after a tightening, the price puzzle. The macroeconomic effect of policy on inflation is the one question in this section the data set cannot answer cleanly, because the system lacks the variable that would let it separate the Federal Reserve's information from its action. ### Model Choice The regression form suits the 30-year yield and does not suit inflation. The Federal Reserve does not set the policy rate in response to the 30-year yield, so treating the policy rate as exogenous to it is a reasonable approximation at monthly frequency, and the estimated pass-through is the kind of quantity a regression can deliver. The Federal Reserve does set the policy rate in response to inflation, and a regression of inflation on the policy rate estimates that response as much as it estimates transmission. The twelve-month lag makes the regressor predetermined but cannot make it exogenous, since the rate a year ago responded to inflation a year ago, which is correlated with inflation today. For questions of effect, the VAR is the right instrument in this data set. Two of the tools applied to the VARs, the Granger-causality test and the impulse response function, go beyond the vector-autoregression material covered in the course, which ends at lag selection, residual diagnostics and forecasting; both are explained where they appear, and the impulse responses rest on an ordering assumption that is stated with them. The two VARs are estimated in levels with a constant, following the course's treatment of a VAR on non-stationary data. Both lag selections settled on two lags, and in both the residuals carry large contemporaneous cross-correlations and either fail the Portmanteau test or lose to a simpler model out of sample, which is the usual position of a small VAR on fifty years of data that include several policy regimes. The macroeconomic VAR would be improved by a measure of expected inflation, which is the variable most often added when the price response comes out positive. ### Forecasting None of the four models forecasts well far ahead, and none should be expected to. The inflation model forecasts 4.5% for the year from September 2023, carrying forward the 2021–2023 persistence, against 3.6% in the year before. The 30-year model has no information beyond two months once the policy rate is held fixed. The yield-curve VAR beats a random walk by 6% one step ahead and 2% at a year, entirely at the short end. The macroeconomic VAR's policy-rate forecasts were worse than a random walk over 2007–2023, and its eight-quarter projection, a policy rate of 2.9%, unemployment at 5.6% and flat real GDP, should be read as the pattern of a typical easing cycle in the sample rather than a prediction. ### Limitations The findings rest on Granger causality and orthogonalised impulse responses, both of which identify effects from the ordering of variables in time and within the period, not from any measure of the Federal Reserve's intentions. The yield-curve VAR uses the effective rate, which moved with the target during the 1982–2008 period and sat at the floor of the target range from 2009 to 2015 and 2020 to 2022; a VAR on that sample treats a rate the Federal Reserve could not lower as an ordinary observation. The short-maturity yields are reconstructed from spreads and differ from the published constant-maturity yields by the averaging convention. The macroeconomic VAR has five variables and 209 quarters; with two lags each equation has 11 coefficients, which is adequate, but the bootstrap bands are wide at long horizons, and the cross-validation on the post-2007 period shows the estimated dynamics do not hold at the zero bound.