---
title: "2.1.6 Kalman filter for dynamic hedge ratios"
description: "A Kalman filter from scratch in NumPy for a changing hedge ratio, versus rolling OLS and the true β; noise tuning and a check with statsmodels."
url: "https://pythonparatrading.com/en/python-trading-course/level-2/kalman-filter-dynamic-hedge-ratio/"
lang: "en"
image: "https://pythonparatrading.com/og/en/python-trading-course/level-2/kalman-filter-dynamic-hedge-ratio.png"
---

# 2.1.6 · Kalman filter for dynamic hedge ratios

> When the relationship between two prices changes over time, a fixed regression goes stale and a rolling regression forces you to choose a window. The Kalman filter treats the hedge ratio as a hidden state that is updated with every session. Here it is coded from scratch and compared with the truth, which in this synthetic data is known day by day.

- Course: [Level 2: Systematic strategies and rigorous backtesting](https://pythonparatrading.com/en/python-trading-course/level-2/)
- Module: 2.1 · Statistics and econometrics for trading
- Duration: 65 min
- Prerequisites: [0.7.3 · Intuitive linear algebra: vectors and matrices as groundwork for portfolios](https://pythonparatrading.com/en/python-trading-course/level-0/linear-algebra-for-portfolios/), [1.6.1 · Returns, volatility, rolling correlations and drawdown](https://pythonparatrading.com/en/python-trading-course/level-1/returns-volatility-correlation-and-drawdown/)

## By the end of this lesson

- Set up a time-varying hedge ratio as a state space model: state equation, observation equation and their two variances.
- Code in NumPy the prediction and update steps of the Kalman filter and use only the estimate available the previous day.
- Compare the filter with the rolling regression against the true β and choose the state noise by looking at the bias-variance tradeoff.
- Estimate the variances by maximum likelihood with `statsmodels` and recognize the look-ahead of smoothing.

## A hedge ratio that does not stay still

In lesson 2.1.5 pair D passed the Engle-Granger test with a fixed β, but the residual of that β reverted four times
more slowly than the true spread: its relationship is $\ln y_t = a_t + \beta_t \ln x_t + u_t$ with a $\beta_t$ that
**changes every day**. A regression with a fixed β for twelve years describes an average relationship that never
existed. It is a common real-world situation: the sensitivity of a stock to its
sector, that of an ETF to its index or that of one bond to another of a different maturity change with the structure
of the company, with interest rates or with the composition of the index.

The hedge ratio matters for two reasons. With the relationship in logarithms, β is an elasticity: to **hedge** a
position of 1 USD in `y` you have to sell about $\beta$ USD of `x`, and a wrong β leaves risk unhedged. To
**trade a pair** (lesson 2.3.3), the spread $\ln y - \beta \ln x$ is the signal: a wrong β puts moves into the spread
that are not reversion.

The data for this lesson is pair D from `pares2.parquet`, synthetic, with the series of the true β for each session in
`pares2_beta_real.csv`: a reflected random walk with a daily step of standard deviation 0.006, and an
Ornstein-Uhlenbeck spread with a half-life of 10 and a standard deviation of 0.025 (2.5%). The relationship pivots
around the first price of `x`, $x_0$: $\ln y_t = \alpha + \ln x_0 + \beta_t \ln(x_t / x_0) + u_t$, so in the usual form
the constant $a_t = \alpha + (1 - \beta_t) \ln x_0$ also moves with β. Everything is done in logarithms. Knowing the
truth makes it possible to measure the error of each method,
something that can never be done with real data.

```python
import numpy as np
import pandas as pd
import statsmodels.api as sm

pares = pd.read_parquet("pares2.parquet")
precios = pares.pivot(index="fecha", columns="activo",
                      values="precio")
y, x = np.log(precios["PZDY"]), np.log(precios["PZDX"])
beta = pd.read_csv("pares2_beta_real.csv", index_col="fecha",
                   parse_dates=True)["beta_real"]
print(f"true beta: from {beta.min():.3f} to {beta.max():.3f}, "
      f"daily step {beta.diff().std():.4f}")

fijo = sm.OLS(y, sm.add_constant(x)).fit().params
error = np.sqrt(((fijo.iloc[1] - beta) ** 2).mean())
print(f"Full-sample OLS: alpha {fijo.iloc[0]:.2f}, "
      f"beta {fijo.iloc[1]:.3f}, root mean squared error {error:.3f}")
```
```output
true beta: from 0.976 to 1.489, daily step 0.0061
Full-sample OLS: alpha -1.03, beta 1.337, root mean squared error 0.166
```

The true β ranges from 0.98 to 1.49. The full-sample regression gives a β of 1.337 and a constant of −1.03: the
regression fits the average relationship of the prices, it does not estimate that of any particular day. Its root
mean squared error against the true β is 0.166. Everything that follows is measured with that yardstick: the root
mean squared error (RMSE) of the estimated β against the true one.

## The rolling regression and its window

The usual answer is a regression over a rolling window: a β is estimated with the last $n$ days, and the next day the
window shifts. There are two ways to do it, on prices or on price changes, and they behave very
differently:

```python
def rmse(est):
    """Error of the beta available the previous day versus the true one."""
    e = (est.shift(1) - beta).iloc[500:]
    return np.sqrt((e ** 2).mean())

dy, dx = y.diff(), x.diff()
filas = {}
for n in [10, 20, 40, 60, 120, 250, 500]:
    en_precios = y.rolling(n).cov(x) / x.rolling(n).var()
    en_variac = ((dy * dx).rolling(n).sum()
                 / (dx ** 2).rolling(n).sum())
    filas[n] = {"precios": rmse(en_precios),
                "variaciones": rmse(en_variac)}
pd.DataFrame(filas).T.rename_axis("ventana").round(4)
```
```output
         precios  variaciones
ventana
10        0.5309       0.2841
20        0.4458       0.1845
40        0.3591       0.1336
60        0.3189       0.1060
120       0.2379       0.0753
250       0.1500       0.0704
500       0.1495       0.0900
```

The error is measured from session 500 onward, so that all windows have data, and always with the previous day's β
(`shift(1)`): today's hedge is decided with what was known yesterday.

On prices (their logarithms, with a constant), the regression fails: the error goes from 0.53 with 10 sessions to
0.15 with 500, never better than a fixed β. The reason is the constant. In a 120-session window the log price of
`PZDX` moves with a standard deviation of about 5% around a level of 4, and in that short stretch many combinations of
α and β fit almost equally well: a higher β with a lower α. The two estimates
offset each other and each one on its own is noise:

```python
b120 = y.rolling(120).cov(x) / x.rolling(120).var()
a120 = y.rolling(120).mean() - b120 * x.rolling(120).mean()
print(f"correlation between rolling alpha and beta: {a120.corr(b120):.2f}")
print(f"rolling beta from {b120.min():.2f} to {b120.max():.2f}; "
      f"alpha from {a120.min():.0f} to {a120.max():.0f}")
```
```output
correlation between rolling alpha and beta: -1.00
rolling beta from 0.45 to 1.96; alpha from -4 to 3
```

On changes (the regression of $\Delta y$ on $\Delta x$, without a constant) the constant disappears and the error
falls a lot, with the U shape typical of any rolling estimate. With short windows β tracks the changes well but is
estimated with little data: **variance**. With long windows it is estimated with a lot of data but mixes βs from
months ago: **bias**. The minimum, 0.070 with 250 sessions (0.075 with 120), is the best tradeoff for this pair;
another pair with a different pace of change would have a different one.

The underlying flaw of the window is that it treats the data point from 249 days ago the same as yesterday's, and
gives no weight at all to the one from 251 days ago. The Kalman filter replaces that rectangular window with an update
whose weights come from a model.

## The state space model

A state space model separates what is not seen (the **state**, here $\beta_t$) from what is observed (the changes in
the log prices, that is, the returns). For the model on changes:

$$
\begin{aligned}
\beta_t &= \beta_{t-1} + \eta_t, & \eta_t &\sim N(0, q) \quad \text{(state)}\\
\Delta y_t &= \beta_t \,\Delta x_t + \varepsilon_t, & \varepsilon_t &\sim N(0, r) \quad \text{(observation)}
\end{aligned}
$$

The first equation says that today's β is yesterday's plus a random change; the second, that the change in `y` is β
times that of `x` plus noise. The two variances are the only decisions: $q$, how much β can move from one day to the
next; $r$, how much noise each observation has. In this data the truth is $q = 0.006^2 = 3.6 \times 10^{-5}$; the
observation noise is not exactly white (it includes the change in the spread and the effect of the change in β on the
price), so $r$ is an approximation: the model is useful even if it is not exact. As an order of magnitude, $r$ is the
daily variance of what the hedge does not explain: with a spread of deviation 0.025 and a half-life of 10, on the
order of $10^{-4}$, which is the value taken below.

The Kalman filter computes, for each day, the best estimate of $\beta_t$ with the data up to that day and its variance
$P_t$. It does so in two steps:

```mermaid
flowchart TB
  A["Yesterday: beta, P"] --> B[Prediction]
  B --> C{"Price arrives"}
  C --> D[Innovation e]
  D --> E[Gain K]
  E --> F["Today: beta, P"]
  F --> B
```

- **Prediction.** Without new data, the best estimate of $\beta_t$ is yesterday's, but with more uncertainty:
  $\beta_{t|t-1} = \beta_{t-1|t-1}$ and $P_{t|t-1} = P_{t-1|t-1} + q$.
- **Update.** $\Delta y_t$ arrives. The **innovation** is the surprise, $e_t = \Delta y_t - \Delta x_t\,\beta_{t|t-1}$,
  with variance $S_t = \Delta x_t^2 P_{t|t-1} + r$. The **gain** $K_t = P_{t|t-1}\,\Delta x_t / S_t$ says what part of
  the surprise is attributed to a change in β, and $\beta_{t|t} = \beta_{t|t-1} + K_t e_t$, $P_{t|t} = (1 - K_t \Delta x_t)
  P_{t|t-1}$.

The gain sums it all up: if the state is very uncertain ($P$ large) or the observation very informative ($r$ small,
$\Delta x_t$ large), β moves a lot with the surprise; if not, little. With a very uncertain initial estimate, what
determines the behavior is the ratio $q/r$.

## The filter from scratch

The filter fits in a loop. This version accepts a state of any dimension, because later it will be used with two
($\alpha$ and $\beta$):

```python
def kalman(z, H, q, r, p0=1e4):
    """Kalman filter with a random walk state.

    z: observations (n,); H: design row for each day (n, k);
    q: state variances (k,); r: observation variance.
    Returns the filtered state (n, k), the innovation and its variance.
    """
    n, k = H.shape
    est, P, Q = np.zeros(k), np.eye(k) * p0, np.diag(q)
    filtrado = np.full((n, k), np.nan)
    innov, var_innov = np.full(n, np.nan), np.full(n, np.nan)
    for t in range(n):
        P = P + Q                       # prediction
        if not np.isfinite(z[t]):
            continue
        h = H[t]
        S = h @ P @ h + r
        K = P @ h / S                   # gain
        e = z[t] - h @ est              # innovation
        est = est + K * e               # update
        P = P - np.outer(K, h) @ P
        filtrado[t], innov[t], var_innov[t] = est, e, S
    return filtrado, innov, var_innov

dz, h = dy.to_numpy(), dx.to_numpy()[:, None]
r = 1e-4
filas = {}
for q in [1e-7, 1e-6, 1e-5, 3.6e-5, 1e-4, 3e-4, 1e-3, 1e-2]:
    f, _, _ = kalman(dz, np.nan_to_num(h), [q], r)
    filas[f"{q:g}"] = rmse(pd.Series(f[:, 0], index=y.index))
pd.Series(filas, name="rmse").rename_axis("q").round(4)
```
```output
q
1e-07      0.1164
1e-06      0.0984
1e-05      0.0699
3.6e-05    0.0614
0.0001     0.0641
0.0003     0.0782
0.001      0.1057
0.01       0.1912
Name: rmse, dtype: float64
```

*Fig. — β error of the Kalman filter on price changes against the state variance q (r = 1e-4): a U between bias and variance.* Curve on a log scale of q, from 1e-7 to 1e-2. The error falls from 0.116 at q = 1e-7 to a minimum of 0.0614 at the true q, 3.6e-5, highlighted, and rises to 0.191 at q = 1e-2. A dashed line marks the error of the best rolling regression, 0.0704, and a vertical one the true q.

The error as a function of q has the same U shape as that of the rolling window, and for the same reason. With a very
small q ($10^{-7}$) the filter barely lets β move and always lags behind the truth (0.116): bias. With a very large q
($10^{-2}$) every price surprise moves β and the estimate is noise (0.191): variance. The minimum, 0.0614, falls at the
true q, $3.6 \times 10^{-5}$, and improves on the best rolling window (0.0704); with $q = 10^{-4}$ the error barely
rises (0.0641): the bottom of the U is flat. The filter is optimal when the model and its variances are the true ones,
and here, with a reasonable r, they almost are.

<aside class="callout" data-kind="nota">

In the pairs trading literature it is common to write $q = \delta / (1 - \delta)$ with $\delta$ between $10^{-5}$ and
$10^{-3}$, and to take $r$ as a small constant. It is the same idea under another name: $\delta$ sets the ratio $q/r$.
There is no universal value; it depends on how fast the relationship changes and on the scale of the prices. The
honest way to choose it is the one that follows: by likelihood or out of sample, never by looking at the backtest
result.

</aside>

## Bias and variance, in plain sight

The table sums up the error, but it does not show how each filter goes wrong. One year of the β estimated with a q
about thirty times smaller and another about thirty times larger than the one at the minimum shows it:

```python
tramo = slice("2023-01-01", "2023-12-31")
ver = pd.DataFrame({"verdad": beta[tramo]})
for q in [1e-6, 1e-3]:
    f, _, _ = kalman(dz, np.nan_to_num(h), [q], r)
    est = pd.Series(f[:, 0], index=y.index).shift(1)[tramo]
    ver[f"q={q:g}"] = est
cambio = ver.diff().abs().mean()
print("Mean daily change of beta:", cambio.round(4).to_dict())
ver.iloc[[0, 60, 120, 180, 240]].round(3)
```
```output
Mean daily change of beta: {'verdad': 0.0047, 'q=1e-06': 0.0006, 'q=0.001': 0.0209}
            verdad  q=1e-06  q=0.001
fecha
2023-01-03   1.342    1.189    1.394
2023-03-30   1.356    1.189    1.284
2023-06-27   1.463    1.200    1.261
2023-09-21   1.435    1.223    1.473
2023-12-15   1.400    1.236    1.438
```

*Fig. — Pair D β in 2023 with two extreme state variances against the true value: top, small q (bias); bottom, large q (variance).* Two panels from January to December 2023, scale from 1.0 to 1.7. The true β, highlighted, goes from 1.34 to 1.46 in June and ends at 1.40. Top, the filter with q = 1e-6 lags below it, from 1.19 to 1.24, barely moving from one day to the next. Bottom, with q = 1e-3, it follows the truth on average but swings every day, between 1.07 and 1.60.

In 2023 the true β goes from 1.34 to 1.46 in June and ends the year at 1.40. With $q = 10^{-6}$ the filter lags far
behind and below (1.19 in January, 1.24 in December) and barely moves from one day to the next: a mean daily change of
0.0006, versus 0.0047 for the truth. With $q = 10^{-3}$ it follows the truth on average, but moves 0.0209 a day, more
than four times as much, and can be at 1.26 when the truth is at 1.46. A hedge ratio that jumps every day forces you to
redo the hedge every
day and pays commissions and bid-ask spread on every adjustment (lesson 2.4.7). The q that is optimal for the error is
not necessarily optimal for the P&L: a little more bias can turn out cheaper.

## Maximum likelihood with statsmodels

With real data there is no true β to measure the error against, so q and r cannot be chosen by looking at the RMSE.
What there is is the **likelihood**: the filter produces an innovation every day with its variance $S_t$, and if the
model and its variances are good, the innovations behave like independent normals with those variances. Maximizing
that likelihood estimates q and r from the data itself.

`statsmodels` has a state space framework with filter, smoother and maximum likelihood. For a time-varying β a small
class is enough: the “design matrix” changes every day (it is $\Delta x_t$) and the two variances are the parameters.

```python
class BetaVariable(sm.tsa.statespace.MLEModel):
    """dy_t = beta_t * dx_t + e_t; beta_t a random walk."""

    def __init__(self, dy, dx):
        super().__init__(dy, k_states=1, initialization="diffuse")
        self["design"] = np.asarray(dx, float).reshape(1, 1, -1)
        self["transition"] = np.eye(1)
        self["selection"] = np.eye(1)

    @property
    def param_names(self):
        return ["r", "q"]

    @property
    def start_params(self):
        return [1e-4, 1e-5]

    def transform_params(self, u):      # positive variances
        return u ** 2

    def untransform_params(self, c):
        return c ** 0.5

    def update(self, params, **kwargs):
        params = super().update(params, **kwargs)
        self["obs_cov", 0, 0] = params[0]
        self["state_cov", 0, 0] = params[1]

# arrays without dates: the session index has no fixed frequency
# and statsmodels would warn (ValueWarning) without anything changing
modelo = BetaVariable(dy.iloc[1:].to_numpy(), dx.iloc[1:].to_numpy())
ajuste = modelo.fit(disp=False)
r_mv, q_mv = ajuste.params
print(f"r = {r_mv:.2e}  q = {q_mv:.2e}  "
      f"(true q {0.006 ** 2:.1e})")
```
```output
r = 8.92e-05  q = 5.12e-05  (true q 3.6e-05)
```

Maximum likelihood gives $q = 5.1 \times 10^{-5}$, close to the true one, and $r = 8.9 \times 10^{-5}$, close to the
$10^{-4}$ that had been taken as an order of magnitude. With real data it would be the reasonable starting point,
without having looked at the truth.

Two more checks: that our own filter and that of `statsmodels` give the same result with the same parameters, and
what happens with **smoothing**:

```python
f_mv, _, _ = kalman(dz, np.nan_to_num(h), [q_mv], r_mv)
propio = pd.Series(f_mv[:, 0], index=y.index)
sm_filt = pd.Series(ajuste.filtered_state[0], index=y.index[1:])
sm_suav = pd.Series(ajuste.smoothed_state[0], index=y.index[1:])

dif = (propio.iloc[500:] - sm_filt.iloc[499:]).abs().max()
print(f"max difference own - statsmodels: {dif:.1e}")
print(f"RMSE filtered (previous day): {rmse(sm_filt):.4f}")
e_suav = (sm_suav - beta).iloc[500:]
print(f"RMSE smoothed (same day):     "
      f"{np.sqrt((e_suav ** 2).mean()):.4f}")
```
```output
max difference own - statsmodels: 1.3e-08
RMSE filtered (previous day): 0.0614
RMSE smoothed (same day):     0.0382
```

Our own filter and that of `statsmodels` agree to the eighth decimal place (the difference comes from the start:
`statsmodels` uses an exact diffuse initialization and ours, an initial variance of 10,000), and the β filtered with
the maximum likelihood variances has an error of 0.061, like the best in the table, without having looked at the
truth. The smoothed one drops to 0.038, and it is a trap: the Kalman smoother (Rauch-Tung-Striebel) runs backward
through the sample from the end, so each day's smoothed β uses all the **later** prices. It is useful for studying
what the relationship was like; in a backtest it is look-ahead bias, the same as the smoothed probabilities of an HMM
(lesson 2.1.7).

## On prices: an α that also moves and a risk for the spread

The model on changes throws away the information in the levels. If the log prices are assumed to be linked by
$\ln y_t = \alpha + \beta_t \ln x_t + u_t$, the state can be the vector $(\alpha, \beta_t)$ and the observation, the
log price:

$$
\ln y_t = \begin{pmatrix} 1 & \ln x_t \end{pmatrix} \begin{pmatrix} \alpha_t \\ \beta_t \end{pmatrix} + \varepsilon_t .
$$

It is the most widespread setup in pairs trading. The `kalman` function works as is with a two-column design and, as
observation noise, the variance of the spread ($0.025^2$). With a fixed α ($q_\alpha = 0$), which is what is usually
done, and several values of $q_\beta$, we measure the error of β and something that matters for trading the spread:
the autocorrelation of the innovation, which is the spread the filter sees.

```python
Hp = np.column_stack([np.ones(len(x)), x.to_numpy()])
filas = {}
for qb in [1e-7, 1e-6, 1e-5, 1e-4, 1e-3, 1e-2]:
    f, e, S = kalman(y.to_numpy(), Hp, [0.0, qb], 0.025 ** 2)
    inn = pd.Series(e, index=y.index).iloc[500:]
    filas[f"{qb:g}"] = {
        "rmse_beta": rmse(pd.Series(f[:, 1], index=y.index)),
                 "alfa_final": f[-1, 0],
                 "ac1_innovacion": inn.autocorr()}
# a_t = alfa + (1 - beta_t) * ln x0, with alfa = 0.3 (pares2_verdad.csv)
print("true a on the last day:",
      round(0.3 + (1 - beta.iloc[-1]) * x.iloc[0], 3))
pd.DataFrame(filas).T.rename_axis("q_beta").round(3)
```
```output
true a on the last day: -1.063
        rmse_beta  alfa_final  ac1_innovacion
q_beta
1e-07       0.131      -0.458           0.891
1e-06       0.136      -0.439           0.791
1e-05       0.120      -0.475           0.545
0.0001      0.112      -0.498           0.182
0.001       0.110      -0.503           0.009
0.01        0.110      -0.503          -0.018
```

*Fig. — Pair D hedge ratio error with each method against the true β (previous day's estimate except for the smoothed one).* Horizontal bars of root mean squared error: Full-sample OLS, 0.166; Rolling OLS on prices (250), 0.150; Rolling OLS on changes (250), 0.070; Kalman on changes (ML), 0.061; Smoothed: look-ahead, 0.038; Kalman on prices (q 1e-4), 0.112. The Kalman filter on changes with maximum-likelihood variances, highlighted, is the best look-ahead-free method; the smoothed one improves because it looks into the future, and the price model with a fixed α, misspecified for this pair, does worse than the methods on changes.

In this data the price model estimates β **worse** than the model on changes: the error does not go below 0.110,
versus 0.061 for the filter on changes. The reason is specification. Here the true constant $a_t$ moves with β (on
the last day it is −1.063) and the model assumes it fixed; the filter estimates it at −0.44 to −0.50 and offsets the
error with β, just like the rolling regression. A model in levels is only more accurate if its form is the right one.

The last column shows the other risk, which exists even when the model is well specified. With a small q the
innovation has an autocorrelation of 0.89: it keeps the memory of the spread, which is what reverts. With
$q_\beta = 10^{-4}$ it already falls to 0.18, and with more, to zero or slightly below: the filter explains every
deviation of the spread as a change in β, and the residual stops reverting. **The filter eats the spread.** For hedging
it may work; for trading the reversion, the signal is gone.

That is why in a pair q is not chosen only by the error of β (which, moreover, is not known with real data): it is
chosen by also checking whether the resulting spread keeps a useful half-life (lesson 2.1.5), checking that the form of
the model is reasonable (here, letting α move, or using the model on changes, which does not need it) and validating
out of sample.

```mermaid
flowchart TB
  Q{"State q"} -- small --> L[Slow beta: bias]
  Q -- large --> R[Noisy beta]
  R --> C[Higher costs]
  R --> S[Spread absorbed]
  L --> P[Late hedge]
```

## Common mistakes with the Kalman filter

- **Using the state updated with today's price to trade today.** The usable β is yesterday's. In the price model, the
  innovation $e_t$ is already computed with yesterday's state: it is the look-ahead-free version of the spread.
- **Using the smoothed estimate.** It is always more accurate and it is always look-ahead.
- **Choosing q and r by the backtest result.** It is overfitting under another name. Use maximum likelihood or
  out-of-sample validation (lesson 2.5.2).
- **Forgetting the start-up.** With a large initial variance, the first few dozen sessions are noise. Discard them, as
  the first 500 are discarded here when measuring the error.
- **Trusting the model.** The observation noise in this data is not white and the filter works anyway, but its
  variances $S_t$ are not exact: a z-score $e_t / \sqrt{S_t}$ used as a signal (lesson 2.3.3) inherits that error.
- **Prices on different scales.** q depends on the units: if `x` is multiplied by 100, β is divided by 100 and its q
  by 10,000. Rescaling the prices changes the appropriate q.

## Exercises

1. Repeat the table of the error as a function of q for the model on changes with $r = 2.5 \times 10^{-5}$ and with
   $r = 4 \times 10^{-4}$. Does the minimum shift?
   By how much?

**Solution.**

```python
for r_alt in [2.5e-5, 4e-4]:
    errores = {}
    for q in [2.25e-6, 9e-6, 3.6e-5, 1.44e-4, 5.76e-4]:
        f, _, _ = kalman(dz, np.nan_to_num(h), [q], r_alt)
        errores[q] = rmse(pd.Series(f[:, 0], index=y.index))
    mejor = min(errores, key=errores.get)
    print(f"r = {r_alt}: best q = {mejor:g}, "
          f"rmse {errores[mejor]:.4f}")
```
```output
r = 2.5e-05: best q = 9e-06, rmse 0.0614
r = 0.0004: best q = 0.000144, rmse 0.0614
```

When r is divided by 4 or multiplied by 4, the q at the minimum shifts in the same proportion (from $3.6 \times
10^{-5}$ with $r = 10^{-4}$ to $9 \times 10^{-6}$ with $r = 2.5 \times 10^{-5}$ and to $1.44 \times 10^{-4}$ with
$r = 4 \times 10^{-4}$) and the minimum error is the same, 0.0614:
what the filter uses is the ratio $q/r$, not each variance separately. That is why a single parameter, like the
$\delta$ of the literature, is enough.

2. With the model on changes and the maximum likelihood variances, compute the z-score of the innovation,
   $e_t / \sqrt{S_t}$. If the model were correct, what mean, standard deviation and autocorrelation should it have?
   Does it have them?

**Solution.**

```python
_, e_mv, S_mv = kalman(dz, np.nan_to_num(h), [q_mv], r_mv)
z = pd.Series(e_mv / np.sqrt(S_mv), index=y.index).iloc[500:]
print(f"mean {z.mean():.3f}  std {z.std():.3f}  "
      f"autocorrelation {z.autocorr():.3f}")
print(f"excess kurtosis {z.kurt():.2f}  autocorrelation "
      f"of z squared {(z ** 2).autocorr():.3f}")
```
```output
mean 0.004  std 1.010  autocorrelation -0.025
excess kurtosis 0.10  autocorrelation of z squared 0.003
```

With a correct model the standardized innovations would be normal white noise: mean 0, standard deviation 1, no
autocorrelation in either the level or the square, and no excess kurtosis. The five figures come out very close
(0.004; 1.01; −0.025; 0.10; 0.003): maximum likelihood fits the variances for that, and in logarithms the observation
noise does not change size with the price. And yet the model is not the true one: the noise contains the daily change
of a mean-reverting spread, whose theoretical autocorrelation is $-(1 - \varphi)/2 \approx -0.033$ with a half-life of
10, on the order of the measured −0.025. Diagnostic tests have little power against small specification errors:
passing them does not prove the model is correct.

3. A colleague proposes choosing the $\delta$ of the filter for their pair by trying twenty values and keeping the one
   that gives the best Sharpe in the 2015 to 2024 backtest. What problems does this have and what would you do
   instead?

**Solution.**

Choosing the parameter by the backtest result is a search with twenty attempts over the same sample: the best Sharpe
includes the luck of that sample (lessons 2.1.8 and 2.5.4) and the chosen δ need not be good afterward. Instead:
estimate q and r by maximum likelihood with data prior to the test period; check that the resulting spread keeps a
useful half-life and that the hedge ratio does not jump so much that rebalancing costs eat the edge; and, if you want
to compare several δ, do it in a walk-forward scheme (lesson 2.5.2), declaring in advance how many are tried.

## Self-assessment

1. In the Kalman filter for a β that follows a random walk, what does the state variance q control relative to the observation noise variance r?

   - A. How much the filter trusts each new observation: with a high q/r, β moves a lot with each data point; with a low q/r, little
   - B. The number of sessions discarded at the start
   - C. The mean of the spread
   - D. Nothing if both are positive

   *Answer A: How much the filter trusts each new observation: with a high q/r, β moves a lot with each data point; with a low q/r, little.* The Kalman gain grows with the uncertainty of the state, which accumulates q at each prediction, and decreases with r. With a diffuse initial estimate, what matters is the ratio: multiplying q and r by the same constant leaves the filtered β unchanged.

2. The rolling regression of the log prices of D (with a constant) gives an error against the true β of 0.15 to 0.53 depending on the window. Why so much?

   - A. Because the data has errors
   - B. Because pandas computes the rolling covariance incorrectly
   - C. Because the true β is constant
   - D. Because in a short window `x` varies little relative to its level, and many combinations of α and β fit almost equally well: the two estimates offset each other

   *Answer D: Because in a short window `x` varies little relative to its level, and many combinations of α and β fit almost equally well: the two estimates offset each other.* Over 120 sessions the log price of PZDX moves with a standard deviation of about 5%, around a level of 4: a higher β with a lower α explains almost the same. In the data, the rolling estimates of α and β have a correlation of −1.00 and β ranges from 0.45 to 1.96.

3. Why is the error of β measured with the previous day's estimate (`shift(1)`)?

   - A. By statistical convention
   - B. Because today's hedge is decided with what was known at yesterday's close: the estimate that already includes today's price was not available when trading
   - C. So that the error comes out smaller
   - D. Because the filter lags by one day

   *Answer B: Because today's hedge is decided with what was known at yesterday's close: the estimate that already includes today's price was not available when trading.* It is the same look-ahead rule used throughout the course. The estimate updated with today's price is for tomorrow.

4. The smoothed β from `statsmodels` has an error of 0.038 and the filtered one, 0.061. Which would you use in a backtest?

   - A. The smoothed one, because it is more accurate
   - B. Either, because both come from the same model
   - C. The filtered one: each day's smoothed estimate uses all the later observations, which did not exist on that day
   - D. The average of the two

   *Answer C: The filtered one: each day's smoothed estimate uses all the later observations, which did not exist on that day.* The smoother runs backward through the sample from the end. Its accuracy advantage is precisely the future information: in a backtest it is look-ahead bias.

5. In the price model (state α and β), raising q barely changes the error of β (from 0.13 to 0.11), but the autocorrelation of the innovation falls from 0.89 to almost 0. What does it imply for a pair?

   - A. That the model is perfect
   - B. That the filter attributes the moves of the spread to β: the residual stops reverting and the mean-reversion signal disappears
   - C. That you have to trade more
   - D. That the spread is a random walk

   *Answer B: That the filter attributes the moves of the spread to β: the residual stops reverting and the mean-reversion signal disappears.* With a very flexible β, any deviation of the spread is explained as a change in β. For hedging a position it may work; for trading the reversion of the spread (lesson 2.3.3), the filter would have eaten the signal.

## Key points

- A time-varying hedge ratio is a hidden state: β_t = β_(t−1) + noise with variance q, and the price observes it with noise of variance r.
- The Kalman filter alternates prediction (uncertainty grows by q) and update (the innovation corrects β in proportion to the gain).
- In pair D, the rolling regression on prices fails (error of 0.15 to 0.53 against the true β); on changes, the best window gives 0.070 and the well-tuned filter, 0.061.
- Choosing q is choosing between bias (slow β, small q) and variance (noisy β, large q); maximum likelihood gives a q close to the best one.
- Only the previous day's filtered β is valid. The smoothed one looks into the future; the price model with a fixed α is misspecified here, and with a large q it eats the spread.

## Resources

### Research

- [Stochastic Processes - Random Walks](https://pythonparatrading.com/en/stochastic-processes/) · Research-0005 · 11 min

### Documentation

- [statsmodels: time series analysis by state space methods (statespace)](https://www.statsmodels.org/stable/statespace.html) · statsmodels.org
- [statsmodels: custom state space models (MLEModel)](https://www.statsmodels.org/stable/examples/notebooks/generated/statespace_custom_models.html) · statsmodels.org

### Books

- [Kalman, R. E. (1960): A New Approach to Linear Filtering and Prediction Problems. Journal of Basic Engineering 82(1), 35-45](https://doi.org/10.1115/1.3662552) · doi.org

### Data

- [Level 2 data: true β of pair D, session by session (CSV)](https://pythonparatrading.com/curso/datos/nivel-2/pares2_beta_real.csv) · Download · pares2_beta_real.csv
- [Level 2 data: eight synthetic pairs, 2014-2025 (Parquet)](https://pythonparatrading.com/curso/datos/nivel-2/pares2.parquet) · Download · pares2.parquet

Previous lesson: [2.1.5 · Cointegration: Engle-Granger, Johansen and the half-life of mean reversion](https://pythonparatrading.com/en/python-trading-course/level-2/cointegration-engle-granger-johansen-half-life/)

Next lesson: [2.1.7 · Regime detection: hidden Markov models and structural breaks](https://pythonparatrading.com/en/python-trading-course/level-2/regime-detection-hmm-structural-breaks/)