Granger causality asks whether the past of one variable improves predictions of another variable, beyond the information in that variable’s own history. However, a useful forecasting relationship does not automatically establish a causal mechanism.
This guide develops that distinction through 14 exercises covering housing starts, consumption, industrial production, stock indexes, and temperature. Along the way, you will work through autoregressive models, moving averages, unit root tests, cointegration, lag selection, and joint significance tests.
Source and calculation note: Exercise numbers 11.20–11.33 refer to the supplied textbook photographs, including visible pages 261–265. Explanations here are original. “Printed result” means the result shown in those pages; “recalculated” means a calculation from the visible table or stated summary statistics. Historical data in the photographs are not current economic estimates. Full book title and edition were not available in the supplied material.
Granger causality: meaning, notation, and model choice
Suppose you forecast household consumption next month. Last month’s consumption provides a natural starting point. Next, you add earlier housing starts and ask whether those extra predictors help explain consumption after accounting for its own lags.
In a linear model, the Granger causality test examines a group of lag coefficients together. Therefore, the relevant comparison involves two models with the same dependent variable, estimation dates, and retained controls.
Here, C denotes consumption and H denotes housing starts in an appropriate stationary form. Under the null hypothesis, every β coefficient equals zero. By contrast, the alternative requires at least one nonzero coefficient; it does not require every lag to matter separately.
Implementation reference: the statsmodels Granger causality documentation describes the joint null, F-based tests, and input direction.
| Symbol | Meaning | Practical implication |
|---|---|---|
| N | Original number of observations | Count observations before transformations. |
| n | Usable regression observations | Differencing and lags reduce this number. |
| p | Number of lags of each variable | A bivariate equation with a constant uses 1 + 2p coefficients. |
| SSER | Restricted residual sum of squares | Measures fit without the tested predictor lags. |
| SSEU | Unrestricted residual sum of squares | Measures fit with those lags included. |
| kU | Unrestricted regression coefficients | Includes the intercept and any extra deterministic terms. |
Choose the data form before testing prediction
Stationarity concerns the stability of a process’s statistical properties, rather than whether its plotted values look flat. For example, a stationary series can fluctuate substantially. Conversely, two smooth upward trends can produce a misleading regression.
| Time-series finding | Typical starting model | Important qualification |
|---|---|---|
| Both variables are I(0) | VAR in levels | Check remaining serial correlation and deterministic terms. |
| Both variables are I(1), without cointegration | VAR in first differences | Interpret predictive relationships between changes. |
| Both variables are I(1), with cointegration | Vector error-correction model | Retain the long-run equilibrium term. |
| Stationarity around a deterministic trend | Model the trend explicitly | Differencing is not an automatic requirement. |
| Uncertain integration or mixed orders | Investigate further before using a standard test | The simple decision table does not cover every system. |
For a first difference, ΔXt = Xt − Xt−1. However, an ordinary difference is not a percentage change. Likewise, a change in a stock index is not automatically a stock return.
Background: Hyndman and Athanasopoulos, Forecasting: Principles and Practice, ARIMA models.
Exercise 11.20: how Granger causality differs from causal effects
A variable Granger-causes another variable when its past adds predictive information, relative to the chosen information set. Nevertheless, an intervention asks a different question: what would happen if someone changed the first variable while holding the relevant causal structure in view?
Consider mortgage demand and furniture purchases. A rise in expected household income could influence both series at different times. Consequently, mortgage demand might predict furniture spending even when a mortgage-demand intervention would not reproduce the same relationship.
Another possibility involves anticipation. For example, firms might change orders because they expect future consumer spending to rise. In that setting, earlier orders can forecast later spending even though expectations about spending drive the ordering decision.
Also, significance does not reveal the duration of a response. A predictive relationship may be temporary or persistent. Therefore, separate the existence of predictive content from its direction, size, stability, and economic mechanism.
Overstatement: “Increasing H will necessarily increase C.”
Exercise 11.21: should housing and consumption enter in levels or differences?
The photographed housing example covers January 1997 through December 1999, giving 36 monthly observations. It labels housing starts Y and real personal consumption X. However, X is the dependent variable in the later predictive regression, so keep the direction clear.
The historical consumption values use billions of 1996 dollars. By comparison, the current FRED PCEC96 page reports billions of chained 2017 dollars. Therefore, downloading today’s series will not necessarily reproduce an old textbook table.
Official metadata: FRED real personal consumption expenditures and FRED housing starts. Both current series pages identify monthly observations expressed at seasonally adjusted annual rates; frequency and annualization describe different features.
Step 1: interpret the unit root null correctly
The unit root null sets ρ = 0 in this parameterization. A sufficiently negative test statistic favors stationarity under the specified deterministic terms. However, failing to reject does not prove that a unit root exists.
| Series tested | Printed statistic | Printed comparison value | Decision using that comparison |
|---|---|---|---|
| Housing level Y | −2.48 | −3.33 | Do not reject the unit root null. |
| Consumption level X | 1.72 | −3.33 | Do not reject the unit root null. |
| Housing difference ΔY | −8.92 | −3.33 | Reject the unit root null. |
| Consumption difference ΔX | −5.84 | −3.33 | Reject the unit root null. |
The last row exposes a wording error in the photographed explanation. Since −5.84 is below −3.33, that comparison supports rejection, despite the printed sentence saying otherwise. Thus, the stated numerical evidence supports treating both variables as I(1) for the exercise.
In a fresh analysis, choose the constant, trend, and augmentation lags deliberately. Moreover, use critical values appropriate to that specification and sample size. A single cutoff cannot serve every Dickey–Fuller regression.
Specification reference: statsmodels ADF documentation distinguishes constant-only, constant-plus-trend, and other deterministic specifications.
Step 2: investigate cointegration before discarding levels
Next, estimate a candidate long-run relationship and test whether its residual behaves as a stationary series. The page reports a residual-test statistic around −2.46 and concludes that it lacks evidence of cointegration. Consequently, the exercise proceeds with first differences and no error-correction term.
However, a residual cointegration test needs its own reference distribution. Do not run an ordinary ADF routine on estimated residuals and mechanically treat the usual ADF p-value as an Engle–Granger p-value. In particular, the first-stage fitted relationship changes the inference.
Method reference: statsmodels augmented Engle–Granger test, which uses the null of no cointegration and dedicated approximate p-values and critical values.
Step 3: retain the correct units
For this example, ΔX measures changes in a real consumption series, while ΔY measures changes in housing starts. Therefore, a coefficient linking them has a specific unit ratio. It is not an elasticity unless the transformation and model justify that interpretation.
Exercise 11.22: calculate AIC for one through six lags
The supplied exercise compares equations with one through six lags of both differenced variables. With N = 36 original monthly observations, first differencing leaves 35 observations. Then p lags leave n = 35 − p usable rows.
This is the exercise’s normalized regression criterion. Standard likelihood-based software may include additional constants and count the estimated innovation variance as a parameter. Therefore, document the convention before comparing printed numbers across packages.
| Lags p | Usable n | Coefficients k | Printed SSE | Recalculated criterion |
|---|---|---|---|---|
| 1 | 34 | 3 | 9,297.932 | 5.787657 |
| 2 | 33 | 5 | 7,797.001 | 5.768017 |
| 3 | 32 | 7 | 7,354.929 | 5.874890 |
| 4 | 31 | 9 | 4,617.587 | 5.584286 |
| 5 | 30 | 11 | 3,742.738 | 5.559709 |
| 6 | 29 | 13 | 3,085.670 | 5.563780 |
Using the printed inputs, five lags gives the lowest criterion. However, the difference between five and six lags is only 0.004071. Rounding both values too early can hide this distinction.
The deeper issue: the estimation window changes
Each row uses a different number of observations. Consequently, the table changes both model complexity and the dates contributing to fit. For a clean conditional-regression comparison, fit every candidate over the same dependent-variable dates.
Here, six lags set the latest common starting point. Therefore, each candidate can use 29 observations, even when it includes fewer than six lags. Later in this guide, a reproducible comparison shows what that adjustment does to the ranking.
Why the sample matters: Rob J Hyndman, Facts and fallacies of the AIC, explains why comparable likelihoods must use the same data and conditioning.
Exercise 11.23: test whether housing starts predict consumption
The printed solution chooses five lags and reports SSER = 5,648.53 and SSEU = 3,742.74. Since the unrestricted equation contains 11 coefficients and uses 30 observations, its residual degrees of freedom equal 19. Meanwhile, the null imposes five restrictions on the housing-lag coefficients.
The 5% critical value for F(5, 19) is approximately 2.7401. Also, the upper-tail p-value is approximately 0.1356. Therefore, the printed-summary calculation does not reject the null at 5%.
A careful conclusion stays local: these five housing lags do not show sufficient incremental predictive content in this fitted consumption equation. However, the result does not prove that housing has no economic connection to spending. Nor does it establish the reverse predictive direction.
Why individual t-tests cannot replace this F-test
Several lagged predictors may overlap strongly. As a result, one coefficient can look individually weak while the group contributes useful information. Conversely, selecting the most favorable t-statistic from many lags changes the question and can encourage selective reporting.
The test’s conventional reference distribution also rests on modeling assumptions. In particular, remaining serial correlation, variance instability, or a poorly specified dynamic equation can undermine the nominal p-value. Therefore, inspect residuals before making a substantive claim.
Exercise 11.24: compute an AR(2) process step by step
An autoregressive model builds today’s value from earlier values of the same series. For this exercise, the coefficients equal 0.4 and −0.3. To reproduce the printed answers, use zero presample values and no intercept.
The first supplied disturbances are 0.1291, 0.6910, 0.1348, and 0.1510. Therefore, the recursion starts with y1 = 0.1291. Next, y2 = 0.4(0.1291) + 0.6910 = 0.74264.
| Requested observation | Recalculated value | Rounded answer |
|---|---|---|
| y4 | 0.08545840 | 0.0855 |
| y8 | 0.36178363 | 0.3618 |
| y20 | 0.76253177 | 0.7625 |
| y30 | 0.51881726 | 0.5188 |
Unlike an MA recursion, this process carries earlier disturbances forward through previous y values. Consequently, you must compute the intermediate observations before reaching y20 or y30. Keep full numerical precision during that recursion and round only the displayed answers.
A small modeling detail with a large effect
The table lists positive numerical inputs for a calculation exercise. However, their positive appearance alone does not establish a mean-zero innovation process. A theoretical AR model usually specifies the innovation distribution and handles the mean through centering or an intercept.
Exercise 11.25: compute an MA(2) process and check the sign convention
By comparison, a moving-average process combines current and earlier disturbances directly. Here, the photographed answers correspond to a minus-sign coefficient convention with θ1 = 0.2 and θ2 = −0.5. Consequently, the second lag contributes positively.
| Observation | Recalculated value | Rounded answer |
|---|---|---|
| y4 | 0.46954 | 0.4695 |
| y8 | 0.60900 | 0.6090 |
| y20 | 1.02536 | 1.0254 |
| y30 | 0.47374 | 0.4737 |
Unlike the AR(2) case, computing y30 only requires disturbances at dates 30, 29, and 28. Nevertheless, software may write MA coefficients using a plus-sign convention. Always compare the actual equation before transferring coefficient signs between a textbook and a package.
Also, an MA(2) stochastic process differs from a two-period moving average used to smooth observed data. The first describes how innovations generate a series. By contrast, the second applies an averaging operation to existing observations.
Exercise 11.26: calculate six temperature autocorrelations
The supplied temperature table contains 31 annual observations for New York’s Central Park, from 1969 through 1999. Autocorrelation compares a series with a lagged version of itself. Therefore, lag one compares adjacent years, while lag six compares years six periods apart.
| Lag | Printed ACF | Recalculated ACF using the displayed formula |
|---|---|---|
| 1 | -0.0051 | -0.0052 |
| 2 | -0.0013 | -0.0046 |
| 3 | -0.2007 | -0.1855 |
| 4 | 0.2448 | 0.2112 |
| 5 | -0.1598 | -0.1370 |
| 6 | 0.1023 | 0.0792 |
The table deliberately separates the answer key from the transparent recalculation. Although the general pattern looks similar, the figures do not match exactly. Differences in normalization can affect sample autocorrelations, but this calculation does not attribute every discrepancy to that explanation without further evidence.
As a rough visual reference, ±1.96/√31 gives bands near ±0.352. None of the six displayed correlations exceeds those bands. However, such pointwise bands do not replace a joint test across all six lags.
Why this is not a climate-trend verdict
A short local temperature example can illustrate statistical methods. Nevertheless, failure to detect autocorrelation does not establish the absence of a deterministic trend, structural change, or a broader climate signal. Those questions require appropriate datasets, specifications, and scientific context.
Exercise 11.27: distinguish Box–Pierce Q from Ljung–Box Q
The photographed answer reports Q ≈ 4.22. Using its six printed autocorrelations, the Box–Pierce formula reproduces that number. However, the Ljung–Box formula applies an additional finite-sample adjustment and gives a different result.
| Inputs and method | Q statistic | Interpretation at 5% |
|---|---|---|
| Printed ACF; Box–Pierce | 4.2233 | Do not reject zero autocorrelation through lag six. |
| Printed ACF; Ljung–Box | 5.1762 | Do not reject under the same six-degree reference. |
| Recalculated ACF; Ljung–Box | 3.9437 | p ≈ 0.6843; do not reject. |
For a raw-series illustration with six lags and no fitted ARMA terms, the usual reference uses six degrees of freedom. Its 5% chi-square cutoff is approximately 12.592. Therefore, all three calculations lead to the same broad decision, despite their different numerical values.
After fitting an ARMA model, account for the degrees of freedom that model estimation consumes. Moreover, use the same autocorrelation definition consistently throughout a replication. Mixing one formula’s input conventions with another formula’s label makes results difficult to audit.
Reference: statsmodels Ljung–Box documentation, including its Box–Pierce option and model degrees-of-freedom adjustment.
Exercise 11.28: test the temperature series without a trend
This exercise uses a constant but no deterministic trend and no lagged differences. Therefore, the regression reduces to the basic Dickey–Fuller form. The zero-lag specification does not add extra terms to absorb serial correlation.
Recalculating the coefficient t-statistic from the visible temperature table gives approximately −5.0908. This matches the printed −5.09 after rounding. Consequently, the exercise’s rejection of the unit root null is consistent with the reported statistic under its constant-only specification.
However, the coefficient t-ratio does not follow an ordinary Student t distribution under the unit root null. Use the relevant Dickey–Fuller reference values. Also, distinguish rejection of a unit root from proof that every possible model feature is stable.
Exercise 11.29: add a trend and interpret the reported F-statistic
Adding a deterministic trend changes the test regression and its reference distribution. For example, the equation can include both an intercept and a time index. Therefore, the constant-only exercise and the trend exercise do not test identical alternatives.
The answer key reports an F-statistic of 15.85 and a rejection of the unit root hypothesis. However, an F-statistic in a Dickey–Fuller setting requires a clearly stated joint restriction and nonstandard reference values. The brief question does not fully identify those restrictions.
For transparency, an OLS regression of the displayed trend equation gives a t-ratio on the lagged temperature level of approximately −5.6109. This is a separate statistic, not a numerical reproduction of the printed F = 15.85. Accordingly, retain the printed F as an answer-key result rather than pretending that the two measures are interchangeable.
Exercise 11.30: examine the DJIA and S&P 500 for cointegration
The stock-index table covers 20 trading observations from January 2 through January 30, 2001. It labels the DJIA Y, the S&P 500 X, and the Toronto index Z. Unlike the monthly macroeconomic examples, these observations follow a trading-day calendar.
The printed answers describe both Y and X as possessing unit roots and report no cointegration between them. In statistical reporting, translate the first statement into a test decision: the specified tests did not reject the unit root null. Likewise, the second means that the procedure did not establish cointegration at its chosen significance level.
Why this result needs unusually careful handling
Only 20 observations support the exercise. Consequently, lag selection and deterministic terms can consume a substantial share of the available information. Short samples also make it difficult to distinguish persistent stationary behavior from a unit root.
Moreover, the photographed S&P row appears to contain an unusually large value of 2360.4 on January 24, surrounded by observations near 1300. This is a data-quality flag, not permission to replace it with a guessed number. Therefore, verify the historical source before claiming an independently replicated cointegration result.
This guide reports the supplied answer-key conclusion without presenting an unverified recalculation as settled evidence. A reliable extension would obtain the original data vintage, document any correction, and rerun the tests with explicit settings.
Exercise 11.31: examine the Toronto index and the S&P 500
The answer key reports a unit root for Z and cointegration between X and Z. Under that exercise’s assumptions, a stationary linear combination links two individually nonstationary series. Consequently, modeling only their first differences would omit information about deviations from that relationship.
The error-correction term tracks the prior deviation from the estimated equilibrium. Meanwhile, lagged differences describe shorter-run movements. Its economic interpretation depends on the normalization, adjustment coefficient, and system specification.
For example, a negative adjustment coefficient can describe movement back toward equilibrium under the displayed normalization. However, no numerical adjustment speed follows from a simple “cointegrated” answer. Estimate that coefficient before discussing how quickly the relationship corrects.
The same small-sample and data-quality limitations apply here. Therefore, treat the printed X–Z conclusion as a pedagogical result until a clean replication resolves the suspicious input and model settings.
Historical context: the Royal Swedish Academy of Sciences’ 2003 announcement recognizes Granger’s work on common trends and cointegration. That recognition concerns the treatment of nonstationary economic series, not a guarantee that every predictive test identifies an intervention effect.
Exercise 11.32: one-lag prediction of industrial production
The second causality example uses 35 monthly observations from February 1998 through December 2000. Importantly, the table already contains first differences. Therefore, do not difference those numbers again when reproducing the exercise.
Let y denote the table’s industrial-production change and x denote its S&P 500 change. With one lag and an intercept, the unrestricted equation uses 34 observations and three coefficients. Recalculation gives the following fitted relationship:
The test of the x-lag coefficient gives F(1, 31) ≈ 4.1802 and p ≈ 0.04947. Therefore, the conventional test just rejects at 5%. However, the result sits close to the threshold and deserves a more cautious interpretation than a strong binary claim.
Interpret the short-run coefficient in its actual units
The coefficient −0.016609 measures the next-period predictive association with a one-point increase in the prior stock-index change, conditional on the lagged production change. For example, a 100-point larger prior stock-index change corresponds to roughly −1.661 units of predicted production change in this fitted equation. This is not a percentage elasticity or an identified policy effect.
Calculate the dynamic multiplier
The printed answer rounds this quantity to approximately −0.0118. Yet the interpretation depends on the variables’ transformation. For a permanent one-unit increase in x, the formula describes the steady-state change in y under the stable fitted recursion.
Alternatively, a one-period unit pulse in x generates a sequence of y responses whose infinite sum equals the same expression. Because these variables already represent changes, neither interpretation automatically becomes a cointegrating long-run relationship between the original index levels. Thus, “long-run causality” would overstate what this simple calculation establishes.
Exercise 11.33: what changes when the model uses six lags?
Six lags leave 29 usable observations from the 35-row differenced table. The unrestricted model estimates 13 coefficients, leaving only 16 residual degrees of freedom. Meanwhile, the null imposes six restrictions on the stock-index lags.
At the 5% level, this specification does not reject the null. However, it would reject at 10%. The printed F ≈ 2.60 therefore agrees with the recalculation, while the “No” answer depends on the stated 5% decision rule.
| Specification | Usable observations | Residual degrees of freedom | F statistic | p-value | 5% decision |
|---|---|---|---|---|---|
| One lag | 34 | 31 | 4.1802 | 0.04947 | Reject |
| Six lags | 29 | 16 | 2.5960 | 0.05957 | Do not reject |
These models test different sets of coefficients with different samples and precision. Consequently, their conclusions need not match. More lags can capture richer dynamics, but they also use more parameters and reduce available observations.
Why a lower AIC does not guarantee a smaller causality p-value
Information criteria assess a fit-and-complexity tradeoff for the whole equation. By contrast, the Granger test assesses a specific restriction on a group of coefficients. Therefore, selecting a model by AIC and obtaining a nonsignificant Granger test are entirely compatible outcomes.
To choose among lag orders, first use a comparable estimation window. Next, examine an information criterion alongside residual diagnostics and the application’s timing. Finally, report sensitivity across a reasonable prespecified range rather than selecting whichever lag produces the desired p-value.
Original audit: what a standard answer key leaves out
A worked example becomes more valuable when readers can see where the inputs, calculations, and interpretation separate. Accordingly, this section distinguishes arithmetic replication from estimation using the visible raw data. That distinction prevents a polished table from implying more certainty than the evidence supports.
Housing: printed summaries versus visible data
Re-estimating the housing regressions from the 36 visible observations does not exactly reproduce the printed residual sums of squares. For example, the five-lag regression produces a normalized criterion near 5.695973, rather than 5.559709 from the printed SSE. Therefore, the printed solution and the visible data do not form a fully consistent replication package under the stated OLS setup.
The cause of that discrepancy remains unresolved. Possible explanations include a different underlying data vintage, an unreported input change, or a numerical error in the source. However, this guide does not choose one explanation without evidence.
| Lags | Criterion from printed SSE | Refit: changing sample | Refit: common 29-row sample |
|---|---|---|---|
| 1 | 5.787657 | 5.790469 | 5.596197 |
| 2 | 5.768017 | 5.760763 | 5.511987 |
| 3 | 5.874890 | 5.875067 | 5.611986 |
| 4 | 5.584286 | 5.638391 | 5.702113 |
| 5 | 5.559709 | 5.695973 | 5.748940 |
| 6 | 5.563780 | 5.578743 | 5.578743 |
Using the visible data with changing samples, six lags minimizes the displayed criterion. By comparison, fixing the sample at 29 observations selects two lags among these six candidates. This is an original sensitivity calculation, not a claim that two lags universally provides the best economic model.
How much predictive fit do the extra lags add?
Using the printed housing sums of squares, the extra housing lags reduce residual variation by about 33.7% relative to the restricted equation. Nevertheless, the joint p-value remains about 0.136 because the model tests five restrictions with limited residual degrees of freedom. Thus, a sizeable in-sample fit change does not automatically imply strong statistical evidence.
For the industrial-production example, the corresponding partial R² is approximately 0.1188 with one lag and 0.4933 with six lags. Yet the six-lag model has the larger p-value. This comparison illustrates why effect size, model flexibility, and statistical uncertainty deserve separate attention.
The main audit findings in one place
| Finding | Why it matters | Better practice |
|---|---|---|
| The printed ΔX unit root sentence contradicts its numerical comparison | It reverses the stated decision. | Follow the stated null, statistic, and correct critical values. |
| Five- and six-lag printed AIC values differ by only 0.004071 | Early rounding can conceal the ranking. | Preserve precision and report sensitivity. |
| Changing lag order changes the sample in the original table | The comparison combines fit changes with date changes. | Hold the estimation window fixed. |
| Raw housing refits differ from the printed summaries | A numerical answer can be internally reproducible yet not match raw data. | Publish inputs and distinguish the two calculations. |
| Q = 4.22 follows Box–Pierce with the printed ACF | Calling every Q statistic Ljung–Box causes confusion. | State the actual formula. |
| Industrial-production p-values straddle 0.05 | A binary conclusion hides close results. | Report p-values, lag order, and sample size together. |
| The stock-index table contains a suspicious observation | An outlier or typo may dominate a very short sample. | Verify the source before correcting or interpreting. |
A practical workflow for a defensible Granger causality analysis
Start with the economic question and the data calendar
Define the dependent variable and forecast horizon before fitting models. For example, predicting next month’s consumption differs from predicting its quarterly growth rate. Also, align observations to a consistent calendar and verify the measurement units.
Release timing matters as much as observation timing. A variable dated January might become available after the information needed for a February forecast. Therefore, an operational forecasting study should use the information actually available at each forecast origin.
Check transformations and model assumptions
Plot the levels and candidate transformations, then investigate trends, seasonality, breaks, and extreme observations. Next, assess integration and cointegration with explicit deterministic terms. Avoid selecting transformations solely because they produce a favorable Granger p-value.
After fitting the equation, examine residual autocorrelation and stability. Moreover, consider whether common shocks or additional predictors could explain the apparent directional relationship. A bivariate result is conditional on a much narrower information set than a larger system.
Separate exploratory lag searches from confirmatory tests
Choose a plausible maximum lag using the data frequency and sample size. Then compare candidate specifications consistently. If you search many lag orders, directions, variables, and significance thresholds, disclose that exploration rather than presenting the most favorable result as a single planned test.
For small samples, inspect a corrected information criterion as well as the ordinary AIC. However, preserve consistent parameter counting and a common data window. No penalty formula repairs a comparison that mixes incompatible likelihoods.
Further reading: Forecasting: Principles and Practice, time series regression, including predictor selection and the small-sample motivation for AICc.
Validate the forecasting claim outside the estimation sample
A rejection in the estimation sample does not guarantee improved future forecasts. Therefore, compare the restricted and unrestricted models in a chronological holdout or rolling-origin evaluation. Match that evaluation to the intended horizon and available information.
For example, compare root mean squared forecast errors across the same forecast dates. Also, inspect whether a small average improvement hides poor performance during important periods. This extends the analysis beyond a single significance threshold.
Write a result that another analyst can reproduce
Reproduce the calculations and inspect the source data
The following code uses NumPy and SciPy to estimate the two nested OLS equations directly. It accepts already-stationary or appropriately transformed inputs. Therefore, it reproduces the regression arithmetic without silently applying a unit root test or choosing a transformation.
For the common-sample housing comparison, start every regression at index six of the differenced data. Conversely, use the default start to reproduce the changing-sample specification. The functions return the residual sums of squares, criterion, F-statistic, and p-value.
import numpy as np
from scipy.stats import f
def granger_ols(target, predictor, p, start=None):
target = np.asarray(target, dtype=float)
predictor = np.asarray(predictor, dtype=float)
assert len(target) == len(predictor)
assert np.isfinite(target).all() and np.isfinite(predictor).all()
start = p if start is None else start
assert start >= p >= 1
stop = len(target)
y = target[start:]
own = [target[start-j:stop-j] for j in range(1, p+1)]
other = [predictor[start-j:stop-j] for j in range(1, p+1)]
restricted = np.column_stack([np.ones(len(y))] + own)
unrestricted = np.column_stack([np.ones(len(y))] + own + other)
n, k = unrestricted.shape
assert n > k
assert np.linalg.matrix_rank(unrestricted) == k
def fit(z):
b = np.linalg.lstsq(z, y, rcond=None)[0]
residual = y - z @ b
return b, float(residual @ residual)
b, sse_u = fit(unrestricted)
_, sse_r = fit(restricted)
F = ((sse_r - sse_u) / p) / (sse_u / (n - k))
return dict(n=n, k=k, coefficients=b, sse_u=sse_u,
sse_r=sse_r, F=F, p_value=f.sf(F, p, n-k),
criterion=np.log(sse_u/n) + 2*k/n)
# Housing arrays below contain LEVELS: difference once.
for p in range(1, 7):
print(p, granger_ols(np.diff(consumption),
np.diff(housing), p, start=6))
# Industrial arrays below already contain FIRST DIFFERENCES.
for p in (1, 6):
print(p, granger_ols(production_change, stock_change, p))
In statsmodels, testing whether housing predicts consumption requires consumption in the first column and housing in the second column. Also, an integer maximum lag requests tests through that order, whereas a list can request selected orders. Confusing either convention can produce a valid computation that answers the wrong question.
Open the exact arrays used in the recalculation
These numbers transcribe the visible tables. The housing refits do not reproduce all printed regression summaries, as documented above. Meanwhile, the industrial-production calculations reproduce the key reported answers.
housing = np.array([
1355, 1486, 1457, 1492, 1442, 1494, 1437, 1390,
1546, 1520, 1510, 1566, 1525, 1584, 1567, 1540,
1536, 1641, 1598, 1614, 1582, 1715, 1660, 1792,
1804, 1738, 1737, 1561, 1649, 1562, 1704, 1657,
1628, 1636, 1663, 1769
])
consumption = np.array([
5342.1, 5351.2, 5358.7, 5368.2, 5361.5, 5397.4, 5454.0, 5464.9,
5467.3, 5484.8, 5506.5, 5530.0, 5540.8, 5573.0, 5603.5, 5609.8,
5658.4, 5686.4, 5685.9, 5708.7, 5738.4, 5758.3, 5771.5, 5809.5,
5817.9, 5854.5, 5908.4, 5915.8, 5928.4, 5976.6, 5987.1, 6020.4,
6033.9, 6062.1, 6090.8, 6150.0
])
production_change = np.array([
1.87, 1.26, -2.0, 0.68, 2.62, -4.56, 7.95, 0.53,
-0.12, -3.48, -1.75, 1.13, 1.88, 2.12, -2.09, 0.41,
4.46, -4.99, 7.12, 0.75, 0.61, -2.49, -1.39, 1.37,
2.16, 2.03, -1.34, 0.2, 5.43, -6.76, 7.66, 1.65,
-1.57, -3.01, -2.19
])
stock_change = np.array([
69.05, 52.41, 10.0, -20.93, 43.02, -13.16, -163.39, 59.73,
81.66, 64.96, 65.59, 50.41, -41.31, 48.04, 48.81, -33.34,
70.87, -43.99, -8.3, -37.7, 80.22, 26.13, 80.18, -74.79,
-28.04, 132.16, -46.14, -31.83, 34.0, -23.77, 86.85, -81.17,
-7.1, -114.45, 5.32
])
temperature = np.array([
12.71, 12.33, 12.58, 12.24, 13.34, 12.59, 12.91, 12.15,
12.13, 11.5, 13.06, 13.16, 12.7, 12.41, 13.68, 12.63,
13.53, 12.68, 12.83, 12.84, 12.67, 13.23, 14.04, 12.27,
13.07, 12.66, 13.39, 11.66, 12.51, 13.73, 13.75
])
innovations = np.array([
0.1291, 0.691, 0.1348, 0.151, 0.3869, 0.7318, 0.4515, 0.3334,
0.8943, 0.0773, 0.4303, 0.5805, 0.925, 0.0408, 0.9621, 0.6577,
0.8292, 0.5996, 0.4197, 0.8095, 0.4661, 0.2208, 0.2334, 0.5894,
0.8296, 0.4352, 0.1958, 0.6074, 0.7228, 0.3146
])
Open the housing and consumption data table: 1997–1999
| Month | Housing starts Y | Real consumption X |
|---|---|---|
| 1997-01 | 1355 | 5342.1 |
| 1997-02 | 1486 | 5351.2 |
| 1997-03 | 1457 | 5358.7 |
| 1997-04 | 1492 | 5368.2 |
| 1997-05 | 1442 | 5361.5 |
| 1997-06 | 1494 | 5397.4 |
| 1997-07 | 1437 | 5454.0 |
| 1997-08 | 1390 | 5464.9 |
| 1997-09 | 1546 | 5467.3 |
| 1997-10 | 1520 | 5484.8 |
| 1997-11 | 1510 | 5506.5 |
| 1997-12 | 1566 | 5530.0 |
| 1998-01 | 1525 | 5540.8 |
| 1998-02 | 1584 | 5573.0 |
| 1998-03 | 1567 | 5603.5 |
| 1998-04 | 1540 | 5609.8 |
| 1998-05 | 1536 | 5658.4 |
| 1998-06 | 1641 | 5686.4 |
| 1998-07 | 1598 | 5685.9 |
| 1998-08 | 1614 | 5708.7 |
| 1998-09 | 1582 | 5738.4 |
| 1998-10 | 1715 | 5758.3 |
| 1998-11 | 1660 | 5771.5 |
| 1998-12 | 1792 | 5809.5 |
| 1999-01 | 1804 | 5817.9 |
| 1999-02 | 1738 | 5854.5 |
| 1999-03 | 1737 | 5908.4 |
| 1999-04 | 1561 | 5915.8 |
| 1999-05 | 1649 | 5928.4 |
| 1999-06 | 1562 | 5976.6 |
| 1999-07 | 1704 | 5987.1 |
| 1999-08 | 1657 | 6020.4 |
| 1999-09 | 1628 | 6033.9 |
| 1999-10 | 1636 | 6062.1 |
| 1999-11 | 1663 | 6090.8 |
| 1999-12 | 1769 | 6150.0 |
Frequently asked questions about Granger causality
Does Granger causality prove causation?
No. It establishes evidence about incremental prediction under a specified information set and model. However, causal intervention claims require additional assumptions or a suitable research design.
Can both variables Granger-cause each other?
Yes, both directions can show predictive content. Nevertheless, test each direction separately because swapping the dependent variable changes the equation and null hypothesis.
What does a p-value above 0.05 mean?
It means the test does not reject its null at the 5% level under the specified assumptions. Therefore, report insufficient evidence rather than claiming proof of no relationship.
Should I always difference economic time series?
No. Stationary levels may need no differencing, while cointegrated I(1) variables call for an error-correction treatment. Also, a deterministic trend requires a different analysis from a stochastic trend.
How many lags should I use?
There is no universal lag count. Instead, consider the frequency, sample size, plausible timing, comparable information criteria, and residual diagnostics. Report sensitivity when conclusions depend on the lag choice.
Why does my software give a different AIC?
Packages may use different constants, parameter counts, samples, or likelihood conditioning. Consequently, verify the underlying definition before interpreting numerical differences as estimation errors.
Is a significant negative coefficient a negative causal effect?
No. It describes a conditional predictive association in that equation. Moreover, several lag coefficients can combine into a response that changes sign over time.
Can I test levels when the variables have unit roots?
A standard stationary-model F-test should not automatically carry over to an arbitrary nonstationary levels regression. Instead, investigate integration, cointegration, and an inference method appropriate to the system before interpreting its p-value.
Conclusion: treat the model, data, and interpretation as one argument
Granger causality provides a structured way to ask whether one series helps predict another. However, the answer depends on transformations, lag order, sample dates, and the information set. A correct F calculation cannot compensate for a poorly defined model.
The exercises show why careful replication matters. In particular, the housing example exposes inconsistencies between printed summaries and visible inputs, while the industrial-production example shows how close p-values can support different threshold decisions. Meanwhile, the AR, MA, and temperature exercises build the mechanics needed to audit those models.
Before publishing a result, identify the exact data, document the specification, preserve numerical precision, and explain the limits of the conclusion. That approach turns a textbook answer into an analysis that readers can check and use.
