Better filter

docs/api/dists/index.md

@@ -5,6 +5,8 @@ return types of the
 [StateSpaceModel](../ta/kalman.md#quantflow.ta.kalman.StateSpaceModel)
 distribution-level interface.
 
+::: quantflow.dists.MeanAndCov
+
 ::: quantflow.dists.Distribution
 
 ::: quantflow.dists.MvNormal

docs/api/rates/index.md

@@ -6,7 +6,7 @@ The central concept is the [discount factor](../../glossary.md#discount-factor)
 
 **[Rate](interest_rate.md)** represents a spot or forward interest rate with a chosen compounding frequency (continuous by default) and day count convention. It supports continuous and periodic compounding and can be bootstrapped directly from a spot/forward pair.
 
-**[YieldCurve](yield_curve.md)** is the abstract base for term-structure models. It defines the interface via `discount_factor` and `instantaneous_forward_rate`, with the two quantities linked by
+**[YieldCurve](yield_curve.md)** is the abstract base for term-structure models. It defines the interface via [discount_factor][quantflow.rates.yield_curve.YieldCurve.discount_factor] and [instantaneous_forward_rate][quantflow.rates.yield_curve.YieldCurve.instantaneous_forward_rate], with the two quantities linked by
 
 \begin{equation}
     f(\tau) = -\frac{\partial \ln D_\tau}{\partial \tau}
@@ -18,4 +18,4 @@ The central concept is the [discount factor](../../glossary.md#discount-factor)
 
 **[VasicekCurve](vasicek.md)** is a Gaussian mean-reverting short-rate model with analytical formulas for discount factors and instantaneous forward rates.
 
-**[Options Discounting](options.md)** provides `YieldCurveCalibration`, the base class for fitting a yield curve to discount factors, and `OptionsDiscountingCalibration`, which bootstraps asset and quote curves from put-call parity observations.
+**[Calibration](calibration.md)** provides [YieldCurveCalibration][quantflow.rates.calibration.YieldCurveCalibration], the base class for fitting a yield curve to discount factors, and [OptionsDiscountingCalibration][quantflow.rates.calibration.OptionsDiscountingCalibration], which bootstraps asset and quote curves from put-call parity observations.

docs/api/ta/kalman.md

@@ -15,8 +15,6 @@ A state-space model separates a time series into two layers:
 See the [kalman filter](../../theory/kalman.md) theory page for more details
 on the algorithms and their applications.
 
-::: quantflow.ta.kalman.MeanAndCov
-
 ::: quantflow.ta.kalman.StateSpaceModel
 
 ::: quantflow.ta.kalman.LinearGaussianModel

docs/examples/rates_kalman.py

@@ -1,5 +1,6 @@
 import asyncio
 from datetime import timedelta
+from decimal import Decimal
 
 import numpy as np
 import pandas as pd
@@ -9,6 +10,7 @@
 from docs.examples._utils import assets_path, cached_df
 from quantflow.data.fed import FederalReserve
 from quantflow.rates.calibration import tenor_to_years
+from quantflow.rates.cir import CIRCurve
 from quantflow.rates.vasicek import VasicekCurve
 
 
@@ -25,50 +27,56 @@ async def fetch() -> pd.DataFrame:
 df = fed_yield_curves()
 # weekly panel (uniform 7-day step) using the average yield over each week
 weekly = df.resample("W-WED").mean().dropna()
+ttm = np.array([tenor_to_years(c) for c in weekly.columns])
 
-# calibrate the Vasicek short-rate model to the panel by Kalman-filter MLE
-calibrator = VasicekCurve().calibrator()
-curve = calibrator.calibrate_historical_rates_dataframe(weekly, frequency=2)
-print(curve.model_dump_json(indent=2, exclude={"ref_date"}))
+# Vasicek: linear-Gaussian, fitted with the exact Kalman filter
+vasicek_cal = VasicekCurve().calibrator()
+vasicek = vasicek_cal.calibrate_historical_rates_dataframe(weekly, frequency=2)
+print("Vasicek (Kalman filter)")
+print(vasicek.model_dump_json(indent=2, exclude={"ref_date"}))
 
-# rebuild model-implied yields from the filtered short rate path: y = (B r - A) / tau
-ttm = np.array([tenor_to_years(c) for c in weekly.columns])
-a, b = curve.affine_coefficients(ttm)
-A, B = np.asarray(a, dtype=float), np.asarray(b, dtype=float)
-short_rate = calibrator.filtered_short_rate  # one value per observation date
+# CIR: state-dependent variance, fitted with the unscented Kalman filter
+cir_cal = CIRCurve().calibrator()
+cir = cir_cal.calibrate_historical_rates_dataframe(weekly, frequency=2)
+print("CIR (unscented Kalman filter)")
+print(cir.model_dump_json(indent=2, exclude={"ref_date"}))
+
+# model-implied semi-annual rates at each date from the filtered short rate paths
+va_short = vasicek_cal.filtered_short_rate
+cir_short = cir_cal.filtered_short_rate
+va_model = np.zeros((len(va_short), len(ttm)))
+cir_model = np.zeros((len(cir_short), len(ttm)))
+for t in range(len(va_short)):
+    vasicek.rate = Decimal(str(float(va_short[t])))
+    va_model[t] = np.asarray(vasicek.rates(ttm), dtype=float)
+    cir.rate = Decimal(str(float(cir_short[t])))
+    cir_model[t] = np.asarray(cir.rates(ttm), dtype=float)
 
-# observed (par -> continuous) and model yields per tenor, over time
+# observed (par -> continuous) and both model yields per tenor, over time
 tenors = ["1Y", "2Y", "5Y", "10Y"]
+colours = dict(observed="#636efa", Vasicek="#ef553b", CIR="#00cc96")
 fig = make_subplots(rows=2, cols=2, subplot_titles=tenors)
 for k, tenor in enumerate(tenors):
     i = weekly.columns.get_loc(tenor)
-    observed = 2.0 * np.log1p(weekly[tenor].to_numpy() / 2.0) * 100
-    model = (B[i] * short_rate - A[i]) / ttm[i] * 100
     row, col = k // 2 + 1, k % 2 + 1
-    fig.add_trace(
-        go.Scatter(
-            x=weekly.index,
-            y=observed,
-            name="observed",
-            legendgroup="observed",
-            showlegend=k == 0,
-            line=dict(color="#636efa"),
-        ),
-        row=row,
-        col=col,
-    )
-    fig.add_trace(
-        go.Scatter(
-            x=weekly.index,
-            y=model,
-            name="model",
-            legendgroup="model",
-            showlegend=k == 0,
-            line=dict(color="#ef553b"),
-        ),
-        row=row,
-        col=col,
-    )
-fig.update_layout(title="Observed vs Vasicek model yields")
+    series = {
+        "observed": 2.0 * np.log1p(weekly[tenor].to_numpy() / 2.0) * 100,
+        "Vasicek": va_model[:, i] * 100,
+        "CIR": cir_model[:, i] * 100,
+    }
+    for name, values in series.items():
+        fig.add_trace(
+            go.Scatter(
+                x=weekly.index,
+                y=values,
+                name=name,
+                legendgroup=name,
+                showlegend=k == 0,
+                line=dict(color=colours[name]),
+            ),
+            row=row,
+            col=col,
+        )
+fig.update_layout(title="Observed vs Vasicek and CIR model yields")
 fig.update_yaxes(title_text="yield (%)")
 fig.write_image(assets_path("rates_kalman.png"), width=1600, height=800)

docs/theory/kalman.md

@@ -239,9 +239,11 @@ without user intervention.
 When the transition $f(x, \Delta t)$ is non-linear the Kalman predict step
 is no longer exact. The unscented Kalman filter (UKF), introduced by
 [Julier & Uhlmann (1997)](../bibliography.md#julier_uhlmann), replaces it
-with a **sigma-point** propagation. The observation update is unchanged: it
-reuses the Kalman update step described above, so the UKF inherits the
-Sherman-Morrison optimisation when applicable.
+with a **sigma-point** propagation. The observation update follows the same
+structure, but builds the innovation covariance $S_t$ and the cross covariance
+$C_t$ from the propagated sigma points and solves the gain $K_t = C_t S_t^{-1}$
+with a dense Cholesky factorisation. The Sherman-Morrison fast path is specific
+to the exact linear filter and is not used here.
 
 ### Sigma-Point Predict
 

docs/tutorials/index.md

@@ -9,4 +9,4 @@ Step-by-step guides for common quantflow workflows.
 | [Heston Volatility Model](heston_calibration.md) | Calibrate the Heston and Heston-jump-diffusion models to an implied volatility surface |
 | [SPX Volatility Surface](spx_vol_surface.md) | Build a 3D implied volatility surface for the S&P 500 from a Yahoo Finance option chain |
 | [BNS Volatility Model](bns_calibration.md) | Calibrate the Barndorff-Nielsen and Shephard stochastic-volatility model to an implied volatility surface |
-| [Vasicek Calibration from Rates](rates_kalman.md) | Fit the Vasicek short-rate model to historical Treasury rates by maximum likelihood with a Kalman filter |
+| [Yield Curve Calibration from Rates](rates_kalman.md) | Fit the Vasicek (Kalman filter) and CIR (unscented Kalman filter) short-rate models to historical Treasury rates by maximum likelihood |

docs/tutorials/rates_kalman.md

@@ -1,27 +1,28 @@
-# Vasicek Calibration from Rates
-
-This tutorial calibrates the [VasicekCurve][quantflow.rates.vasicek.VasicekCurve]
-short-rate model to a panel of historical US Treasury yields by maximum
-likelihood, using a [KalmanFilter][quantflow.ta.kalman.KalmanFilter] to evaluate
-the likelihood. It shows the full workflow: pulling the data from the Federal
-Reserve, reshaping it into a uniform panel, fitting the model, and comparing the
-fitted curve with the observed yields.
-
-For the mechanics of the filter itself, see the
+# Yield Curve Calibration from Rates
+
+This tutorial calibrates two short-rate models, the
+[VasicekCurve][quantflow.rates.vasicek.VasicekCurve] and the
+[CIRCurve][quantflow.rates.cir.CIRCurve], to a panel of historical US Treasury
+yields by maximum likelihood. Both treat the short rate as a latent state and
+the yields as noisy observations of it, but they need different filters: Vasicek
+is linear-Gaussian and uses the exact
+[KalmanFilter][quantflow.ta.kalman.KalmanFilter], while CIR has a state-dependent
+variance and uses the
+[UnscentedKalmanFilter][quantflow.ta.kalman.UnscentedKalmanFilter].
+
+For the mechanics of the filters themselves, see the
 [Kalman Filter](../theory/kalman.md) theory page.
 
 ## The idea
 
-The Vasicek short rate is an Ornstein-Uhlenbeck process. We never observe it
-directly: what we observe is a cross section of yields at each date. Treating
-the short rate as a latent state and the yields as noisy linear observations of
-it turns calibration into a state-space estimation problem.
+A short-rate model never observes its state directly: what we observe at each
+date is a cross section of yields. Treating the short rate as a latent state and
+the yields as noisy observations of it turns calibration into a state-space
+estimation problem.
 
-Over a uniform time step the dynamics reduce to a Gaussian AR(1) and each yield
-is affine in the short rate. The
-[calibrate_historical_rates][quantflow.rates.vasicek.VasicekCurveCalibration.calibrate_historical_rates]
-method documents these equations. The Kalman filter computes the exact Gaussian
-log-likelihood of the observed panel, and the calibrator maximises it over
+At each date every yield is affine in the short rate, so the observation is
+linear in both models. The filter computes the Gaussian log-likelihood of the
+observed panel, and the calibrator maximises it over
 $(\kappa, \theta, \sigma, h)$, where $h$ is the observation noise standard
 deviation.
 
@@ -36,15 +37,50 @@ The calibration assumes a uniform time step, while the raw data is sampled on
 business days. Resampling to weekly Wednesdays with the average yield over each
 week gives an evenly spaced panel.
 
+## Vasicek: the exact Kalman filter
+
+The Vasicek short rate is an Ornstein-Uhlenbeck process. Over a uniform time step
+its dynamics reduce to a Gaussian AR(1) with a constant innovation variance, and
+each yield is affine in the short rate. The
+[calibrate_historical_rates][quantflow.rates.vasicek.VasicekCurveCalibration.calibrate_historical_rates]
+method documents these equations.
+
+Because the model is fully linear-Gaussian, the exact Kalman filter gives the
+exact log-likelihood. One full Kalman pass over the panel is performed per
+optimiser iteration.
+
+## CIR: why the unscented filter
+
+The CIR short rate is a square-root diffusion,
+$dr_t = \kappa(\theta - r_t)\,dt + \sigma\sqrt{r_t}\,dW_t$. Its conditional mean
+is still linear in the previous state, but its conditional variance depends on
+the state:
+
+\begin{equation}
+    \text{Var}[r_t \mid r_{t-1}] =
+        r_{t-1}\,\frac{\sigma^2}{\kappa}\left(\phi - \phi^2\right)
+        + \theta\,\frac{\sigma^2}{2\kappa}\left(1 - \phi\right)^2,
+        \quad \phi = e^{-\kappa \Delta t}.
+\end{equation}
+
+The exact Kalman filter assumes a constant process-noise covariance, so it
+cannot represent this heteroskedasticity. The unscented Kalman filter only needs
+the conditional mean and covariance of the transition, which it propagates
+through sigma points, so the state-dependent variance drops straight in. The
+[CIRStateSpaceModel][quantflow.rates.cir.CIRStateSpaceModel] supplies those
+moments and the affine observation, and
+[calibrate_historical_rates][quantflow.rates.cir.CIRCurveCalibration.calibrate_historical_rates]
+runs the unscented filter inside the same maximum-likelihood loop.
+
 ## Calibrating
 
+For both models,
 [calibrate_historical_rates_dataframe][quantflow.rates.calibration.YieldCurveCalibration.calibrate_historical_rates_dataframe]
 parses the tenor columns into times to maturity, infers the time step from the
 index, converts the par yields to continuously compounded rates (here
 `frequency=2` for semiannual compounding), and runs the maximum-likelihood fit.
-One full Kalman pass over the panel is performed per optimiser iteration.
 
-The fitted parameters and the final filtered short rate are returned on the
+The fitted parameters and the final filtered short rate are returned on each
 calibrated curve:
 
 ```
@@ -53,18 +89,20 @@ calibrated curve:
 
 ## Model versus observed through time
 
-The fit is a time-series fit, so the right check is whether the model tracks the
-history of each tenor. The calibrator exposes the
-[filtered_short_rate][quantflow.rates.vasicek.VasicekCurveCalibration.filtered_short_rate]
-path, from which each tenor's model-implied yield is reconstructed through the
-affine yield relation and plotted against its observed history.
+The fit is a time-series fit, so the right check is whether each model tracks the
+history of every tenor. Both calibrators expose a `filtered_short_rate` path
+([Vasicek][quantflow.rates.vasicek.VasicekCurveCalibration.filtered_short_rate],
+[CIR][quantflow.rates.cir.CIRCurveCalibration.filtered_short_rate]), from which
+each tenor's model-implied yield is reconstructed through the affine yield
+relation and plotted against its observed history.
 
-The single factor tracks the short and intermediate tenors (1Y, 2Y) closely, but
-the model yields are smoother than the data at the long end (5Y, 10Y): one
-mean-reverting factor cannot capture the independent variation of the long end,
-the expected limitation of the one-factor Vasicek model.
+Both single-factor models track the short and intermediate tenors (1Y, 2Y)
+closely and are smoother than the data at the long end (5Y, 10Y): one
+mean-reverting factor cannot capture the independent variation of the long end.
+The two fits stay close in this rate regime, the difference being the CIR
+volatility that scales with the level of rates.
 
-[![Observed vs Vasicek model yields](../assets/examples/rates_kalman.png)](../assets/examples/rates_kalman.png){target="_blank"}
+[![Observed vs Vasicek and CIR model yields](../assets/examples/rates_kalman.png)](../assets/examples/rates_kalman.png){target="_blank"}
 
 ## Code
 

mkdocs.yml

@@ -124,8 +124,8 @@ nav:
       - Option Pricing: tutorials/option_pricing.md
       - Pricing Method Comparison: tutorials/pricing_method_comparison.md
       - SPX Volatility Surface: tutorials/spx_vol_surface.md
-      - Vasicek Calibration from Rates: tutorials/rates_kalman.md
       - Volatility Surface: tutorials/volatility_surface.md
+      - Yield Curve Calibration from Rates: tutorials/rates_kalman.md
   - Theory:
       - theory/index.md
       - Characteristic Function: theory/characteristic.md