Time Series Econometrics
xβ own and xβ cross response Β· A = [0.5, aββ; 0.6, 0.5], Ξ£ = I Β· largest |eigenvalue| ${modulus.toFixed(3)} Β· roots are ${disc < 0 ? "complex: oscillation" : "real: monotone decay"}${modulus >= 1 ? " Β· not stable" : ""}
` ``` --- ## βοΈ Worked Example: IRF by Hand ::: {style="font-size: 28px"} VAR(1) for $(\Delta m_t,\Delta p_t)'$: $A=\begin{pmatrix}0.5&0.2\\0.3&0.4\end{pmatrix}$, $\Sigma=\begin{pmatrix}1&0.5\\0.5&2\end{pmatrix}$. ::: {.fragment} **Factor.** $L=\begin{pmatrix}1&0\\0.5&1\end{pmatrix}$, $D=\operatorname{diag}(1,\ 2-0.5^2)=\operatorname{diag}(1,1.75)$. Check: $LDL'=\Sigma$. ::: ::: {.fragment} **Responses.** $\Psi_0=L$; $\ \Psi_1=AL=\begin{pmatrix}0.5+0.1&0.2\\0.3+0.2&0.4\end{pmatrix}=\begin{pmatrix}0.6&0.2\\0.5&0.4\end{pmatrix}$. A money shock raises inflation by 0.5 on impact and 0.5 a quarter later; a price shock leaves $\Delta m$ unchanged on impact and raises it by 0.2 after one quarter. ::: ::: {.fragment} **Granger.** The (1,2) element of $A$ is $0.2\neq0$, so $\Delta p$ Granger-causes $\Delta m$; and $0.3\neq0$, so $\Delta m$ Granger-causes $\Delta p$. ::: ::: --- ## π Variance Decomposition (FEVD) ::: {style="font-size: 28px"} The $h$-step forecast error is $x_{t+h}-E_tx_{t+h}=\sum_{s=0}^{h-1}\Psi_sv_{t+h-s}$. With orthogonal shocks, $$\operatorname{var}_t(x_{j,t+h})=\sum_{i=1}^m f_{ij}(h)\,\sigma_i^2,\qquad f_{ij}(h)=\sum_{s=0}^{h-1}(\Psi_s)_{ji}^2 .$$ These are the diagonal elements of $\operatorname{var}_t x_{t+h}$. The **share of shock $i$** in variable $j$: $$\omega_{ij}(h)=\frac{f_{ij}(h)\,\sigma_i^2}{\sum_{k=1}^m f_{kj}(h)\,\sigma_k^2}\times100\%,\qquad \sum_i\omega_{ij}(h)=100\% .$$ ::: {.fragment} **Example above, $h=2$, $j=\Delta m$:** $f_{11}=1^2+0.6^2=1.36$, $f_{21}=0^2+0.2^2=0.04$. Share of the money shock $=\dfrac{1.36}{1.36+0.04\times1.75}=\dfrac{1.36}{1.43}=95.1\%$. ::: ::: --- ## π» The Example in R ```{r} #| label: worked-var A <- matrix(c(0.5, 0.2, 0.3, 0.4), 2, byrow = TRUE) Sigma <- matrix(c(1, 0.5, 0.5, 2), 2) C <- t(chol(Sigma)) # Sigma = C C', C lower triangular L <- C %*% diag(1 / diag(C)) # unit diagonal D <- diag(diag(C)^2) # Sigma = L D L' Psi <- list(L, A %*% L) # Psi_0, Psi_1 f <- Reduce(`+`, lapply(Psi, function(P) P^2)) # f_ij(2), rows = variable j share <- 100 * (f %*% D) / rowSums(f %*% D) dimnames(share) <- list(c("dm", "dp"), c("v_m", "v_p")) Psi_01 <- cbind(Psi[[1]], Psi[[2]]); colnames(Psi_01) <- c("Psi0_m", "Psi0_p", "Psi1_m", "Psi1_p") list(D = diag(D), Psi_0_and_Psi_1 = Psi_01, FEVD_h2_percent = round(share, 1)) ``` --- ## π Granger Causality ::: {style="font-size: 29px"} ::: {.callout-note} ## Definition $\Delta m$ does **not** Granger-cause $\Delta p$ ($\Delta m\overset{GC}{\not\to}\Delta p$) if $$E(\Delta p_{t+1}\mid \Delta p_t,\Delta p_{t-1},\dots,\Delta m_t,\Delta m_{t-1},\dots)=E(\Delta p_{t+1}\mid \Delta p_t,\Delta p_{t-1},\dots).$$ ::: In a VAR: all $p$ coefficients on lagged $\Delta m$ in the $\Delta p$ equation are zero. Test this with an $F$ (or Wald) statistic, $F\sim F(p,\,T-mp-1)$ approximately. ::: {.fragment} It is a **statistical** notion, not an economic one. Stock prices Granger-cause dividends, because prices are forward-looking; economically, it is expected dividends that determine prices. ::: ::: --- ## π§ Think-Pair-Share ```{r} #| label: tps-timer #| echo: false countdown::countdown(minutes = 3, seconds = 0, top = 0, right = 0, font_size = "2em", warn_when = 30) ``` ::: {style="font-size: 30px"} A VAR(1) has $A=\begin{pmatrix}0.7&0\\0.5&0.2\end{pmatrix}$ and $\Sigma=\begin{pmatrix}1&0.3\\0.3&1\end{pmatrix}$, with $x_1$ ordered first. 1. Is the VAR stationary? 2. Which variable Granger-causes which? 3. Compute $\Psi_0$ and $\Psi_1$. 4. What share of the $h$-step forecast error variance of $x_1$ comes from $v_2$, for any $h$? ::: --- ## β
Think-Pair-Share: Solution ::: {style="font-size: 28px"} 1. $A$ is triangular: eigenvalues $0.7$ and $0.2$, both inside the unit circle. **Stationary.** 2. $a_{21}=0.5\neq0$: $x_1$ Granger-causes $x_2$. $a_{12}=0$: $x_2$ does **not** Granger-cause $x_1$. 3. $L=\begin{pmatrix}1&0\\0.3&1\end{pmatrix}$, $D=\operatorname{diag}(1,0.91)$. $\Psi_0=L$, $\Psi_1=AL=\begin{pmatrix}0.7&0\\0.56&0.2\end{pmatrix}$. ::: {.fragment} 4. **Zero.** $x_1$ does not respond to $v_2$ on impact (it is ordered first) and $A$ never feeds $x_2$ back into $x_1$, so $(\Psi_s)_{12}=0$ for all $s$. Its forecast errors are entirely its own shock. ::: ::: --- ## β οΈ Spurious Regression ::: {style="font-size: 29px"} Let $x_t,y_t\sim I(1)$ and run $y_t=\beta x_t+\varepsilon_t$. There are two possibilities. 1. $\varepsilon_t\sim I(1)$ **for whatever $\beta$**: a **spurious regression**. OLS reports a large $R^2$ and a "significant" $t$, though there is no relation at all. Granger and Newbold's example: rabbits in Australia against US GNP. 2. $\varepsilon_t\sim I(0)$ **for a certain $\beta$**: a **cointegrating regression**. ::: {.fragment} In case 1 the $t$-statistic does not converge to $N(0,1)$; it diverges at rate $\sqrt T$, and $R^2$ converges to a random variable instead of to zero (Phillips, 1986). More data makes the illusion **stronger**. ::: ::: --- ## π» Two Independent Random Walks ```{r} #| label: spurious-sim spurious <- function(n, reps = 1000) { out <- replicate(reps, { x <- cumsum(rnorm(n)); y <- cumsum(rnorm(n)) # independent I(1) s_lev <- summary(lm(y ~ x)); s_dif <- summary(lm(diff(y) ~ diff(x))) c(t_lev = s_lev$coefficients[2, 3], R2_lev = s_lev$r.squared, t_dif = s_dif$coefficients[2, 3]) }) c(T = n, reject_levels = mean(abs(out["t_lev", ]) > 1.96), median_R2 = median(out["R2_lev", ]), median_abs_t = median(abs(out["t_lev", ])), reject_differences = mean(abs(out["t_dif", ]) > 1.96)) } sims <- t(sapply(c(50, 100, 200, 400, 800), spurious)) round(sims, 3) ``` ::: {style="font-size: 28px"} A nominal 5% test rejects two times in three at $T=50$ and nine times in ten at $T=800$; the median $|t|$ grows like $\sqrt T$. In first differences it rejects about 5% of the time, as it should. ::: --- ## π The $t$-Statistic Explodes ```{r} #| label: spurious-fig #| echo: false #| fig-width: 11 #| fig-height: 4.2 tstats <- do.call(rbind, lapply(c(100, 800), function(n) data.frame( T = paste("T =", n), t = replicate(1000, { x <- cumsum(rnorm(n)); y <- cumsum(rnorm(n)) summary(lm(y ~ x))$coefficients[2, 3] })))) tstats$T <- factor(tstats$T, levels = c("T = 100", "T = 800")) ggplot(tstats, aes(t)) + geom_histogram(aes(y = after_stat(density)), bins = 60, fill = "#8ba3c7") + stat_function(fun = dnorm, n = 1001, colour = "#8b2635", linewidth = 1.2) + geom_vline(xintercept = c(-1.96, 1.96), linetype = "dashed") + facet_wrap(~ T) + coord_cartesian(xlim = c(-60, 60)) + labs(x = "t-statistic on x in the levels regression", y = "density") + theme(strip.text = element_text(size = 20, face = "bold")) ``` ::: {style="font-size: 28px"} The red curve is the $N(0,1)$ the textbook test assumes, squeezed inside the dashed lines. The actual spread is many times wider, and grows by about $\sqrt8\approx2.8$ from $T=100$ to $T=800$. ::: --- ## π The Drunk and Her Dog ```{r} #| label: drunk-dog #| echo: false #| fig-width: 11 #| fig-height: 3.8 set.seed(7) n <- 300; drunk <- cumsum(rnorm(n)) gap <- as.numeric(arima.sim(list(ar = 0.7), n, sd = 0.8)) dog <- drunk + gap walk <- data.frame(t = rep(1:n, 3), value = c(drunk, dog, gap - 12), who = rep(c("drunk x_t", "dog y_t", "leash y_t - x_t (shifted down)"), each = n)) ggplot(walk, aes(t, value, colour = who)) + geom_line(linewidth = 0.9) + scale_colour_manual(values = c("#2f4b7c", "#c0392b", "grey40"), name = NULL) + labs(x = "t", y = NULL) + theme(legend.position = "top") ``` ::: {style="font-size: 28px"} Each wanders like a random walk, but the leash keeps the distance $y_t-x_t$ stationary: $x_t,y_t\sim I(1)$ while $y_t-\beta x_t\sim I(0)$ with $\beta=1$. That is **cointegration**. ::: --- ## π Cointegration in Economics ::: {style="font-size: 29px"} **Known $\beta$**, implied by theory: - a constant consumptionβincome ratio: $\log c_t-\log y_t\sim I(0)$, so the regression of $\log c$ on $\log y$ should be cointegrating - purchasing power parity: the real exchange rate $e_t+p_t^*-p_t\sim I(0)$ - long and short interest rates: $R_t-r_t\sim I(0)$ - spot and futures exchange rates ::: {.fragment} **Unknown $\beta$**, to be estimated. Money demand: $$m_t-p_t=\alpha_1+\alpha_2y_t+\alpha_3r_t+I(0).$$ ::: ::: --- ## π° Case Study: A Monetary Policy VAR :::: columns ::: {.column width="50%"} ::: {.callout-note} ## π The question A VAR in inflation, unemployment and the federal funds rate, in the spirit of Stock and Watson (2001). 1. Which lag order do the criteria choose? 2. Does the policy rate Granger-cause inflation, and inflation the rate? 3. How do the three variables respond to a policy-rate shock? ::: ::: ::: {.column width="50%"} ::: {.callout-tip} ## π The data **Series**: CPI (`CPIAUCSL`), unemployment rate (`UNRATE`), effective fed funds rate (`FEDFUNDS`), FRED **Period**: 1960Q1β2007Q4, quarterly averages; ends before the zero lower bound $\pi_t=400\,\Delta\log\text{CPI}_t$; ordering $(\pi_t,u_t,R_t)$ ::: ::: :::: ```{.r .display-only} base <- "https://fred.stlouisfed.org/graph/fredgraph.csv?id=" raw <- lapply(c("CPIAUCSL", "UNRATE", "FEDFUNDS"), function(i) read.csv(paste0(base, i))) ``` --- ## π° Case Study: Estimation in R ```{r} #| label: case-var #| code-summary: "π¦ Data, lag order, OLS, Granger tests" qavg <- function(id) { z <- read.csv(paste0("data/", id, ".csv")); d <- as.Date(z$observation_date) tapply(z[[id]], paste(format(d, "%Y"), (as.numeric(format(d, "%m")) - 1) %/% 3 + 1), mean) } cpi <- qavg("CPIAUCSL"); un <- qavg("UNRATE"); ff <- qavg("FEDFUNDS") Y <- cbind(infl = 400 * diff(log(cpi)), u = un[-1], R = ff[-1]) # 1960Q1-2007Q4 m <- ncol(Y) var_ols <- function(Y, p) { # OLS equation by equation n <- nrow(Y); Z <- cbind(1, do.call(cbind, lapply(1:p, function(j) Y[(p + 1 - j):(n - j), ]))) X <- Y[(p + 1):n, ]; B <- solve(crossprod(Z), crossprod(Z, X)); E <- X - Z %*% B list(B = B, E = E, Z = Z, X = X, Sigma = crossprod(E) / nrow(E), p = p) } crit <- sapply(1:6, function(p) { f <- var_ols(Y[(7 - p):nrow(Y), ], p); Tn <- nrow(f$E) ld <- log(det(f$Sigma)); c(AIC = ld + 2 * p * m^2 / Tn, BIC = ld + log(Tn) * p * m^2 / Tn, HQ = ld + 2 * log(log(Tn)) * p * m^2 / Tn) }) p <- which.min(crit["BIC", ]); fit <- var_ols(Y, p) granger <- function(fit, cause, eq) { # F test: lags of `cause` out of `eq` drop <- 1 + (0:(fit$p - 1)) * m + cause; y <- fit$X[, eq] rss <- function(Z) sum(lm.fit(Z, y)$residuals^2); df2 <- nrow(fit$Z) - ncol(fit$Z) F <- ((rss(fit$Z[, -drop]) - rss(fit$Z)) / fit$p) / (rss(fit$Z) / df2) pf(F, fit$p, df2, lower.tail = FALSE) } list(order_chosen = c(AIC = which.min(crit["AIC", ]), BIC = p, HQ = which.min(crit["HQ", ])), granger_p_values = c(R_to_infl = granger(fit, 3, 1), infl_to_R = granger(fit, 1, 3), u_to_R = granger(fit, 2, 3)) |> signif(2)) ``` ::: {style="font-size: 28px"} BIC and HQ choose $p=3$; the VAR(3) is used below. Small $p$-values reject Granger non-causality in all three directions. ::: --- ## π° Case Study: A Rate Shock ```{r} #| label: case-irf #| echo: false #| fig-width: 10 #| fig-height: 3.9 irf_var <- function(B, Sigma, p, H = 20) { A <- lapply(1:p, function(j) t(B[1 + (j - 1) * m + 1:m, ])) Phi <- list(diag(m)) for (h in 1:H) Phi[[h + 1]] <- Reduce(`+`, lapply(1:min(h, p), function(j) A[[j]] %*% Phi[[h + 1 - j]])) P <- t(chol(Sigma)); lapply(Phi, function(F) F %*% P) # one-s.d. orthogonal shocks } point <- sapply(irf_var(fit$B, fit$Sigma, p), function(M) M[, 3]) # responses to shock in R boot <- replicate(300, { # residual bootstrap Yb <- Y; E <- fit$E[sample(nrow(fit$E), replace = TRUE), ] for (t in (p + 1):nrow(Y)) Yb[t, ] <- c(1, as.vector(t(Yb[(t - 1):(t - p), ]))) %*% fit$B + E[t - p, ] fb <- var_ols(Yb, p); sapply(irf_var(fb$B, fb$Sigma, p), function(M) M[, 3]) }) lab <- c("inflation (pp)", "unemployment (pp)", "fed funds rate (pp)"); lab2 <- c("inflation", "unemployment", "fed funds rate") df <- do.call(rbind, lapply(1:3, function(k) data.frame(h = 0:20, var = lab[k], irf = point[k, ], lo = apply(boot[k, , ], 1, quantile, 0.05), hi = apply(boot[k, , ], 1, quantile, 0.95)))) df$var <- factor(df$var, levels = lab) ggplot(df, aes(h, irf)) + geom_hline(yintercept = 0, colour = "grey50") + geom_ribbon(aes(ymin = lo, ymax = hi, fill = "90% bootstrap band"), alpha = 0.5) + geom_line(aes(colour = "Cholesky IRF"), linewidth = 1.2) + facet_wrap(~ var, scales = "free_y") + scale_fill_manual(values = "#8ba3c7", name = NULL) + scale_colour_manual(values = "#8b2635", name = NULL) + labs(x = "quarters after the shock", y = NULL) + theme_minimal(base_size = 20) + theme(legend.position = "top") ``` ::: {style="font-size: 28px"} **What the data say:** after a one-s.d. rate shock the rate decays over about two years and unemployment rises, peaking after two years, as theory predicts. Inflation **rises** at first and turns negative only after six quarters: the "price puzzle". A three-variable Cholesky VAR omits information the Fed reacted to (commodity prices, expectations), so part of the "shock" anticipates inflation. ::: --- ## π° Case Study: Variance Decomposition ```{r} #| label: case-fevd #| echo: false #| fig-width: 10 #| fig-height: 3.6 Psi <- irf_var(fit$B, fit$Sigma, p, H = 24) fevd <- function(h) { f <- Reduce(`+`, lapply(Psi[1:h], function(M) M^2)); 100 * f / rowSums(f) } fd <- do.call(rbind, lapply(1:24, function(h) { s <- fevd(h) data.frame(h = h, var = rep(lab2, 3), shock = rep(c("inflation shock", "unemployment shock", "rate shock"), each = 3), share = as.vector(s)) })) fd$var <- factor(fd$var, levels = lab2) fd$shock <- factor(fd$shock, levels = c("inflation shock", "unemployment shock", "rate shock")) ggplot(fd, aes(h, share, fill = shock)) + geom_area() + facet_wrap(~ var) + scale_x_continuous(breaks = c(1, 8, 16, 24)) + scale_fill_manual(values = c("#8b2635", "#d9a441", "#2f4b7c"), name = NULL) + labs(x = "horizon h (quarters)", y = "share of FEV (%)") + theme_minimal(base_size = 20) + theme(legend.position = "top", panel.spacing = unit(2, "lines")) ``` ::: {style="font-size: 28px"} Each panel stacks $\omega_{ij}(h)$ over the three shocks to 100%. Inflation is driven mostly by its own shock; beyond a year, **unemployment** shocks carry the largest share of the rate's forecast errors: the systematic reaction of policy to the real economy. ::: --- ## π Quiz #1: Granger Causality {.quiz-question} In a bivariate VAR(2) for $(y_t,z_t)$, which restriction says $z$ does **not** Granger-cause $y$? - [The coefficients on $z_{t-1}$ and $z_{t-2}$ in the $y$ equation are both zero]{.correct data-explanation="β
z does not Granger-cause y when past z adds nothing to the forecast of y given past y: both lag coefficients of z in the y equation are zero. Test with F(2, T β 5)."} - The coefficients on $y_{t-1}$ and $y_{t-2}$ in the $z$ equation are both zero - The correlation between the two innovations is zero - The impact response of $y$ to the $z$ shock is zero --- ## π Quiz #2: Variance Decomposition {.quiz-question} For variable $j$ at horizon $h$, $f_{1j}(h)\sigma_1^2=3$ and $f_{2j}(h)\sigma_2^2=1$. What is the share of shock 2 in its forecast error variance? - [$25\%$]{.correct data-explanation="β
Ο_2j(h) = 1 / (3 + 1) = 25%. The shares over the shocks sum to 100%."} - $33\%$ - $75\%$ - $1\%$ --- ## π Quiz #3: Spurious or Not? {.quiz-question} $x_t$ and $y_t$ are independent random walks, $T=1000$. You regress $y_t$ on $x_t$ in levels. What do you expect? - [A large $|t|$-statistic and a sizeable $R^2$, although there is no relation]{.correct data-explanation="β
This is the spurious regression: the error is I(1) for every Ξ², the t-statistic diverges at rate βT and RΒ² has a non-degenerate limit. Regress first differences, or test for cointegration, instead."} - A $t$-statistic close to $N(0,1)$, rejecting about 5% of the time - $R^2$ close to zero, since the series are independent - A cointegrating regression with $\beta=0$ --- ## π Key Formulas ::: {style="font-size: 28px"} | | | |---|---| | VAR($p$) | $x_t=a_0+\sum_{j=1}^pA_jx_{t-j}+\varepsilon_t$, $\operatorname{var}\varepsilon_t=\Sigma$ | | Stationarity | $\det(I_m-A_1z-\dots-A_pz^p)\neq0$ for $\lvert z\rvert\le1$ | | Lag choice | $\text{BIC}(p)=\log\det\hat\Sigma_p+pm^2\log T/T$ | | Orthogonal shocks | $\Sigma=LDL'$, $v_t=L^{-1}\varepsilon_t$, $\Psi_s=\Phi_sL$ | | IRF | $\mathrm{IRF}(h)=\partial x_{t+h}/\partial v_t'=\Psi_h$ | | FEVD | $\omega_{ij}(h)=f_{ij}(h)\sigma_i^2/\sum_kf_{kj}(h)\sigma_k^2$, $f_{ij}(h)=\sum_{s